Skip to content

Refactor LogitNormal #4703

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Merged
merged 6 commits into from
May 17, 2021
Merged
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
62 changes: 26 additions & 36 deletions pymc3/distributions/continuous.py
Original file line number Diff line number Diff line change
@@ -47,6 +47,7 @@
from aesara.tensor.var import TensorVariable
from scipy import stats
from scipy.interpolate import InterpolatedUnivariateSpline
from scipy.special import expit

from pymc3.aesaraf import floatX
from pymc3.distributions import logp_transform, transforms
@@ -66,7 +67,7 @@
)
from pymc3.distributions.distribution import Continuous
from pymc3.distributions.special import log_i0
from pymc3.math import invlogit, log1mexp, log1pexp, logdiffexp, logit
from pymc3.math import log1mexp, log1pexp, logdiffexp, logit
from pymc3.util import UNSET

__all__ = [
@@ -3727,6 +3728,21 @@ def logcdf(value, mu, s):
)


class LogitNormalRV(RandomVariable):
name = "logit_normal"
ndim_supp = 0
ndims_params = [0, 0]
dtype = "floatX"
_print_name = ("logitNormal", "\\operatorname{logitNormal}")

@classmethod
def rng_fn(cls, rng, mu, sigma, size=None):
return expit(stats.norm.rvs(loc=mu, scale=sigma, size=size, random_state=rng))


logit_normal = LogitNormalRV()


class LogitNormal(UnitContinuous):
r"""
Logit-Normal log-likelihood.
@@ -3771,44 +3787,22 @@ class LogitNormal(UnitContinuous):
tau: float
Scale parameter (tau > 0).
"""
rv_op = logit_normal

def __init__(self, mu=0, sigma=None, tau=None, sd=None, **kwargs):
@classmethod
def dist(cls, mu=0, sigma=None, tau=None, sd=None, **kwargs):
if sd is not None:
sigma = sd
self.mu = mu = at.as_tensor_variable(floatX(mu))
mu = at.as_tensor_variable(floatX(mu))
tau, sigma = get_tau_sigma(tau=tau, sigma=sigma)
self.sigma = self.sd = at.as_tensor_variable(sigma)
self.tau = tau = at.as_tensor_variable(tau)

self.median = invlogit(mu)
sigma = sd = at.as_tensor_variable(sigma)
tau = at.as_tensor_variable(tau)
assert_negative_support(sigma, "sigma", "LogitNormal")
assert_negative_support(tau, "tau", "LogitNormal")

super().__init__(**kwargs)

def random(self, point=None, size=None):
"""
Draw random values from LogitNormal distribution.
Parameters
----------
point: dict, optional
Dict of variable values on which random values are to be
conditioned (uses default point if not specified).
size: int, optional
Desired size of random sample (returns one sample if not
specified).
Returns
-------
array
"""
# mu, _, sigma = draw_values([self.mu, self.tau, self.sigma], point=point, size=size)
# return expit(
# generate_samples(stats.norm.rvs, loc=mu, scale=sigma, dist_shape=self.shape, size=size)
# )
return super().dist([mu, sigma], **kwargs)

def logp(self, value):
def logp(value, mu, sigma):
"""
Calculate log-probability of LogitNormal distribution at specified value.
@@ -3822,8 +3816,7 @@ def logp(self, value):
-------
TensorVariable
"""
mu = self.mu
tau = self.tau
tau, sigma = get_tau_sigma(sigma=sigma)
return bound(
-0.5 * tau * (logit(value) - mu) ** 2
+ 0.5 * at.log(tau / (2.0 * np.pi))
@@ -3833,9 +3826,6 @@ def logp(self, value):
tau > 0,
)

def _distr_parameters_for_repr(self):
return ["mu", "sigma"]


class Interpolated(BoundedContinuous):
r"""
1 change: 0 additions & 1 deletion pymc3/tests/test_distributions.py
Original file line number Diff line number Diff line change
@@ -2515,7 +2515,6 @@ def test_logistic(self):
decimal=select_by_precision(float64=6, float32=1),
)

@pytest.mark.xfail(reason="Distribution not refactored yet")
def test_logitnormal(self):
self.check_logp(
LogitNormal,
41 changes: 30 additions & 11 deletions pymc3/tests/test_distributions_random.py
Original file line number Diff line number Diff line change
@@ -158,7 +158,7 @@ def setup_method(self, *args, **kwargs):
self.model = pm.Model()

def get_random_variable(self, shape, with_vector_params=False, name=None):
""" Creates a RandomVariable of the parametrized distribution. """
"""Creates a RandomVariable of the parametrized distribution."""
if with_vector_params:
params = {
key: value * np.ones(self.shape, dtype=np.dtype(type(value)))
@@ -187,7 +187,7 @@ def get_random_variable(self, shape, with_vector_params=False, name=None):

@staticmethod
def sample_random_variable(random_variable, size):
""" Draws samples from a RandomVariable using its .random() method. """
"""Draws samples from a RandomVariable using its .random() method."""
if size is None:
return random_variable.eval()
else:
@@ -196,7 +196,7 @@ def sample_random_variable(random_variable, size):
@pytest.mark.parametrize("size", [None, (), 1, (1,), 5, (4, 5)], ids=str)
@pytest.mark.parametrize("shape", [None, ()], ids=str)
def test_scalar_distribution_shape(self, shape, size):
""" Draws samples of different [size] from a scalar [shape] RV. """
"""Draws samples of different [size] from a scalar [shape] RV."""
rv = self.get_random_variable(shape)
exp_shape = self.default_shape if shape is None else tuple(np.atleast_1d(shape))
exp_size = self.default_size if size is None else tuple(np.atleast_1d(size))
@@ -216,7 +216,7 @@ def test_scalar_distribution_shape(self, shape, size):
"shape", [None, (), (1,), (1, 1), (1, 2), (10, 11, 1), (9, 10, 2)], ids=str
)
def test_scalar_sample_shape(self, shape, size):
""" Draws samples of scalar [size] from a [shape] RV. """
"""Draws samples of scalar [size] from a [shape] RV."""
rv = self.get_random_variable(shape)
exp_shape = self.default_shape if shape is None else tuple(np.atleast_1d(shape))
exp_size = self.default_size if size is None else tuple(np.atleast_1d(size))
@@ -301,12 +301,6 @@ class TestExGaussian(BaseTestCases.BaseTestCase):
params = {"mu": 0.0, "sigma": 1.0, "nu": 1.0}


@pytest.mark.xfail(reason="This distribution has not been refactored for v4")
class TestLogitNormal(BaseTestCases.BaseTestCase):
distribution = pm.LogitNormal
params = {"mu": 0.0, "sigma": 1.0}


@pytest.mark.xfail(reason="This distribution has not been refactored for v4")
class TestZeroInflatedNegativeBinomial(BaseTestCases.BaseTestCase):
distribution = pm.ZeroInflatedNegativeBinomial
@@ -515,6 +509,32 @@ class TestNormal(BaseTestDistribution):
]


class TestLogitNormal(BaseTestDistribution):
def logit_normal_rng_fn(self, rng, size, loc, scale):
return expit(st.norm.rvs(loc=loc, scale=scale, size=size, random_state=rng))

pymc_dist = pm.LogitNormal
pymc_dist_params = {"mu": 5.0, "sigma": 10.0}
expected_rv_op_params = {"mu": 5.0, "sigma": 10.0}
reference_dist_params = {"loc": 5.0, "scale": 10.0}
reference_dist = lambda self: functools.partial(
self.logit_normal_rng_fn, rng=self.get_random_state()
)
tests_to_run = [
"check_pymc_params_match_rv_op",
"check_pymc_draws_match_reference",
"check_rv_size",
]


class TestLogitNormalTau(BaseTestDistribution):
pymc_dist = pm.LogitNormal
tau, sigma = get_tau_sigma(tau=25.0)
pymc_dist_params = {"mu": 1.0, "tau": tau}
expected_rv_op_params = {"mu": 1.0, "sigma": sigma}
tests_to_run = ["check_pymc_params_match_rv_op"]


class TestNormalTau(BaseTestDistribution):
pymc_dist = pm.Normal
tau, sigma = get_tau_sigma(tau=25.0)
@@ -1384,7 +1404,6 @@ def test_dirichlet_multinomial_dist_ShapeError(self, n, a, shape, expectation):
with expectation:
m.random()

@pytest.mark.xfail(reason="This distribution has not been refactored for v4")
def test_logitnormal(self):
def ref_rand(size, mu, sigma):
return expit(st.norm.rvs(loc=mu, scale=sigma, size=size))