kumaraswamy.py 3.2 KB

12345678910111213141516171819202122232425262728293031323334353637383940414243444546474849505152535455565758596061626364656667686970717273747576
  1. import torch
  2. from torch._six import nan
  3. from torch.distributions import constraints
  4. from torch.distributions.uniform import Uniform
  5. from torch.distributions.transformed_distribution import TransformedDistribution
  6. from torch.distributions.transforms import AffineTransform, PowerTransform
  7. from torch.distributions.utils import broadcast_all, euler_constant
  8. def _moments(a, b, n):
  9. """
  10. Computes nth moment of Kumaraswamy using using torch.lgamma
  11. """
  12. arg1 = 1 + n / a
  13. log_value = torch.lgamma(arg1) + torch.lgamma(b) - torch.lgamma(arg1 + b)
  14. return b * torch.exp(log_value)
  15. class Kumaraswamy(TransformedDistribution):
  16. r"""
  17. Samples from a Kumaraswamy distribution.
  18. Example::
  19. >>> m = Kumaraswamy(torch.tensor([1.0]), torch.tensor([1.0]))
  20. >>> m.sample() # sample from a Kumaraswamy distribution with concentration alpha=1 and beta=1
  21. tensor([ 0.1729])
  22. Args:
  23. concentration1 (float or Tensor): 1st concentration parameter of the distribution
  24. (often referred to as alpha)
  25. concentration0 (float or Tensor): 2nd concentration parameter of the distribution
  26. (often referred to as beta)
  27. """
  28. arg_constraints = {'concentration1': constraints.positive, 'concentration0': constraints.positive}
  29. support = constraints.unit_interval
  30. has_rsample = True
  31. def __init__(self, concentration1, concentration0, validate_args=None):
  32. self.concentration1, self.concentration0 = broadcast_all(concentration1, concentration0)
  33. finfo = torch.finfo(self.concentration0.dtype)
  34. base_dist = Uniform(torch.full_like(self.concentration0, 0),
  35. torch.full_like(self.concentration0, 1),
  36. validate_args=validate_args)
  37. transforms = [PowerTransform(exponent=self.concentration0.reciprocal()),
  38. AffineTransform(loc=1., scale=-1.),
  39. PowerTransform(exponent=self.concentration1.reciprocal())]
  40. super(Kumaraswamy, self).__init__(base_dist, transforms, validate_args=validate_args)
  41. def expand(self, batch_shape, _instance=None):
  42. new = self._get_checked_instance(Kumaraswamy, _instance)
  43. new.concentration1 = self.concentration1.expand(batch_shape)
  44. new.concentration0 = self.concentration0.expand(batch_shape)
  45. return super(Kumaraswamy, self).expand(batch_shape, _instance=new)
  46. @property
  47. def mean(self):
  48. return _moments(self.concentration1, self.concentration0, 1)
  49. @property
  50. def mode(self):
  51. # Evaluate in log-space for numerical stability.
  52. log_mode = self.concentration0.reciprocal() * \
  53. (-self.concentration0).log1p() - (-self.concentration0 * self.concentration1).log1p()
  54. log_mode[(self.concentration0 < 1) | (self.concentration1 < 1)] = nan
  55. return log_mode.exp()
  56. @property
  57. def variance(self):
  58. return _moments(self.concentration1, self.concentration0, 2) - torch.pow(self.mean, 2)
  59. def entropy(self):
  60. t1 = (1 - self.concentration1.reciprocal())
  61. t0 = (1 - self.concentration0.reciprocal())
  62. H0 = torch.digamma(self.concentration0 + 1) + euler_constant
  63. return t0 + t1 * H0 - torch.log(self.concentration1) - torch.log(self.concentration0)