Added ppimport.py

Small bugfix to estimation.py and distributions.py
master
Per.Andreas.Brodtkorb 15 years ago
parent 652bf4981d
commit da96ab184a

@ -1,5 +1,5 @@
#-------------------------------------------------------------------------------
# Name: module1
# Purpose:
#
@ -8,7 +8,7 @@
# Created: 16.09.2008
# Copyright: (c) pab 2008
# Licence: <your licence>
#-------------------------------------------------------------------------------
#!/usr/bin/env python

@ -0,0 +1,441 @@
#!/usr/bin/env python
"""
Postpone module import to future.
Python versions: 1.5.2 - 2.3.x
Author: Pearu Peterson <pearu@cens.ioc.ee>
Created: March 2003
$Revision: 922 $
$Date: 2004-11-27 14:23:27 -0700 (Sat, 27 Nov 2004) $
"""
__all__ = ['ppimport','ppimport_attr','ppresolve']
import os
import sys
import types
import traceback
DEBUG=0
_ppimport_is_enabled = 1
def enable():
""" Enable postponed importing."""
global _ppimport_is_enabled
_ppimport_is_enabled = 1
def disable():
""" Disable postponed importing."""
global _ppimport_is_enabled
_ppimport_is_enabled = 0
class PPImportError(ImportError):
pass
def _get_so_ext(_cache={}):
so_ext = _cache.get('so_ext')
if so_ext is None:
if sys.platform[:5]=='linux':
so_ext = '.so'
else:
try:
# if possible, avoid expensive get_config_vars call
from distutils.sysconfig import get_config_vars
so_ext = get_config_vars('SO')[0] or ''
except ImportError:
#XXX: implement hooks for .sl, .dll to fully support
# Python 1.5.x
so_ext = '.so'
_cache['so_ext'] = so_ext
return so_ext
def _get_frame(level=0):
try:
return sys._getframe(level+1)
except AttributeError:
# Python<=2.0 support
frame = sys.exc_info()[2].tb_frame
for i in range(level+1):
frame = frame.f_back
return frame
def ppimport_attr(module, name):
""" ppimport(module, name) is 'postponed' getattr(module, name)
"""
global _ppimport_is_enabled
if _ppimport_is_enabled and isinstance(module, _ModuleLoader):
return _AttrLoader(module, name)
return getattr(module, name)
class _AttrLoader(object):
def __init__(self, module, name):
self.__dict__['_ppimport_attr_module'] = module
self.__dict__['_ppimport_attr_name'] = name
def _ppimport_attr_getter(self):
module = self.__dict__['_ppimport_attr_module']
if isinstance(module, _ModuleLoader):
# in case pp module was loaded by other means
module = sys.modules[module.__name__]
attr = getattr(module,
self.__dict__['_ppimport_attr_name'])
try:
d = attr.__dict__
if d is not None:
self.__dict__ = d
except AttributeError:
pass
self.__dict__['_ppimport_attr'] = attr
return attr
def __nonzero__(self):
return 1
def __getattr__(self, name):
try:
attr = self.__dict__['_ppimport_attr']
except KeyError:
attr = self._ppimport_attr_getter()
if name=='_ppimport_attr':
return attr
return getattr(attr, name)
def __repr__(self):
if '_ppimport_attr' in self.__dict__:
return repr(self._ppimport_attr)
module = self.__dict__['_ppimport_attr_module']
name = self.__dict__['_ppimport_attr_name']
return "<attribute %s of %s>" % (`name`,`module`)
__str__ = __repr__
# For function and class attributes.
def __call__(self, *args, **kwds):
return self._ppimport_attr(*args,**kwds)
def _is_local_module(p_dir,name,suffices):
base = os.path.join(p_dir,name)
for suffix in suffices:
if os.path.isfile(base+suffix):
if p_dir:
return base+suffix
return name+suffix
def ppimport(name):
""" ppimport(name) -> module or module wrapper
If name has been imported before, return module. Otherwise
return ModuleLoader instance that transparently postpones
module import until the first attempt to access module name
attributes.
"""
global _ppimport_is_enabled
level = 1
parent_frame = p_frame = _get_frame(level)
while '__name__' not in p_frame.f_locals:
level = level + 1
p_frame = _get_frame(level)
p_name = p_frame.f_locals['__name__']
if p_name=='__main__':
p_dir = ''
fullname = name
elif '__path__' in p_frame.f_locals:
# python package
p_path = p_frame.f_locals['__path__']
p_dir = p_path[0]
fullname = p_name + '.' + name
else:
# python module
p_file = p_frame.f_locals['__file__']
p_dir = os.path.dirname(p_file)
fullname = p_name + '.' + name
# module may be imported already
module = sys.modules.get(fullname)
if module is not None:
if _ppimport_is_enabled or isinstance(module, types.ModuleType):
return module
return module._ppimport_importer()
so_ext = _get_so_ext()
py_exts = ('.py','.pyc','.pyo')
so_exts = (so_ext,'module'+so_ext)
for d,n,fn,e in [\
# name is local python module or local extension module
(p_dir, name, fullname, py_exts+so_exts),
# name is local package
(os.path.join(p_dir, name), '__init__', fullname, py_exts),
# name is package in parent directory (scipy specific)
(os.path.join(os.path.dirname(p_dir), name), '__init__', name, py_exts),
]:
location = _is_local_module(d, n, e)
if location is not None:
fullname = fn
break
if location is None:
# name is to be looked in python sys.path.
fullname = name
location = 'sys.path'
# Try once more if module is imported.
# This covers the case when importing from python module
module = sys.modules.get(fullname)
if module is not None:
if _ppimport_is_enabled or isinstance(module,types.ModuleType):
return module
return module._ppimport_importer()
# It is OK if name does not exists. The ImportError is
# postponed until trying to use the module.
loader = _ModuleLoader(fullname,location,p_frame=parent_frame)
if _ppimport_is_enabled:
return loader
return loader._ppimport_importer()
def _get_frame_code(frame):
filename = frame.f_code.co_filename
lineno = frame.f_lineno
result = '%s in %s:\n' % (filename,frame.f_code.co_name)
if not os.path.isfile(filename):
return result
f = open(filename)
i = 1
line = f.readline()
while line:
line = f.readline()
i = i + 1
if (abs(i-lineno)<2):
result += '#%d: %s\n' % (i,line.rstrip())
if i>lineno+3:
break
f.close()
return result
def frame_traceback(frame):
if not frame:
return
blocks = []
f = frame
while f:
blocks.insert(0,_get_frame_code(f))
f = f.f_back
print '='*50
print '\n'.join(blocks)
print '='*50
class _ModuleLoader(object):
# Don't use it directly. Use ppimport instead.
def __init__(self,name,location,p_frame=None):
# set attributes, avoid calling __setattr__
self.__dict__['__name__'] = name
self.__dict__['__file__'] = location
self.__dict__['_ppimport_p_frame'] = p_frame
if location != 'sys.path':
from numpy.testing import Tester
self.__dict__['test'] = Tester(os.path.dirname(location)).test
# install loader
sys.modules[name] = self
def _ppimport_importer(self):
name = self.__name__
try:
module = sys.modules[name]
except KeyError:
raise ImportError,self.__dict__.get('_ppimport_exc_info')[1]
if module is not self:
exc_info = self.__dict__.get('_ppimport_exc_info')
if exc_info is not None:
raise PPImportError,\
''.join(traceback.format_exception(*exc_info))
else:
assert module is self,`(module, self)`
# uninstall loader
del sys.modules[name]
if DEBUG:
print 'Executing postponed import for %s' %(name)
try:
module = __import__(name,None,None,['*'])
except Exception,msg: # ImportError:
if DEBUG:
p_frame = self.__dict__.get('_ppimport_p_frame',None)
frame_traceback(p_frame)
self.__dict__['_ppimport_exc_info'] = sys.exc_info()
raise
assert isinstance(module,types.ModuleType),`module`
#self.__dict__ = module.__dict__
self.__dict__.update(module.__dict__)
self.__dict__['_ppimport_module'] = module
# XXX: Should we check the existence of module.test? Warn?
from numpy.testing import Tester
test = Tester(os.path.dirname(module.__file__)).test
return module
def __setattr__(self, name, value):
try:
module = self.__dict__['_ppimport_module']
except KeyError:
module = self._ppimport_importer()
return setattr(module, name, value)
def __getattr__(self, name):
try:
module = self.__dict__['_ppimport_module']
except KeyError:
module = self._ppimport_importer()
return getattr(module, name)
def __repr__(self):
global _ppimport_is_enabled
if not _ppimport_is_enabled:
try:
module = self.__dict__['_ppimport_module']
except KeyError:
module = self._ppimport_importer()
return module.__repr__()
if '_ppimport_module' in self.__dict__:
status = 'imported'
elif '_ppimport_exc_info' in self.__dict__:
status = 'import error'
else:
status = 'import postponed'
return '<module %s from %s [%s]>' \
% (`self.__name__`,`self.__file__`, status)
__str__ = __repr__
def ppresolve(a,ignore_failure=None):
""" Return resolved object a.
a can be module name, postponed module, postponed modules
attribute, string representing module attribute, or any
Python object.
"""
global _ppimport_is_enabled
if _ppimport_is_enabled:
disable()
a = ppresolve(a,ignore_failure=ignore_failure)
enable()
return a
if type(a) is type(''):
ns = a.split('.')
if ignore_failure:
try:
a = ppimport(ns[0])
except:
return a
else:
a = ppimport(ns[0])
b = [ns[0]]
del ns[0]
while ns:
if hasattr(a,'_ppimport_importer') or \
hasattr(a,'_ppimport_module'):
a = getattr(a,'_ppimport_module',a)
if hasattr(a,'_ppimport_attr'):
a = a._ppimport_attr
b.append(ns[0])
del ns[0]
if ignore_failure and not hasattr(a, b[-1]):
a = '.'.join(ns+b)
b = '.'.join(b)
if b in sys.modules and sys.modules[b] is None:
del sys.modules[b]
return a
a = getattr(a,b[-1])
if hasattr(a,'_ppimport_importer') or \
hasattr(a,'_ppimport_module'):
a = getattr(a,'_ppimport_module',a)
if hasattr(a,'_ppimport_attr'):
a = a._ppimport_attr
return a
def _ppresolve_ignore_failure(a):
return ppresolve(a,ignore_failure=1)
try:
import pydoc as _pydoc
except ImportError:
_pydoc = None
if _pydoc is not None:
# Redefine __call__ method of help.__class__ to
# support ppimport.
import new as _new
_old_pydoc_help_call = _pydoc.help.__class__.__call__
def _ppimport_pydoc_help_call(self,*args,**kwds):
return _old_pydoc_help_call(self, *map(_ppresolve_ignore_failure,args),
**kwds)
_ppimport_pydoc_help_call.__doc__ = _old_pydoc_help_call.__doc__
_pydoc.help.__class__.__call__ = _new.instancemethod(_ppimport_pydoc_help_call,
None,
_pydoc.help.__class__)
_old_pydoc_Doc_document = _pydoc.Doc.document
def _ppimport_pydoc_Doc_document(self,*args,**kwds):
args = (_ppresolve_ignore_failure(args[0]),) + args[1:]
return _old_pydoc_Doc_document(self,*args,**kwds)
_ppimport_pydoc_Doc_document.__doc__ = _old_pydoc_Doc_document.__doc__
_pydoc.Doc.document = _new.instancemethod(_ppimport_pydoc_Doc_document,
None,
_pydoc.Doc)
_old_pydoc_describe = _pydoc.describe
def _ppimport_pydoc_describe(object):
return _old_pydoc_describe(_ppresolve_ignore_failure(object))
_ppimport_pydoc_describe.__doc__ = _old_pydoc_describe.__doc__
_pydoc.describe = _ppimport_pydoc_describe
import inspect as _inspect
_old_inspect_getfile = _inspect.getfile
def _ppimport_inspect_getfile(object):
if isinstance(object,_ModuleLoader):
return object.__dict__['__file__']
return _old_inspect_getfile(_ppresolve_ignore_failure(object))
_ppimport_inspect_getfile.__doc__ = _old_inspect_getfile.__doc__
_inspect.getfile = _ppimport_inspect_getfile
_old_inspect_getdoc = _inspect.getdoc
def _ppimport_inspect_getdoc(object):
return _old_inspect_getdoc(_ppresolve_ignore_failure(object))
_ppimport_inspect_getdoc.__doc__ = _old_inspect_getdoc.__doc__
_inspect.getdoc = _ppimport_inspect_getdoc
_old_inspect_getsource = _inspect.getsource
def _ppimport_inspect_getsource(object):
return _old_inspect_getsource(_ppresolve_ignore_failure(object))
_ppimport_inspect_getsource.__doc__ = _old_inspect_getsource.__doc__
_inspect.getsource = _ppimport_inspect_getsource
import __builtin__ as _builtin
_old_builtin_dir = _builtin.dir
def _ppimport_builtin_dir(*arg):
if not arg:
p_frame = _get_frame(1)
g = p_frame.f_globals
l = p_frame.f_locals
l['_ppimport_old_builtin_dir'] = _old_builtin_dir
r = eval('_ppimport_old_builtin_dir()',g,l)
del r[r.index('_ppimport_old_builtin_dir')]
return r
return _old_builtin_dir(*map(_ppresolve_ignore_failure,arg))
_ppimport_builtin_dir.__doc__ = _old_builtin_dir.__doc__
_builtin.dir = _ppimport_builtin_dir

@ -6,7 +6,7 @@ Statistics package in WAFO Toolbox.
"""
from scipy.stats import *
from wafo.stats.core import *
#from wafo.stats.distributions import *
from wafo.stats.distributions import *
#from wafo.spectrum.core import SpecData1D
#import wafo.spectrum.models
#import wafo.spectrum.dispersion_relation

@ -163,7 +163,7 @@ class rv_frozen(object):
def __init__(self, dist, *args, **kwds):
self.dist = dist
loc0, scale0 = map(kwds.get, ['loc', 'scale'])
if isinstance(dist, rv_continuous):
if hasattr(dist, 'fix_loc_scale'): #isinstance(dist, rv_continuous):
args, loc0, scale0 = dist.fix_loc_scale(args, loc0, scale0)
self.par = args + (loc0, scale0)
else: # rv_discrete
@ -569,7 +569,7 @@ class rv_generic(object):
and rv_continuous.
"""
def fix_loc_scale(self, args, loc, scale=1):
def _fix_loc_scale(self, args, loc, scale=1):
N = len(args)
if N > self.numargs:
if N == self.numargs + 1 and loc is None:
@ -585,8 +585,8 @@ class rv_generic(object):
loc = 0.0
return args, loc, scale
def fix_loc(self, args, loc):
args, loc, scale = self.fix_loc_scale(args, loc)
def _fix_loc(self, args, loc):
args, loc, scale = self._fix_loc_scale(args, loc)
return args, loc
# These are actually called, and should not be overwritten if you
@ -617,7 +617,7 @@ class rv_generic(object):
loc, scale, size, discrete = map(kwds.get, kwd_names,
[None] * len(kwd_names))
args, loc, scale = self.fix_loc_scale(args, loc, scale)
args, loc, scale = self._fix_loc_scale(args, loc, scale)
cond = logical_and(self._argcheck(*args), (scale >= 0))
if not all(cond):
raise ValueError, "Domain error in arguments."
@ -765,6 +765,7 @@ class rv_continuous(rv_generic):
shapes=None, extradoc=None):
rv_generic.__init__(self)
self.fix_loc_scale = self._fix_loc_scale
if badvalue is None:
badvalue = nan
@ -4518,7 +4519,7 @@ class rv_discrete(rv_generic):
shapes=None, extradoc=None):
rv_generic.__init__(self)
self.fix_loc = self._fix_loc
if badvalue is None:
badvalue = nan
self.badvalue = badvalue
@ -4543,9 +4544,9 @@ class rv_discrete(rv_generic):
self.pk = take(ravel(self.pk), indx, 0)
self.a = self.xk[0]
self.b = self.xk[-1]
self.P = make_dict(self.xk, self.pk)
self.P = dict(zip(self.xk, self.pk)) #make_dict(self.xk, self.pk)
self.qvals = numpy.cumsum(self.pk, axis=0)
self.F = make_dict(self.xk, self.qvals)
self.F = dict(zip(self.xk, self.qvals)) #make_dict(self.xk, self.qvals)
self.Finv = reverse_dict(self.F)
self._ppf = new.instancemethod(sgf(_drv_ppf, otypes='d'),
self, rv_discrete)

@ -1,26 +1,23 @@
# Functions to implement several important functions for
# various Continous and Discrete Probability Distributions
#
# Author: Travis Oliphant 2002-2003
'''
Contains FitDistribution and Profile class, which are
important classes for fitting to various Continous and Discrete Probability Distributions
Author: Per A. Brodtkorb 2008
'''
from __future__ import division
from wafo.plotbackend import plotbackend
from wafo.misc import ecross, findcross, JITImport
from wafo.misc import ecross, findcross
from scipy.misc.ppimport import ppimport
import numdifftools
from scipy import special
from scipy.linalg import pinv2
from scipy import optimize
# Trick to avoid error due to circular import
#_WAFOCOV = ppimport('wafo.covariance')
_WAFODIST = ppimport('wafo.stats.distributions')
#_WAFODIST = JITImport('wafo.stats.distributions')
from numpy import alltrue, arange, \
ravel, ones, sum, \
zeros, log, sqrt, exp
@ -43,6 +40,14 @@ arr = asarray
all = alltrue
def chi2isf(p,df):
return special.chdtri(df, p)
def chi2sf(x, df):
return special.chdtrc(df, x)
def norm_ppf(q):
return special.ndtri(q)
def valarray(shape, value=nan, typecode=None):
"""Return an array of all value.
"""
@ -87,7 +92,7 @@ class rv_frozen(object):
def __init__(self, dist, *args, **kwds):
self.dist = dist
loc0, scale0 = map(kwds.get, ['loc', 'scale'])
if isinstance(dist, _WAFODIST.rv_continuous):
if hasattr(dist, 'fix_loc_scale'): #isinstance(dist, _WAFODIST.rv_continuous):
args, loc0, scale0 = dist.fix_loc_scale(args, loc0, scale0)
self.par = args + (loc0, scale0)
else: # rv_discrete
@ -247,7 +252,7 @@ class Profile(object):
self.i_free = nonzero(isfree)
self.Lmax = Lmax
self.alpha_Lrange = 0.5 * _WAFODIST.chi2.isf(self.alpha, 1)
self.alpha_Lrange = 0.5 * chi2isf(self.alpha, 1) #_WAFODIST.chi2.isf(self.alpha, 1)
self.alpha_cross_level = Lmax - self.alpha_Lrange
lowLevel = self.alpha_cross_level - self.alpha_Lrange / 7.0
@ -344,7 +349,8 @@ class Profile(object):
pcov = self.fit_dist.par_cov[i_notfixed, :][:, i_notfixed]
pvar = sum(numpy.dot(drl, pcov) * drl)
p_crit = _WAFODIST.norm.isf(self.alpha / 2.0) * sqrt(numpy.ravel(pvar)) * 1.5
#p_crit = _WAFODIST.norm.isf(self.alpha / 2.0) * sqrt(numpy.ravel(pvar)) * 1.5
p_crit = -norm_ppf(self.alpha / 2.0) * sqrt(numpy.ravel(pvar)) * 1.5
if self.pmin == None:
self.pmin = p_opt - 5.0 * p_crit
if self.pmax == None:
@ -392,7 +398,7 @@ class Profile(object):
if alpha < self.alpha:
raise ValueError('Unable to return CI with alpha less than %g' % self.alpha)
cross_level = self.Lmax - 0.5 * _WAFODIST.chi2.isf(alpha, 1)
cross_level = self.Lmax - 0.5 * chi2isf(alpha, 1) #_WAFODIST.chi2.isf(alpha, 1)
ind = findcross(self.data, cross_level)
N = len(ind)
if N == 0:
@ -536,7 +542,7 @@ class FitDistribution(rv_frozen):
self.par_cov[:, :] = nan
pvar = numpy.diag(self.par_cov)
zcrit = -_WAFODIST.norm.ppf(self.alpha / 2.0)
zcrit = -norm_ppf(self.alpha / 2.0)#_WAFODIST.norm.ppf(self.alpha / 2.0)
self.par_lower = self.par - zcrit * sqrt(pvar)
self.par_upper = self.par + zcrit * sqrt(pvar)
@ -702,7 +708,7 @@ class FitDistribution(rv_frozen):
Other distribution types will introduce deviations in the plot.
'''
bin, limits = numpy.histogram(self.data, normed=True, new=True)
bin, limits = numpy.histogram(self.data, normed=True) #, new=True)
limits.shape = (-1, 1)
xx = limits.repeat(3, axis=1)
xx.shape = (-1,)
@ -793,7 +799,7 @@ class FitDistribution(rv_frozen):
C1 = m - sqrt(0.5 * n * v)
C2 = sqrt(v / (2.0 * n))
Tn = (T + 0.5 * k * isParUnKnown - C1) / C2 # chi2 with n degrees of freedom
pvalue = _WAFODIST.chi2.sf(Tn, n)
pvalue = chi2sf(Tn, n) #_WAFODIST.chi2.sf(Tn, n)
return pvalue
def nlogps(self, theta, x):
@ -964,6 +970,7 @@ class FitDistribution(rv_frozen):
return - H
def main():
_WAFODIST = ppimport('wafo.stats.distributions')
#nbinom(10, 0.75).rvs(3)
import matplotlib
matplotlib.interactive(True)
@ -986,15 +993,15 @@ def main():
Lx.get_CI(alpha=0.2)
# CI for logSF=log(SF)
Lpr = phat.profile(i=1, logSF=log(R), link=phat.dist.link)
Lpr = phat.profile(i=0, logSF=log(R), link=phat.dist.link)
Lpr.plot()
Lpr.get_CI(alpha=0.075)
dlaplace.stats(0.8, loc=0)
_WAFODIST.dlaplace.stats(0.8, loc=0)
# pass
t = planck(0.51000000000000001)
t = _WAFODIST.planck(0.51000000000000001)
t.ppf(0.5)
t = zipf(2)
t = _WAFODIST.zipf(2)
t.ppf(0.5)
import pylab as plb
_WAFODIST.rice.rvs(1)

Loading…
Cancel
Save