[Scipy-svn] r6963 - in trunk/scipy/stats: . tests
scipy-svn at scipy.org
scipy-svn at scipy.org
Sun Nov 28 08:37:29 EST 2010
Author: rgommers
Date: 2010-11-28 07:37:29 -0600 (Sun, 28 Nov 2010)
New Revision: 6963
Modified:
trunk/scipy/stats/distributions.py
trunk/scipy/stats/info.py
trunk/scipy/stats/tests/test_continuous_basic.py
Log:
DEP: deprecate stats.invnorm, add a copy as invgauss. See #1172.
Modified: trunk/scipy/stats/distributions.py
===================================================================
--- trunk/scipy/stats/distributions.py 2010-11-28 13:37:08 UTC (rev 6962)
+++ trunk/scipy/stats/distributions.py 2010-11-28 13:37:29 UTC (rev 6963)
@@ -6,6 +6,7 @@
#
import math
+import warnings
from copy import copy
from scipy.misc import comb, derivative
@@ -57,10 +58,10 @@
'weibull_max', 'genlogistic', 'genpareto', 'genexpon', 'genextreme',
'gamma', 'gengamma', 'genhalflogistic', 'gompertz', 'gumbel_r',
'gumbel_l', 'halfcauchy', 'halflogistic', 'halfnorm', 'hypsecant',
- 'gausshyper', 'invgamma', 'invnorm', 'invweibull', 'johnsonsb',
- 'johnsonsu', 'laplace', 'levy', 'levy_l', 'levy_stable',
- 'logistic', 'loggamma', 'loglaplace', 'lognorm', 'gilbrat',
- 'maxwell', 'mielke', 'nakagami', 'ncx2', 'ncf', 't',
+ 'gausshyper', 'invgamma', 'invnorm', 'invgauss', 'invweibull',
+ 'johnsonsb', 'johnsonsu', 'laplace', 'levy', 'levy_l',
+ 'levy_stable', 'logistic', 'loggamma', 'loglaplace', 'lognorm',
+ 'gilbrat', 'maxwell', 'mielke', 'nakagami', 'ncx2', 'ncf', 't',
'nct', 'pareto', 'lomax', 'powerlaw', 'powerlognorm', 'powernorm',
'rdist', 'rayleigh', 'reciprocal', 'rice', 'recipinvgauss',
'semicircular', 'triang', 'truncexpon', 'truncnorm',
@@ -3247,30 +3248,67 @@
## Inverse Normal Distribution
# scale is gamma from DATAPLOT and B from Regress
+_invnorm_msg = \
+"""The `invnorm` distribution will be renamed to `invgauss` after scipy 0.9"""
class invnorm_gen(rv_continuous):
def _rvs(self, mu):
+ warnings.warn(_invnorm_msg, DeprecationWarning)
return mtrand.wald(mu, 1.0, size=self._size)
def _pdf(self, x, mu):
+ warnings.warn(_invnorm_msg, DeprecationWarning)
return 1.0/sqrt(2*pi*x**3.0)*exp(-1.0/(2*x)*((x-mu)/mu)**2)
def _logpdf(self, x, mu):
+ warnings.warn(_invnorm_msg, DeprecationWarning)
return -0.5*log(2*pi) - 1.5*log(x) - ((x-mu)/mu)**2/(2*x)
def _cdf(self, x, mu):
+ warnings.warn(_invnorm_msg, DeprecationWarning)
fac = sqrt(1.0/x)
C1 = norm.cdf(fac*(x-mu)/mu)
C1 += exp(2.0/mu)*norm.cdf(-fac*(x+mu)/mu)
return C1
def _stats(self, mu):
+ warnings.warn(_invnorm_msg, DeprecationWarning)
return mu, mu**3.0, 3*sqrt(mu), 15*mu
invnorm = invnorm_gen(a=0.0, name='invnorm', longname="An inverse normal",
shapes="mu",extradoc="""
Inverse normal distribution
+NOTE: `invnorm` will be renamed to `invgauss` after scipy 0.9
+
invnorm.pdf(x,mu) = 1/sqrt(2*pi*x**3) * exp(-(x-mu)**2/(2*x*mu**2))
for x > 0.
"""
)
+## Inverse Gaussian Distribution (used to be called 'invnorm'
+# scale is gamma from DATAPLOT and B from Regress
+
+class invgauss_gen(rv_continuous):
+ def _rvs(self, mu):
+ return mtrand.wald(mu, 1.0, size=self._size)
+ def _pdf(self, x, mu):
+ return 1.0/sqrt(2*pi*x**3.0)*exp(-1.0/(2*x)*((x-mu)/mu)**2)
+ def _logpdf(self, x, mu):
+ return -0.5*log(2*pi) - 1.5*log(x) - ((x-mu)/mu)**2/(2*x)
+ def _cdf(self, x, mu):
+ fac = sqrt(1.0/x)
+ C1 = norm.cdf(fac*(x-mu)/mu)
+ C1 += exp(2.0/mu)*norm.cdf(-fac*(x+mu)/mu)
+ return C1
+ def _stats(self, mu):
+ return mu, mu**3.0, 3*sqrt(mu), 15*mu
+invgauss = invgauss_gen(a=0.0, name='invgauss', longname="An inverse Gaussian",
+ shapes="mu",extradoc="""
+
+Inverse Gaussian distribution
+
+invgauss.pdf(x,mu) = 1/sqrt(2*pi*x**3) * exp(-(x-mu)**2/(2*x*mu**2))
+for x > 0.
+"""
+ )
+
+
## Inverted Weibull
class invweibull_gen(rv_continuous):
@@ -4136,7 +4174,7 @@
# FIXME: PPF does not work.
class recipinvgauss_gen(rv_continuous):
- def _rvs(self, mu): #added, taken from invnorm
+ def _rvs(self, mu): #added, taken from invgauss
return 1.0/mtrand.wald(mu, 1.0, size=self._size)
def _pdf(self, x, mu):
return 1.0/sqrt(2*pi*x)*exp(-(1-mu*x)**2.0 / (2*x*mu**2.0))
@@ -4403,7 +4441,7 @@
## Wald distribution (Inverse Normal with shape parameter mu=1.0)
-class wald_gen(invnorm_gen):
+class wald_gen(invgauss_gen):
"""A Wald continuous random variable.
%(before_notes)s
@@ -4418,11 +4456,11 @@
def _rvs(self):
return mtrand.wald(1.0, 1.0, size=self._size)
def _pdf(self, x):
- return invnorm._pdf(x, 1.0)
+ return invgauss._pdf(x, 1.0)
def _logpdf(self, x):
- return invnorm._logpdf(x, 1.0)
+ return invgauss._logpdf(x, 1.0)
def _cdf(self, x):
- return invnorm._cdf(x, 1.0)
+ return invgauss._cdf(x, 1.0)
def _stats(self):
return 1.0, 1.0, 3.0, 15.0
wald = wald_gen(a=0.0, name="wald", extradoc="""
Modified: trunk/scipy/stats/info.py
===================================================================
--- trunk/scipy/stats/info.py 2010-11-28 13:37:08 UTC (rev 6962)
+++ trunk/scipy/stats/info.py 2010-11-28 13:37:29 UTC (rev 6963)
@@ -81,6 +81,7 @@
hypsecant Hyperbolic Secant
invgamma Inverse Gamma
invnorm Inverse Normal
+ invgauss Inverse Gaussian
invweibull Inverse Weibull
johnsonsb Johnson SB
johnsonsu Johnson SU
Modified: trunk/scipy/stats/tests/test_continuous_basic.py
===================================================================
--- trunk/scipy/stats/tests/test_continuous_basic.py 2010-11-28 13:37:08 UTC (rev 6962)
+++ trunk/scipy/stats/tests/test_continuous_basic.py 2010-11-28 13:37:29 UTC (rev 6963)
@@ -1,3 +1,5 @@
+import warnings
+
import numpy.testing as npt
import numpy as np
import nose
@@ -15,7 +17,7 @@
not for numerically exact results.
-TODO:
+TODO:
* make functioning test for skew and kurtosis
still known failures - skip for now
@@ -70,6 +72,7 @@
['hypsecant', ()],
['invgamma', (2.0668996136993067,)],
['invnorm', (0.14546264555347513,)],
+ ['invgauss', (0.14546264555347513,)],
['invweibull', (10.58,)], # sample mean test fails at(0.58847112119264788,)]
['johnsonsb', (4.3172675099141058, 3.1837781130785063)],
['johnsonsu', (2.554395574161155, 2.2482281679651965)],
@@ -144,7 +147,7 @@
'loglaplace', 'rdist', 'semicircular', 'invweibull', 'ksone',
'cosine', 'kstwobign', 'truncnorm', 'mielke', 'recipinvgauss', 'levy',
'johnsonsu', 'levy_l', 'powernorm', 'wrapcauchy',
- 'johnsonsb', 'truncexpon', 'rice', 'invnorm', 'invgamma',
+ 'johnsonsb', 'truncexpon', 'rice', 'invnorm', 'invgauss', 'invgamma',
'powerlognorm']
distmiss = [[dist,args] for dist,args in distcont if dist in distmissing]
@@ -166,7 +169,7 @@
skurt = stats.kurtosis(rvs)
sskew = stats.skew(rvs)
m,v = distfn.stats(*arg)
-
+
yield check_sample_meanvar_, distfn, arg, m, v, sm, sv, sn, distname + \
'sample mean test'
# the sample skew kurtosis test has known failures, not very good distance measure
@@ -219,14 +222,14 @@
if not np.isinf(m):
npt.assert_almost_equal(m1, m, decimal=10, err_msg= msg + \
' - 1st moment')
- else: # or np.isnan(m1),
- npt.assert_(np.isinf(m1),
+ else: # or np.isnan(m1),
+ npt.assert_(np.isinf(m1),
msg + ' - 1st moment -infinite, m1=%s' % str(m1))
#np.isnan(m1) temporary special treatment for loggamma
if not np.isinf(v):
npt.assert_almost_equal(m2-m1*m1, v, decimal=10, err_msg= msg + \
' - 2ndt moment')
- else: #or np.isnan(m2),
+ else: #or np.isnan(m2),
npt.assert_(np.isinf(m2),
msg + ' - 2nd moment -infinite, m2=%s' % str(m2))
#np.isnan(m2) temporary special treatment for loggamma
@@ -356,6 +359,10 @@
npt.assert_(pval > alpha, "D = " + str(D) + "; pval = " + str(pval) +
"; alpha = " + str(alpha) + "\nargs = " + str(args))
+
+warnings.filterwarnings('ignore', message="The `invnorm` distribution")
+
+
if __name__ == "__main__":
#nose.run(argv=['', __file__])
nose.runmodule(argv=[__file__,'-s'], exit=False)
More information about the Scipy-svn
mailing list