Source Code for Module mvpa.support.stats

 1  # emacs: -*- mode: python; py-indent-offset: 4; indent-tabs-mode: nil -*- 
 2  # vi: set ft=python sts=4 ts=4 sw=4 et: 
 3  ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ## 
 4  # 
 5  #   See COPYING file distributed along with the PyMVPA package for the 
 6  #   copyright and license terms. 
 7  # 
 8  ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ## 
 9  """Fixer for rdist in scipy 
10  """ 
11   
12  __docformat__ = 'restructuredtext' 
13   
14  from mvpa.base import externals 
15   
16  if __debug__: 
17      from mvpa.base import debug 
18   
19  if externals.exists('scipy', raiseException=True): 
20      import scipy 
21      import scipy.stats 
22      import scipy.stats as stats 
23   
24  if not externals.exists('good scipy.stats.rdist'): 
25      if __debug__: 
26          debug("EXT", "Fixing up scipy.stats.rdist") 
27      # Lets fix it up, future imports of scipy.stats should carry fixed 
28      # version, isn't python is \emph{evil} ;-) 
29      import numpy as N 
30   
31      from scipy.stats.distributions import rv_continuous 
32      from scipy import special 
33      import scipy.integrate 
34   
35      # NB: Following function is copied from scipy SVN rev.5236 
36      #     and fixed with pow -> N.power (thanks Josef!) 
37      # FIXME: PPF does not work. 
38 -    class rdist_gen(rv_continuous): 
39 -        def _pdf(self, x, c): 
40              return N.power((1.0-x*x),c/2.0-1) / special.beta(0.5,c/2.0) 
41 -        def _cdf_skip(self, x, c): 
42              #error inspecial.hyp2f1 for some values see tickets 758, 759 
43              return 0.5 + x/special.beta(0.5,c/2.0)* \ 
44                     special.hyp2f1(0.5,1.0-c/2.0,1.5,x*x) 
45 -        def _munp(self, n, c): 
46              return (1-(n % 2))*special.beta((n+1.0)/2,c/2.0) 
47   
48      # Lets try to avoid at least some of the numerical problems by removing points 
49      # around edges 
50      rdist = rdist_gen(a=-1.0, b=1.0, name="rdist", longname="An R-distributed", 
51                        shapes="c", extradoc=""" 
52   
53      R-distribution 
54   
55      rdist.pdf(x,c) = (1-x**2)**(c/2-1) / B(1/2, c/2) 
56      for -1 <= x <= 1, c > 0. 
57      """ 
58                        ) 
59      # Fix up number of arguments for veccdf's vectorize 
60      if rdist.veccdf.nin == 1: 
61          if __debug__: 
62              debug("EXT", "Fixing up veccdf.nin to make 2 for rdist") 
63          rdist.veccdf.nin = 2 
64   
65      scipy.stats.distributions.rdist_gen = scipy.stats.rdist_gen = rdist_gen 
66      scipy.stats.distributions.rdist = scipy.stats.rdist = rdist 
67