pytorch

1
import torch
2
from torch import nan
3
from torch.distributions import constraints
4
from torch.distributions.transformed_distribution import TransformedDistribution
5
from torch.distributions.transforms import AffineTransform, PowerTransform
6
from torch.distributions.uniform import Uniform
7
from torch.distributions.utils import broadcast_all, euler_constant
8

9
__all__ = ["Kumaraswamy"]
10

11

12
def _moments(a, b, n):
13
    """
14
    Computes nth moment of Kumaraswamy using using torch.lgamma
15
    """
16
    arg1 = 1 + n / a
17
    log_value = torch.lgamma(arg1) + torch.lgamma(b) - torch.lgamma(arg1 + b)
18
    return b * torch.exp(log_value)
19

20

21
class Kumaraswamy(TransformedDistribution):
22
    r"""
23
    Samples from a Kumaraswamy distribution.
24

25
    Example::
26

27
        >>> # xdoctest: +IGNORE_WANT("non-deterministic")
28
        >>> m = Kumaraswamy(torch.tensor([1.0]), torch.tensor([1.0]))
29
        >>> m.sample()  # sample from a Kumaraswamy distribution with concentration alpha=1 and beta=1
30
        tensor([ 0.1729])
31

32
    Args:
33
        concentration1 (float or Tensor): 1st concentration parameter of the distribution
34
            (often referred to as alpha)
35
        concentration0 (float or Tensor): 2nd concentration parameter of the distribution
36
            (often referred to as beta)
37
    """
38
    arg_constraints = {
39
        "concentration1": constraints.positive,
40
        "concentration0": constraints.positive,
41
    }
42
    support = constraints.unit_interval
43
    has_rsample = True
44

45
    def __init__(self, concentration1, concentration0, validate_args=None):
46
        self.concentration1, self.concentration0 = broadcast_all(
47
            concentration1, concentration0
48
        )
49
        finfo = torch.finfo(self.concentration0.dtype)
50
        base_dist = Uniform(
51
            torch.full_like(self.concentration0, 0),
52
            torch.full_like(self.concentration0, 1),
53
            validate_args=validate_args,
54
        )
55
        transforms = [
56
            PowerTransform(exponent=self.concentration0.reciprocal()),
57
            AffineTransform(loc=1.0, scale=-1.0),
58
            PowerTransform(exponent=self.concentration1.reciprocal()),
59
        ]
60
        super().__init__(base_dist, transforms, validate_args=validate_args)
61

62
    def expand(self, batch_shape, _instance=None):
63
        new = self._get_checked_instance(Kumaraswamy, _instance)
64
        new.concentration1 = self.concentration1.expand(batch_shape)
65
        new.concentration0 = self.concentration0.expand(batch_shape)
66
        return super().expand(batch_shape, _instance=new)
67

68
    @property
69
    def mean(self):
70
        return _moments(self.concentration1, self.concentration0, 1)
71

72
    @property
73
    def mode(self):
74
        # Evaluate in log-space for numerical stability.
75
        log_mode = (
76
            self.concentration0.reciprocal() * (-self.concentration0).log1p()
77
            - (-self.concentration0 * self.concentration1).log1p()
78
        )
79
        log_mode[(self.concentration0 < 1) | (self.concentration1 < 1)] = nan
80
        return log_mode.exp()
81

82
    @property
83
    def variance(self):
84
        return _moments(self.concentration1, self.concentration0, 2) - torch.pow(
85
            self.mean, 2
86
        )
87

88
    def entropy(self):
89
        t1 = 1 - self.concentration1.reciprocal()
90
        t0 = 1 - self.concentration0.reciprocal()
91
        H0 = torch.digamma(self.concentration0 + 1) + euler_constant
92
        return (
93
            t0
94
            + t1 * H0
95
            - torch.log(self.concentration1)
96
            - torch.log(self.concentration0)
97
        )
98