/usr/share/pyshared/pymc/ScipyDistributions.py is in python-pymc 2.2+ds-1.
This file is owned by root:root, with mode 0o644.
The actual contents of the file can be viewed below.
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 146 147 148 149 150 151 152 153 154 155 156 157 158 159 160 161 162 163 164 165 166 167 168 169 170 171 172 173 174 175 176 177 178 179 180 181 182 183 184 185 186 187 188 189 190 191 192 193 194 195 196 197 198 199 200 201 202 203 204 205 206 207 208 209 210 | import scipy.stats.distributions as sc_dst
import inspect
import numpy as np
from pymc import Stochastic
from copy import copy
from .distributions import *
"""
Wraps a SciPy rv object in a PyMC random variable object.
Needs to wrap the following methods:
generic.rvs(<shape(s)>,loc=0,scale=1)
- random variates
generic.pdf(x,<shape(s)>,loc=0,scale=1)
- probability density function
generic.cdf(x,<shape(s)>,loc=0,scale=1)
- cumulative density function
generic.sf(x,<shape(s)>,loc=0,scale=1)
- survival function (1-cdf --- sometimes more accurate)
generic.ppf(q,<shape(s)>,loc=0,scale=1)
- percent point function (inverse of cdf --- percentiles)
generic.isf(q,<shape(s)>,loc=0,scale=1)
- inverse survival function (inverse of sf)
generic.stats(<shape(s)>,loc=0,scale=1,moments='mv')
- mean('m',axis=0), variance('v'), skew('s'), and/or kurtosis('k')
generic.entropy(<shape(s)>,loc=0,scale=1)
- (differential) entropy of the RV.
"""
__all__ = ['stochastic_from_scipy_dist']
def separate_shape_args(kwds, shape_args):
new_kwds = {}
for name in kwds:
new_kwds[name] = kwds[name]
args = [new_kwds.pop(name) for name in shape_args]
return args, new_kwds
def remove_size_arg(fun):
def wrapper(size=1,*args, **kwds):
return fun(*args, **kwds)
def stochastic_from_scipy_dist(scipy_dist):
"""
Return a Stochastic subclass made from a particular SciPy distribution.
"""
name = scipy_dist.__class__.__name__.replace('_gen','').capitalize()
(args, varargs, varkw, defaults) = inspect.getargspec(scipy_dist._cdf)
shape_args = args[2:]
if isinstance(scipy_dist, sc_dst.rv_continuous):
dtype=float
def logp(value, **kwds):
args, kwds = separate_shape_args(kwds, shape_args)
return np.sum(scipy_dist.logpdf(value,*args,**kwds))
parent_names = shape_args + ['loc', 'scale']
defaults = [None] * (len(parent_names)-2) + [0., 1.]
elif isinstance(scipy_dist, sc_dst.rv_discrete):
dtype=int
def logp(value, **kwds):
args, kwds = separate_shape_args(kwds, shape_args)
return np.sum(scipy_dist.logpmf(value,*args,**kwds))
parent_names = shape_args + ['loc']
defaults = [None] * (len(parent_names)-1) + [0]
else:
return None
parents_default = dict(zip(parent_names, defaults))
def random(shape=None,**kwds):
args, kwds = separate_shape_args(kwds, shape_args)
if shape is None:
return scipy_dist.rvs(*args,**kwds)
else:
return np.reshape(scipy_dist.rvs(*args,**kwds), shape)
# Build docstring from distribution
docstr = name[0]+' = '+name + '(name, '+', '.join(parent_names)+', value=None, shape=None, trace=True, rseed=True, doc=None)\n\n'
docstr += 'Stochastic variable with '+name+' distribution.\nParents are: '+', '.join(parent_names) + '.\n\n'
docstr += """
Methods:
random()
- draws random value
sets value to return value
ppf(q)
- percent point function (inverse of cdf --- percentiles)
sets value to return value
isf(q)
- inverse survival function (inverse of sf)
sets value to return value
stats(moments='mv')
- mean('m',axis=0), variance('v'), skew('s'), and/or kurtosis('k')
Attributes:
logp
- sum(log(pdf())) or sum(log(pmf()))
cdf
- cumulative distribution function
sf
- survival function (1-cdf --- sometimes more accurate)
entropy
- (differential) entropy of the RV.
NOTE: If you encounter difficulties with this object, please try the analogous
computation using the rv objects in scipy.stats.distributions directly before
reporting the bug.
"""
new_class = new_dist_class(dtype, name, parent_names, parents_default, docstr, logp, random, True, None)
class newer_class(new_class):
__doc__ = docstr
rv = scipy_dist
def __init__(self, *args, **kwds):
new_class.__init__(self, *args, **kwds)
self.args, self.kwds = separate_shape_args(self.parents, shape_args)
self.frozen_rv = self.rv(self.args, self.kwds)
def _pymc_dists_to_value(self, args):
"""Replace arguments that are a pymc.Node with their value."""
# This is needed because the scipy rv function transforms
# every input argument which causes new pymc lambda
# functions to be generated. Thus, when calling this many
# many times, excessive amounts of RAM are used.
new_args = []
for arg in args:
if isinstance(arg, pm.Node):
new_args.append(arg.value)
else:
new_args.append(arg)
return new_args
def _cdf(self):
"""
The cumulative distribution function of self conditional on parents
evaluated at self's current value
"""
return self.rv.cdf(self.value, *self._pymc_dists_to_value(self.args), **self.kwds)
cdf = property(_cdf, doc=_cdf.__doc__)
def _sf(self):
"""
The survival function of self conditional on parents
evaluated at self's current value
"""
return self.rv.sf(self.value, *self._pymc_dists_to_value(self.args), **self.kwds)
sf = property(_sf, doc=_sf.__doc__)
def ppf(self, q):
"""
The percentile point function (inverse cdf) of self conditional on parents.
Self's value will be set to the return value.
"""
self.value = self.rv.ppf(q, *self._pymc_dists_to_value(self.args), **self.kwds)
return self.value
def isf(self, q):
"""
The inverse survival function of self conditional on parents.
Self's value will be set to the return value.
"""
self.value = self.rv.isf(q, *self._pymc_dists_to_value(self.args), **self.kwds)
return self.value
def stats(self, moments='mv'):
"""The first few moments of self's distribution conditional on parents"""
return self.rv.stats(moments=moments, *self._pymc_dists_to_value(self.args), **self.kwds)
def _entropy(self):
"""The entropy of self's distribution conditional on its parents"""
return self.rv.entropy(*self._pymc_dists_to_value(self.args), **self.kwds)
entropy = property(_entropy, doc=_entropy.__doc__)
newer_class.__name__ = new_class.__name__
return newer_class
for scipy_dist_name in sc_dst.__all__:
scipy_dist = sc_dst.__dict__[scipy_dist_name]
if isinstance(scipy_dist, sc_dst.rv_continuous) or isinstance(scipy_dist, sc_dst.rv_discrete):
new_dist = stochastic_from_scipy_dist(scipy_dist)
locals()[new_dist.__name__] = new_dist
__all__.append(new_dist.__name__)
|