started splitting up dist. classes; some small testval test changes

pymc3/distributions/continuous.py

@@ -9,35 +9,54 @@
 
 import numpy as np
 import theano.tensor as tt
+from theano.gof.op import get_test_value
 from scipy import stats
 
 from . import transforms
 from .dist_math import bound, logpow, gammaln, betaln, std_cdf, i0, i1
-from .distribution import UnivariateContinuous, draw_values, generate_samples
+from .distribution import Univariate, Continuous, draw_values, generate_samples
 
 __all__ = ['Uniform', 'Flat', 'Normal', 'Beta', 'Exponential', 'Laplace',
            'StudentT', 'Cauchy', 'HalfCauchy', 'Gamma', 'Weibull',
            'Bound', 'StudentTpos', 'Lognormal', 'ChiSquared', 'HalfNormal',
            'Wald', 'Pareto', 'InverseGamma', 'ExGaussian', 'VonMises']
 
 
-class PositiveUnivariateContinuous(UnivariateContinuous):
+class PositiveUnivariateContinuous(Univariate, Continuous):
     """Base class for positive univariate continuous distributions"""
 
     def __init__(self, *args, **kwargs):
         transform = kwargs.get('transform', transforms.log)
         super(PositiveUnivariateContinuous, self).__init__(transform=transform,
                                                            *args, **kwargs)
+        # TODO: is there a better way to use Theano's
+        # existing `...tag.test_value` mechanics?
+        ndim_sum = self.ndim_supp + self.ndim_ind + self.ndim_reps
+        if self.testval is None:
+            if ndim_sum == 0:
+                self.testval = 0.5
+            else:
+                self.testval = get_test_value(tt.alloc(*((0.5,) + self.shape)))
 
 
-class UnitUnivariateContinuous(UnivariateContinuous):
+class UnitUnivariateContinuous(Univariate, Continuous):
     """Base class for univariate continuous distributions in [0,1]"""
 
     def __init__(self, *args, **kwargs):
         transform = kwargs.get('transform', transforms.logodds)
         super(UnitUnivariateContinuous, self).__init__(transform=transform,
                                                        *args, **kwargs)
 
+        # TODO: is there a better way to use Theano's
+        # existing `...tag.test_value` mechanics?
+        ndim_sum = self.ndim_supp + self.ndim_ind + self.ndim_reps
+        if self.testval is None:
+            if ndim_sum == 0:
+                self.testval = 0.5
+            else:
+                self.testval = get_test_value(tt.alloc(*((0.5,) + self.shape)))
+
+
 
 def get_tau_sd(tau=None, sd=None):
     """
@@ -80,7 +99,7 @@ def get_tau_sd(tau=None, sd=None):
     return (tt.as_tensor_variable(tau), tt.as_tensor_variable(sd))
 
 
-class Uniform(UnivariateContinuous):
+class Uniform(Univariate, Continuous):
     R"""
     Continuous uniform log-likelihood.
 
@@ -134,7 +153,7 @@ def logp(self, value):
                      value >= lower, value <= upper)
 
 
-class Flat(UnivariateContinuous):
+class Flat(Univariate, Continuous):
     """
     Uninformative log-likelihood that returns 0 regardless of
     the passed value.
@@ -155,7 +174,7 @@ def logp(self, value):
         return tt.zeros_like(value)
 
 
-class Normal(UnivariateContinuous):
+class Normal(Univariate, Continuous):
     R"""
     Univariate normal log-likelihood.
 
@@ -518,7 +537,7 @@ def logp(self, value):
         return bound(tt.log(lam) - lam * value, value > 0, lam > 0)
 
 
-class Laplace(UnivariateContinuous):
+class Laplace(Univariate, Continuous):
     R"""
     Laplace log-likelihood.
 
@@ -628,7 +647,7 @@ def logp(self, value):
                      tt.log(value), tau > 0)
 
 
-class StudentT(UnivariateContinuous):
+class StudentT(Univariate, Continuous):
     r"""
     Non-central Student's T log-likelihood.
 
@@ -754,7 +773,7 @@ def logp(self, value):
                      value >= m, alpha > 0, m > 0)
 
 
-class Cauchy(UnivariateContinuous):
+class Cauchy(Univariate, Continuous):
     R"""
     Cauchy log-likelihood.
 
@@ -1083,7 +1102,7 @@ def logp(self, value):
                      value >= 0, alpha > 0, beta > 0)
 
 
-class Bounded(UnivariateContinuous):
+class Bounded(Univariate, Continuous):
     """A bounded distribution."""
 
     def __init__(self, distribution, lower, upper, *args, **kwargs):
@@ -1141,8 +1160,7 @@ def dist(self, *args, **kwargs):
 
 StudentTpos = Bound(StudentT, 0)
 
-
-class ExGaussian(UnivariateContinuous):
+class ExGaussian(Univariate, Continuous):
     R"""
     Exponentially modified Gaussian log-likelihood.
 
@@ -1227,7 +1245,7 @@ def logp(self, value):
         return bound(lp, sigma > 0., nu > 0.)
 
 
-class VonMises(UnivariateContinuous):
+class VonMises(Univariate, Continuous):
     R"""
     Univariate VonMises log-likelihood.
     .. math::

pymc3/distributions/discrete.py

@@ -6,14 +6,14 @@
 from scipy import stats
 
 from .dist_math import bound, factln, binomln, betaln, logpow
-from .distribution import UnivariateDiscrete, draw_values, generate_samples
+from .distribution import Univariate, Discrete, draw_values, generate_samples
 
 __all__ = ['Binomial',  'BetaBinomial',  'Bernoulli',  'Poisson',
            'NegativeBinomial', 'ConstantDist', 'ZeroInflatedPoisson',
            'DiscreteUniform', 'Geometric', 'Categorical']
 
 
-class Binomial(UnivariateDiscrete):
+class Binomial(Univariate, Discrete):
     R"""
     Binomial log-likelihood.
 
@@ -63,7 +63,7 @@ def logp(self, value):
             0 <= p, p <= 1)
 
 
-class BetaBinomial(UnivariateDiscrete):
+class BetaBinomial(Univariate, Discrete):
     R"""
     Beta-binomial log-likelihood.
 
@@ -132,7 +132,7 @@ def logp(self, value):
                      alpha > 0, beta > 0)
 
 
-class Bernoulli(UnivariateDiscrete):
+class Bernoulli(Univariate, Discrete):
     R"""Bernoulli log-likelihood
 
     The Bernoulli distribution describes the probability of successes
@@ -175,7 +175,7 @@ def logp(self, value):
             p >= 0, p <= 1)
 
 
-class Poisson(UnivariateDiscrete):
+class Poisson(Univariate, Discrete):
     R"""
     Poisson log-likelihood.
 
@@ -223,7 +223,7 @@ def logp(self, value):
             mu >= 0, value >= 0)
 
 
-class NegativeBinomial(UnivariateDiscrete):
+class NegativeBinomial(Univariate, Discrete):
     R"""
     Negative binomial log-likelihood.
 
@@ -280,7 +280,7 @@ def logp(self, value):
                          negbinom)
 
 
-class Geometric(UnivariateDiscrete):
+class Geometric(Univariate, Discrete):
     R"""
     Geometric log-likelihood.
 
@@ -321,7 +321,7 @@ def logp(self, value):
                      0 <= p, p <= 1, value >= 1)
 
 
-class DiscreteUniform(UnivariateDiscrete):
+class DiscreteUniform(Univariate, Discrete):
     R"""
     Discrete uniform distribution.
 
@@ -371,7 +371,7 @@ def logp(self, value):
                      lower <= value, value <= upper)
 
 
-class Categorical(UnivariateDiscrete):
+class Categorical(Univariate, Discrete):
     R"""
     Categorical log-likelihood.
 
@@ -423,7 +423,7 @@ def logp(self, value):
                      sumto1)
 
 
-class ConstantDist(UnivariateDiscrete):
+class ConstantDist(Univariate, Discrete):
     """
     Constant log-likelihood.
 
@@ -454,7 +454,7 @@ def logp(self, value):
         return bound(0, tt.eq(value, c))
 
 
-class ZeroInflatedPoisson(UnivariateDiscrete):
+class ZeroInflatedPoisson(Univariate, Discrete):
 
     def __init__(self, theta, z, ndim=None, size=None, dtype=None, *args, **kwargs):
         self.theta = tt.as_tensor_variable(theta)

pymc3/distributions/distribution.py

@@ -41,32 +41,6 @@ def _as_tensor_shape_variable(var):
     return res
 
 
-def _as_tensor_shape_variable(var):
-    """ Just a collection of useful shape stuff from
-    `_infer_ndim_bcast` """
-
-    if var is None:
-        return T.constant([], dtype='int64')
-
-    res = var
-    if isinstance(res, (tuple, list)):
-        if len(res) == 0:
-            return T.constant([], dtype='int64')
-        res = T.as_tensor_variable(res, ndim=1)
-
-    else:
-        if res.ndim != 1:
-            raise TypeError("shape must be a vector or list of scalar, got\
-                            '%s'" % res)
-
-    if (not (res.dtype.startswith('int') or
-             res.dtype.startswith('uint'))):
-
-        raise TypeError('shape must be an integer vector or list',
-                        res.dtype)
-    return res
-
-
 class Distribution(object):
     """Statistical distribution"""
     def __new__(cls, name, *args, **kwargs):
@@ -185,23 +159,20 @@ def __init__(self, shape_supp, shape_ind, shape_reps, bcast, dtype,
                 tuple(self.shape_ind) +\
                 tuple(self.shape_supp)
 
-        if testval is None:
-            if ndim_sum == 0:
-                testval = tt.constant(0, dtype=dtype)
-            else:
-                testval = tt.zeros(self.shape)
-
         self.ndim = tt.get_vector_length(self.shape)
-
-        self.testval = testval
         self.defaults = defaults
         self.transform = transform
+
+        if testval is None:
+            testval = self.get_test_value(defaults=self.defaults)
+
+        self.testval = testval
         self.type = tt.TensorType(str(dtype), bcast)
 
     def default(self):
-        return self.get_test_val(self.testval, self.defaults)
+        return self.get_test_value(self.testval, self.defaults)
 
-    def get_test_val(self, val, defaults):
+    def get_test_value(self, val=None, defaults=None):
         if val is None:
             for v in defaults:
                 the_attr = getattr(self, v, None)
@@ -216,11 +187,13 @@ def get_test_val(self, val, defaults):
 
     def getattr_value(self, val):
         """ Attempts to obtain a non-symbolic value for an attribute
-        (potentially in str form)
+        (potentially given in str form)
         """
         if isinstance(val, str):
             val = getattr(self, val)
 
+        # Could use Theano's:
+        # val = theano.gof.op.get_test_value(val)
         if isinstance(val, tt.sharedvar.SharedVariable):
             return val.get_value()
         elif isinstance(val, tt.TensorVariable):
@@ -231,6 +204,10 @@ def getattr_value(self, val):
         return val
 
 
+def TensorType(dtype, shape):
+    return tt.TensorType(str(dtype), np.atleast_1d(shape) == 1)
+
+
 class NoDistribution(Distribution):
 
     def __init__(self, shape_supp, shape_ind, shape_reps, bcast, dtype,
@@ -284,7 +261,7 @@ def __init__(self, logp, ndim_support, ndim, size, bcast, dtype='float64',
         self.logp = logp
 
 
-class UnivariateContinuous(Continuous):
+class Univariate(Distribution):
 
     def __init__(self, dist_params, ndim=None, size=None, dtype=None,
                  bcast=None, *args, **kwargs):
@@ -307,29 +284,15 @@ def __init__(self, dist_params, ndim=None, size=None, dtype=None,
                                    tuple(x.dtype for x in dist_params))
 
         # We just assume
-        super(UnivariateContinuous, self).__init__(
+        super(Univariate, self).__init__(
             tuple(), tuple(), size, bcast, *args, **kwargs)
 
 
-class MultivariateContinuous(Continuous):
-
-    pass
-
-
-class MultivariateDiscrete(Discrete):
+class Multivariate(Distribution):
 
     pass
 
 
-class UnivariateDiscrete(Discrete):
-
-    def __init__(self, ndim, size, bcast, *args, **kwargs):
-        self.shape_supp = ()
-
-        super(UnivariateDiscrete, self).__init__(
-            0, ndim, size, bcast, *args, **kwargs)
-
-
 def draw_values(params, point=None):
     """
     Draw (fix) parameter values. Handles a number of cases: