diff --git a/RELEASE-NOTES.md b/RELEASE-NOTES.md index 35ae968c54..ebe696e824 100644 --- a/RELEASE-NOTES.md +++ b/RELEASE-NOTES.md @@ -29,6 +29,7 @@ It also brings some dreadfully awaited fixes, so be sure to go through the chang - Fixed mathematical formulation in `MvStudentT` random method. (see [#4359](https://github.com/pymc-devs/pymc3/pull/4359)) - Fix issue in `logp` method of `HyperGeometric`. It now returns `-inf` for invalid parameters (see [4367](https://github.com/pymc-devs/pymc3/pull/4367)) - Fixed `MatrixNormal` random method to work with parameters as random variables. (see [#4368](https://github.com/pymc-devs/pymc3/pull/4368)) +- Update the `logcdf` method of several continuous distributions to return -inf for invalid parameters and values, and raise an informative error when multiple values cannot be evaluated in a single call. (see [4393](https://github.com/pymc-devs/pymc3/pull/4393)) ## PyMC3 3.10.0 (7 December 2020) diff --git a/pymc3/distributions/continuous.py b/pymc3/distributions/continuous.py index 08102d2022..5653c2fd15 100644 --- a/pymc3/distributions/continuous.py +++ b/pymc3/distributions/continuous.py @@ -278,7 +278,7 @@ def logcdf(self, value): Parameters ---------- - value: numeric + value: numeric or np.ndarray or theano.tensor Value(s) for which log CDF is calculated. If the log CDF for multiple values are desired the values must be provided in a numpy array or theano tensor. @@ -286,13 +286,16 @@ def logcdf(self, value): ------- TensorVariable """ + lower = self.lower + upper = self.upper + return tt.switch( - tt.or_(tt.lt(value, self.lower), tt.gt(value, self.upper)), + tt.lt(value, lower) | tt.lt(upper, lower), -np.inf, tt.switch( - tt.eq(value, self.upper), + tt.lt(value, upper), + tt.log(value - lower) - tt.log(upper - lower), 0, - tt.log(value - self.lower) - tt.log(self.upper - self.lower), ), ) @@ -344,7 +347,7 @@ def logcdf(self, value): Parameters ---------- - value: numeric + value: numeric or np.ndarray or theano.tensor Value(s) for which log CDF is calculated. If the log CDF for multiple values are desired the values must be provided in a numpy array or theano tensor. @@ -401,7 +404,7 @@ def logcdf(self, value): Parameters ---------- - value: numeric + value: numeric or np.ndarray or theano.tensor Value(s) for which log CDF is calculated. If the log CDF for multiple values are desired the values must be provided in a numpy array or theano tensor. @@ -542,7 +545,7 @@ def logcdf(self, value): Parameters ---------- - value: numeric + value: numeric or np.ndarray or theano.tensor Value(s) for which log CDF is calculated. If the log CDF for multiple values are desired the values must be provided in a numpy array or theano tensor. @@ -900,7 +903,7 @@ def logcdf(self, value): Parameters ---------- - value: numeric + value: numeric or np.ndarray or theano.tensor Value(s) for which log CDF is calculated. If the log CDF for multiple values are desired the values must be provided in a numpy array or theano tensor. @@ -910,10 +913,10 @@ def logcdf(self, value): """ sigma = self.sigma z = zvalue(value, mu=0, sigma=sigma) - return tt.switch( - tt.lt(z, -1.0), - tt.log(tt.erfcx(-z / tt.sqrt(2.0))) - tt.sqr(z), + return bound( tt.log1p(-tt.erfc(z / tt.sqrt(2.0))), + 0 <= value, + 0 < sigma, ) @@ -1106,7 +1109,7 @@ def logcdf(self, value): Parameters ---------- - value: numeric + value: numeric or np.ndarray or theano.tensor Value(s) for which log CDF is calculated. If the log CDF for multiple values are desired the values must be provided in a numpy array or theano tensor. @@ -1297,20 +1300,30 @@ def logcdf(self, value): Parameters ---------- value: numeric - Value(s) for which log CDF is calculated. If the log CDF for multiple - values are desired the values must be provided in a numpy array or theano tensor. + Value(s) for which log CDF is calculated. Returns ------- TensorVariable """ - value = floatX(tt.as_tensor(value)) - a = floatX(tt.as_tensor(self.alpha)) - b = floatX(tt.as_tensor(self.beta)) - return tt.switch( - tt.le(value, 0), - -np.inf, - tt.switch(tt.ge(value, 1), 0, tt.log(incomplete_beta(a, b, value))), + # incomplete_beta function can only handle scalar values (see #4342) + if np.ndim(value): + raise TypeError( + f"Beta.logcdf expects a scalar value but received a {np.ndim(value)}-dimensional object." + ) + + a = self.alpha + b = self.beta + + return bound( + tt.switch( + tt.lt(value, 1), + tt.log(incomplete_beta(a, b, value)), + 0, + ), + 0 <= value, + 0 < a, + 0 < b, ) def _distr_parameters_for_repr(self): @@ -1521,7 +1534,7 @@ def logcdf(self, value): Parameters ---------- - value: numeric + value: numeric or np.ndarray or theano.tensor Value(s) for which log CDF is calculated. If the log CDF for multiple values are desired the values must be provided in a numpy array or theano tensor. @@ -1636,7 +1649,7 @@ def logcdf(self, value): Parameters ---------- - value: numeric + value: numeric or np.ndarray or theano.tensor Value(s) for which log CDF is calculated. If the log CDF for multiple values are desired the values must be provided in a numpy array or theano tensor. @@ -1792,7 +1805,7 @@ def logcdf(self, value): Parameters ---------- - value: numeric + value: numeric or np.ndarray or theano.tensor Value(s) for which log CDF is calculated. If the log CDF for multiple values are desired the values must be provided in a numpy array or theano tensor. @@ -1955,20 +1968,32 @@ def logcdf(self, value): Parameters ---------- value: numeric - Value(s) for which log CDF is calculated. If the log CDF for multiple - values are desired the values must be provided in a numpy array or theano tensor. + Value(s) for which log CDF is calculated. Returns ------- TensorVariable """ + # incomplete_beta function can only handle scalar values (see #4342) + if np.ndim(value): + raise TypeError( + f"StudentT.logcdf expects a scalar value but received a {np.ndim(value)}-dimensional object." + ) + nu = self.nu mu = self.mu sigma = self.sigma + lam = self.lam t = (value - mu) / sigma sqrt_t2_nu = tt.sqrt(t ** 2 + nu) x = (t + sqrt_t2_nu) / (2.0 * sqrt_t2_nu) - return tt.log(incomplete_beta(nu / 2.0, nu / 2.0, x)) + + return bound( + tt.log(incomplete_beta(nu / 2.0, nu / 2.0, x)), + 0 < nu, + 0 < sigma, + 0 < lam, + ) class Pareto(Continuous): @@ -2090,7 +2115,7 @@ def logcdf(self, value): Parameters ---------- - value: numeric + value: numeric or np.ndarray or theano.tensor Value(s) for which log CDF is calculated. If the log CDF for multiple values are desired the values must be provided in a numpy array or theano tensor. @@ -2209,7 +2234,7 @@ def logcdf(self, value): Parameters ---------- - value: numeric + value: numeric or np.ndarray or theano.tensor Value(s) for which log CDF is calculated. If the log CDF for multiple values are desired the values must be provided in a numpy array or theano tensor. @@ -2317,7 +2342,7 @@ def logcdf(self, value): Parameters ---------- - value: numeric + value: numeric or np.ndarray or theano.tensor Value(s) for which log CDF is calculated. If the log CDF for multiple values are desired the values must be provided in a numpy array or theano tensor. @@ -2468,7 +2493,7 @@ def logcdf(self, value): Parameters ---------- - value: numeric + value: numeric or np.ndarray or theano.tensor Value(s) for which log CDF is calculated. If the log CDF for multiple values are desired the values must be provided in a numpy array or theano tensor. @@ -2478,7 +2503,17 @@ def logcdf(self, value): """ alpha = self.alpha beta = self.beta - return bound(tt.log(tt.gammainc(alpha, beta * value)), value >= 0, alpha > 0, beta > 0) + # Avoid C-assertion when the gammainc function is called with invalid values (#4340) + safe_alpha = tt.switch(tt.lt(alpha, 0), 0, alpha) + safe_beta = tt.switch(tt.lt(beta, 0), 0, beta) + safe_value = tt.switch(tt.lt(value, 0), 0, value) + + return bound( + tt.log(tt.gammainc(safe_alpha, safe_beta * safe_value)), + 0 <= value, + 0 < alpha, + 0 < beta, + ) def _distr_parameters_for_repr(self): return ["alpha", "beta"] @@ -2632,7 +2667,7 @@ def logcdf(self, value): Parameters ---------- - value: numeric + value: numeric or np.ndarray or theano.tensor Value(s) for which log CDF is calculated. If the log CDF for multiple values are desired the values must be provided in a numpy array or theano tensor. @@ -2642,11 +2677,16 @@ def logcdf(self, value): """ alpha = self.alpha beta = self.beta + # Avoid C-assertion when the gammaincc function is called with invalid values (#4340) + safe_alpha = tt.switch(tt.lt(alpha, 0), 0, alpha) + safe_beta = tt.switch(tt.lt(beta, 0), 0, beta) + safe_value = tt.switch(tt.lt(value, 0), 0, value) + return bound( - tt.log(tt.gammaincc(alpha, beta / value)), - value >= 0, - alpha > 0, - beta > 0, + tt.log(tt.gammaincc(safe_alpha, safe_beta / safe_value)), + 0 <= value, + 0 < alpha, + 0 < beta, ) @@ -2814,7 +2854,7 @@ def logcdf(self, value): Parameters ---------- - value: numeric + value: numeric or np.ndarray or theano.tensor Value(s) for which log CDF is calculated. If the log CDF for multiple values are desired the values must be provided in a numpy array or theano tensor. @@ -3114,7 +3154,7 @@ def logcdf(self, value): Parameters ---------- - value: numeric + value: numeric or np.ndarray or theano.tensor Value(s) for which log CDF is calculated. If the log CDF for multiple values are desired the values must be provided in a numpy array or theano tensor. @@ -3503,7 +3543,7 @@ def logcdf(self, value): Parameters ---------- - value: numeric + value: numeric or np.ndarray or theano.tensor Value(s) for which log CDF is calculated. If the log CDF for multiple values are desired the values must be provided in a numpy array or theano tensor. @@ -3632,7 +3672,7 @@ def logcdf(self, value): Parameters ---------- - value: numeric + value: numeric or np.ndarray or theano.tensor Value(s) for which log CDF is calculated. If the log CDF for multiple values are desired the values must be provided in a numpy array or theano tensor. @@ -3914,7 +3954,7 @@ def logcdf(self, value): Parameters ---------- - value: numeric + value: numeric or np.ndarray or theano.tensor Value(s) for which log CDF is calculated. If the log CDF for multiple values are desired the values must be provided in a numpy array or theano tensor. @@ -4256,7 +4296,7 @@ def logcdf(self, value): Parameters ---------- - value: numeric + value: numeric or np.ndarray or theano.tensor Value(s) for which log CDF is calculated. If the log CDF for multiple values are desired the values must be provided in a numpy array or theano tensor. diff --git a/pymc3/distributions/discrete.py b/pymc3/distributions/discrete.py index d1dc858f9c..1710b12440 100644 --- a/pymc3/distributions/discrete.py +++ b/pymc3/distributions/discrete.py @@ -166,9 +166,7 @@ def logcdf(self, value): # incomplete_beta function can only handle scalar values (see #4342) if np.ndim(value): raise TypeError( - "Binomial.logcdf expects a scalar value but received a {}-dimensional object.".format( - np.ndim(value) - ) + f"Binomial.logcdf expects a scalar value but received a {np.ndim(value)}-dimensional object." ) n = self.n @@ -336,9 +334,7 @@ def logcdf(self, value): # logcdf can only handle scalar values at the moment if np.ndim(value): raise TypeError( - "BetaBinomial.logcdf expects a scalar value but received a {}-dimensional object.".format( - np.ndim(value) - ) + f"BetaBinomial.logcdf expects a scalar value but received a {np.ndim(value)}-dimensional object." ) alpha = self.alpha @@ -463,7 +459,7 @@ def logcdf(self, value): Parameters ---------- - value: numeric + value: numeric or np.ndarray or theano.tensor Value(s) for which log CDF is calculated. If the log CDF for multiple values are desired the values must be provided in a numpy array or theano tensor. @@ -602,7 +598,7 @@ def logcdf(self, value): Parameters ---------- - value: numeric + value: numeric or np.ndarray or theano.tensor Value(s) for which log CDF is calculated. If the log CDF for multiple values are desired the values must be provided in a numpy array or theano tensor. @@ -717,7 +713,7 @@ def logcdf(self, value): Parameters ---------- - value: numeric + value: numeric or np.ndarray or theano.tensor Value(s) for which log CDF is calculated. If the log CDF for multiple values are desired the values must be provided in a numpy array or theano tensor. @@ -727,7 +723,7 @@ def logcdf(self, value): """ mu = self.mu value = tt.floor(value) - # To avoid gammaincc C-assertion when given invalid values (#4340) + # Avoid C-assertion when the gammaincc function is called with invalid values (#4340) safe_mu = tt.switch(tt.lt(mu, 0), 0, mu) safe_value = tt.switch(tt.lt(value, 0), 0, value) @@ -910,9 +906,7 @@ def logcdf(self, value): # incomplete_beta function can only handle scalar values (see #4342) if np.ndim(value): raise TypeError( - "NegativeBinomial.logcdf expects a scalar value but received a {}-dimensional object.".format( - np.ndim(value) - ) + f"NegativeBinomial.logcdf expects a scalar value but received a {np.ndim(value)}-dimensional object." ) # TODO: avoid `p` recomputation if distribution was defined in terms of `p` @@ -1017,7 +1011,7 @@ def logcdf(self, value): Parameters ---------- - value: numeric + value: numeric or np.ndarray or theano.tensor Value(s) for which log CDF is calculated. If the log CDF for multiple values are desired the values must be provided in a numpy array or theano tensor. @@ -1166,9 +1160,7 @@ def logcdf(self, value): # logcdf can only handle scalar values at the moment if np.ndim(value): raise TypeError( - "BetaBinomial.logcdf expects a scalar value but received a {}-dimensional object.".format( - np.ndim(value) - ) + f"HyperGeometric.logcdf expects a scalar value but received a {np.ndim(value)}-dimensional object." ) # TODO: Use lower upper in locgdf for smarter logsumexp? @@ -1288,7 +1280,7 @@ def logcdf(self, value): Parameters ---------- - value: numeric + value: numeric or np.ndarray or theano.tensor Value(s) for which log CDF is calculated. If the log CDF for multiple values are desired the values must be provided in a numpy array or theano tensor. @@ -1300,7 +1292,11 @@ def logcdf(self, value): lower = self.lower return bound( - tt.log(tt.minimum(tt.floor(value), upper) - lower + 1) - tt.log(upper - lower + 1), + tt.switch( + tt.lt(value, upper), + tt.log(tt.minimum(tt.floor(value), upper) - lower + 1) - tt.log(upper - lower + 1), + 0, + ), lower <= value, lower <= upper, ) @@ -1601,7 +1597,7 @@ def logcdf(self, value): Parameters ---------- - value: numeric + value: numeric or np.ndarray or theano.tensor Value(s) for which log CDF is calculated. If the log CDF for multiple values are desired the values must be provided in a numpy array or theano tensor. @@ -1743,9 +1739,7 @@ def logcdf(self, value): # logcdf can only handle scalar values due to limitation in Binomial.logcdf if np.ndim(value): raise TypeError( - "ZeroInflatedBinomial.logcdf expects a scalar value but received a {}-dimensional object.".format( - np.ndim(value) - ) + f"ZeroInflatedBinomial.logcdf expects a scalar value but received a {np.ndim(value)}-dimensional object." ) psi = self.psi @@ -1913,9 +1907,7 @@ def logcdf(self, value): # logcdf can only handle scalar values due to limitation in NegativeBinomial.logcdf if np.ndim(value): raise TypeError( - "ZeroInflatedNegativeBinomial.logcdf expects a scalar value but received a {}-dimensional object.".format( - np.ndim(value) - ) + f"ZeroInflatedNegativeBinomial.logcdf expects a scalar value but received a {np.ndim(value)}-dimensional object." ) psi = self.psi diff --git a/pymc3/tests/test_distributions.py b/pymc3/tests/test_distributions.py index 9a21e00c15..a2d00a2d60 100644 --- a/pymc3/tests/test_distributions.py +++ b/pymc3/tests/test_distributions.py @@ -575,6 +575,35 @@ def check_logcdf( err_msg=str(pt), ) + # Test that values below domain evaluate to -np.inf + if np.isfinite(domain.lower): + below_domain = domain.lower - 1 + assert_equal( + dist.logcdf(below_domain).tag.test_value, + -np.inf, + err_msg=str(below_domain), + ) + + # Test that values above domain evaluate to 0 + # Natural domains do not have inf as the upper edge, but should also be ignored + not_nat_domain = domain not in (NatSmall, Nat, NatBig, PosNat) + if not_nat_domain and np.isfinite(domain.upper): + above_domain = domain.upper + 1 + assert_equal( + dist.logcdf(above_domain).tag.test_value, + 0, + err_msg=str(above_domain), + ) + + # Test that method works with multiple values or raises informative TypeError + try: + dist.logcdf(np.array([value, value])).tag.test_value + except TypeError as err: + if not str(err).endswith( + ".logcdf expects a scalar value but received a 1-dimensional object." + ): + raise + def check_selfconsistency_discrete_logcdf( self, distribution, domain, paramdomains, decimal=None, n_samples=100 ): @@ -681,7 +710,7 @@ def test_flat(self): with Model(): x = Flat("a") assert_allclose(x.tag.test_value, 0) - self.check_logcdf(Flat, Runif, {}, lambda value: np.log(0.5)) + self.check_logcdf(Flat, R, {}, lambda value: np.log(0.5)) # Check infinite cases individually. assert 0.0 == Flat.dist().logcdf(np.inf).tag.test_value assert -np.inf == Flat.dist().logcdf(-np.inf).tag.test_value @@ -692,7 +721,7 @@ def test_half_flat(self): x = HalfFlat("a", shape=2) assert_allclose(x.tag.test_value, 1) assert x.tag.test_value.shape == (2,) - self.check_logcdf(HalfFlat, Runif, {}, lambda value: -np.inf) + self.check_logcdf(HalfFlat, Rplus, {}, lambda value: -np.inf) # Check infinite cases individually. assert 0.0 == HalfFlat.dist().logcdf(np.inf).tag.test_value assert -np.inf == HalfFlat.dist().logcdf(-np.inf).tag.test_value