Skip to content

refactor Moyal distribution #4704

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 2 commits into from
May 16, 2021
Merged
Show file tree
Hide file tree
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
61 changes: 23 additions & 38 deletions pymc3/distributions/continuous.py
Original file line number Diff line number Diff line change
Expand Up @@ -3959,6 +3959,21 @@ def _distr_parameters_for_repr(self):
return []


class MoyalRV(RandomVariable):
name = "moyal"
ndim_supp = 0
ndims_params = [0, 0]
dtype = "floatX"
_print_name = ("Moyal", "\\operatorname{Moyal}")

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


moyal = MoyalRV()


class Moyal(Continuous):
r"""
Moyal log-likelihood.
Expand Down Expand Up @@ -4006,43 +4021,18 @@ class Moyal(Continuous):
sigma: float
Scale parameter (sigma > 0).
"""
rv_op = moyal

def __init__(self, mu=0, sigma=1.0, *args, **kwargs):
self.mu = at.as_tensor_variable(floatX(mu))
self.sigma = at.as_tensor_variable(floatX(sigma))
@classmethod
def dist(cls, mu=0, sigma=1.0, *args, **kwargs):
mu = at.as_tensor_variable(floatX(mu))
sigma = at.as_tensor_variable(floatX(sigma))

assert_negative_support(sigma, "sigma", "Moyal")

self.mean = self.mu + self.sigma * (np.euler_gamma + at.log(2))
self.median = self.mu - self.sigma * at.log(2 * at.erfcinv(1 / 2) ** 2)
self.mode = self.mu
self.variance = (np.pi ** 2 / 2.0) * self.sigma ** 2

super().__init__(*args, **kwargs)

def random(self, point=None, size=None):
"""
Draw random values from Moyal 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.sigma], point=point, size=size)
# return generate_samples(
# stats.moyal.rvs, loc=mu, scale=sigma, dist_shape=self.shape, size=size
# )
return super().dist([mu, sigma], *args, **kwargs)

def logp(self, value):
def logp(value, mu, sigma):
"""
Calculate log-probability of Moyal distribution at specified value.

Expand All @@ -4056,15 +4046,13 @@ def logp(self, value):
-------
TensorVariable
"""
mu = self.mu
sigma = self.sigma
scaled = (value - mu) / sigma
return bound(
(-(1 / 2) * (scaled + at.exp(-scaled)) - at.log(sigma) - (1 / 2) * at.log(2 * np.pi)),
0 < sigma,
)

def logcdf(self, value):
def logcdf(value, mu, sigma):
"""
Compute the log of the cumulative distribution function for Moyal distribution
at the specified value.
Expand All @@ -4079,9 +4067,6 @@ def logcdf(self, value):
-------
TensorVariable
"""
mu = self.mu
sigma = self.sigma

scaled = (value - mu) / sigma
return bound(
at.log(at.erfc(at.exp(-scaled / 2) * (2 ** -0.5))),
Expand Down
2 changes: 0 additions & 2 deletions pymc3/tests/test_distributions.py
Original file line number Diff line number Diff line change
Expand Up @@ -2547,7 +2547,6 @@ def test_rice(self):
lambda value, b, sigma: sp.rice.logpdf(value, b=b, loc=0, scale=sigma),
)

@pytest.mark.xfail(reason="Distribution not refactored yet")
def test_moyal_logp(self):
# Using a custom domain, because the standard `R` domain undeflows with scipy in float64
value_domain = Domain([-inf, -1.5, -1, -0.01, 0.0, 0.01, 1, 1.5, inf])
Expand All @@ -2558,7 +2557,6 @@ def test_moyal_logp(self):
lambda value, mu, sigma: floatX(sp.moyal.logpdf(value, mu, sigma)),
)

@pytest.mark.xfail(reason="Distribution not refactored yet")
@pytest.mark.xfail(
condition=(aesara.config.floatX == "float32"),
reason="Pymc3 underflows earlier than scipy on float32",
Expand Down
20 changes: 13 additions & 7 deletions pymc3/tests/test_distributions_random.py
Original file line number Diff line number Diff line change
Expand Up @@ -319,12 +319,6 @@ class TestZeroInflatedBinomial(BaseTestCases.BaseTestCase):
params = {"n": 10, "p": 0.6, "psi": 0.3}


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


class BaseTestDistribution(SeededTest):
pymc_dist: Optional[Callable] = None
pymc_dist_params = dict()
Expand Down Expand Up @@ -491,6 +485,19 @@ class TestStudentT(BaseTestDistribution):
]


class TestMoyal(BaseTestDistribution):
pymc_dist = pm.Moyal
pymc_dist_params = {"mu": 0.0, "sigma": 1.0}
expected_rv_op_params = {"mu": 0.0, "sigma": 1.0}
reference_dist_params = {"loc": 0.0, "scale": 1.0}
reference_dist = seeded_scipy_distribution_builder("moyal")
tests_to_run = [
"check_pymc_params_match_rv_op",
"check_pymc_draws_match_reference",
"check_rv_size",
]


class TestStudentTLam(BaseTestDistribution):
pymc_dist = pm.StudentT
lam, sigma = get_tau_sigma(tau=2.0)
Expand Down Expand Up @@ -1391,7 +1398,6 @@ def ref_rand(size, mu, sigma):

pymc3_random(pm.LogitNormal, {"mu": R, "sigma": Rplus}, ref_rand=ref_rand)

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