Mailman 3 November 2008 - SciPy-Dev

symeig integration ready
by Tiziano Zito Nov. 18, 2008

Nov. 18, 2008

Dear devs, I have managed to integrate symeig in scipy, but before doing a commit I would like to get your feedback: - scipy.linalg.eigh gets the following signature: eigh(a, b=None, lower=True, eigvals_only=False, overwrite_a=False, overwrite_b=False, turbo=True, eigvals=None, type=1): this is different from the current signature: eigh(a, lower=True, eigvals_only=False, overwrite_a=False): Is it OK to break backward-compatibility? Code would only break if someone was specifing keyword arguments as positional arguments... Same thing for the eigvalsh routine. - I added the pyf wrappers for the new lapack routines in generic_lapack.pyf . They look somewhat different from all other wrappers (they were written 6 years ago without scipy in mind). In the long run I think those fortran routines should get the same kind of wrappers all others have, but I can not do such a rewriting without the guidance of the guy who wrote the generic_XXX.pyf in the first place, which I assume is pearu. In particular the calc_lwork.f module puzzles me. The new lapack routine get the minimum "lwork" and not the optimal one, which should be delivered by calc_lwork, if I only knew how... - I added 129 tests to test_decomp.py using the marvelous nose parametric tests feature. They run pretty fast, and they excercise all relevant lapack routines with most (but not all) combinations of parameters. Is it OK? Should I flag those tests as "less" important so that they are not run by default? If everything is all right with you I would commit the changes tomorrow. I will also continue to provide symeig as a standalone module until scipy 0.7 is released and a corresponding debian package is available. A new release with all changes (and bug-fixes :-)) I made to incorporate it in scipy will be released soon. Let me know, tiziano

2 1

percentileofscore
by josef.pktd＠gmail.com Nov. 17, 2008

Nov. 17, 2008

What is percentileofscore supposed to do? I did not find any good interpretation what the numbers are supposed to mean. >From statistics, I am used to a definition according to the cdf, i.e. fraction of elements weakly smaller than the "score". Instead a strictly smaller definition could be useful, as used eg. in ranking of schools. The current implementation with histogram, does not give results that I can easily interpret. The proposed implementation, has still one error as mentioned by Stefan. It uses the mean when there are multiple elements presents. I looked at 3 cases: * the score element is uniquely present in array * multiple elements in the array are equal to the score * no element in the array is equal to the score I tried out 5 different definitions percentileofscore_proposed: taken from google review with correction percentileofscore_mean: similar to proposed, give mean rank if multiple present This just adds another correction to the proposed version (start index at one instead of zero) percentileofscore_meaninterp: similar to proposed, interpolate if missing percentileofscore_strict: one liner, Fraction(x<score) percentileofscore_weak one liner, Fraction(x<=score) this is what I get: #unique element >>> percentileofscore_proposed([1,2,3,4,5,6,7,8,9,10],4) 30.0 >>> percentileofscore_mean([1,2,3,4,5,6,7,8,9,10],4) 40.0 >>> percentileofscore_meaninterp([1,2,3,4,5,6,7,8,9,10],4) 40.0 >>> percentileofscore_strict([1,2,3,4,5,6,7,8,9,10],4) 30.0 >>> percentileofscore_weak([1,2,3,4,5,6,7,8,9,10],4) 40.0 #multiple elements >>> percentileofscore_proposed([1,2,3,4,4,5,6,7,8,9],4) 35.0 >>> percentileofscore_mean([1,2,3,4,4,5,6,7,8,9],4) 45.0 >>> percentileofscore_meaninterp([1,2,3,4,4,5,6,7,8,9],4) 45.0 >>> percentileofscore_weak([1,2,3,4,4,5,6,7,8,9],4) 50.0 >>> percentileofscore_strict([1,2,3,4,4,5,6,7,8,9],4) 30.0 #missing elements >>> percentileofscore_proposed([1,2,3,5,6,7,8,9,10,11],4) 30.0 >>> percentileofscore_mean([1,2,3,5,6,7,8,9,10,11],4) 30.0 >>> percentileofscore_meaninterp([1,2,3,5,6,7,8,9,10,11],4) 35.0 >>> percentileofscore_weak([1,2,3,5,6,7,8,9,10,11],4) 30.0 >>> percentileofscore_strict([1,2,3,5,6,7,8,9,10,11],4) 30.0 What's the use case for percentileofscore? I just use Fraction(x<=score) or Fraction(x<score) Josef

2 4

status of stats.distributions
by josef.pktd＠gmail.com Nov. 17, 2008

Nov. 17, 2008

I finished with the basic cleanup of scipy.stats.distribution. All generic methods work now, all distributions (except logser.rvs) pass the basic tests for the given parameter values. Test coverage according to figleaf is about 91%. There are some remaining problems: Entropy and fit test for the continuous rv are not included in the test suite. The entropy integration fails for 6 (out of more than 80) continuous distributions and returns nans, I haven't looked at this in detail. Also the entropy test only checks for nan, I didn't find a quick, general test for the numerical correctness of the entropy calculation. The parameter estimation with fit also does not converge very well for for some distribution with sample size up to 10000, and it takes pretty long to run. Some methods defined in the specific distributions don't work correctly, but I did not find any mistakes or I could not find enough information of the statistical properties of these distributions with googling or the bugs are outside of scipy.stats. I replaced these methods by their generic counterparts which work correctly although maybe slower. The skipped methods were renamed by appending "_skip " to the method name. If someone finds the correction, then any help is appreciated. All my tests are currently for chosen parameter values, but I know of a few cases that are broken for some parameter values that are in the valid (but maybe uncommon) range. I did quite a bit of fuzz testing earlier on, but don't have the time now to go over the remaining cases. Tickets 697, 758, 766 and my ticket 745 can be closed now. ticket 620, I would close as don't fix, but I'm not sure how important users would think this is. I also just fixed 769, which looks correct to me. Enhancement tickets 767, 768 are about including limiting cases in distributions. I don't have a strong opinion about these. Is the speed penalty important in this case or not? Are boundary cases important in applications? There are also possibly some details that I missed, e.g. I didn't check return types, but I am basically finished and waiting to see what a more wide spread testing will bring. Josef

2 1

fftw?
by Neal Becker Nov. 16, 2008

Nov. 16, 2008

So what happens to fftw? Will this be available as an additional package? Has any work started on this?

2 1

Preparing to tag 0.7.0b1 in 1 week
by Jarrod Millman Nov. 16, 2008

Nov. 16, 2008

I would like to tag the 0.7.0b1 release on Saturday (11/23) night and announce the beta on Monday (11/24). If you know of any reason that we should delay the first beta, please let me know ASAP. If there is anything that you were planning to commit before the beta, please let me know. All the tests--that aren't known failures or marked to be skipped--pass on my computer. I don't know of any major feature that are waiting to be merged. It would also be great if everyone could test the trunk. Double check to make sure there aren't any regressions. We will be releasing the beta on sourceforge with binaries for Windows and Mac. I will also announce the release on all the mailing lists. I hope to have very wide-spread testing of this beta, so we need to be fairly certain there aren't any unknown regressions. I have a few remaining tasks I would like to see taken care of before the beta: * review, fix, and apply the percentileofscore rewrite and close ticket 560: http://codereview.appspot.com/7913 * work on release notes http://projects.scipy.org/scipy/scipy/browser/trunk/doc * quickly review and bring up-to-date the scipy tutorial http://docs.scipy.org/doc/scipy/reference/tutorial/index.html Please take a look at the release notes and let me know if you see anything that needs to be changed or updated: http://projects.scipy.org/scipy/scipy/milestone/0.7.0 Once we get the beta out we will focus even more on bugs and further improving documentation. Stefan is organizing another scipy sprint as part of the Cape Town Python Users' Group for Saturday, November 29. I will organize a Berkeley sprint as well. So mark your calendars. Thanks, -- Jarrod Millman Computational Infrastructure for Research Labs 10 Giannini Hall, UC Berkeley phone: 510.643.4014 http://cirl.berkeley.edu/

1 0

scipy doc and python 2.4
by Nils Wagner Nov. 16, 2008

Nov. 16, 2008

Hi all, AFAIK scipy 0.7 requires python 2.4 and higher. "finally" is not available in 2.4. Failed to import 'scipy.fftpack.dftfreq': Failed to import 'scipy.io.wavfile.scipy.io.arff.loadarff': Failed to import 'scipy.signal.medfilt2': Failed to import 'scipy.stats.anova': Failed to import 'scipy.stats.fretcher_l': Failed to import 'scipy.stats.fretchet_r': Failed to import 'scipy.stats.meanwhitneyu': Failed to import 'scipy.stats.paired': Failed to import 'scipy.stats.von_mises': touch build/generate-stamp mkdir -p build/html build/doctrees LANG=C sphinx-build -b html -d build/doctrees source build/html Exception occurred: File "/usr/lib/python2.4/site-packages/Sphinx-0.5dev_20081113-py2.4.egg/sphinx/application.py", line 161, in setup_extension mod = __import__(extension, None, None, ['setup']) File "/home/nwagner/svn/scipy/doc/ext/plot_directive.py", line 64 finally: ^ SyntaxError: invalid syntax Nils

2 1

help with special.nbdtrc
by josef.pktd＠gmail.com Nov. 16, 2008

Nov. 16, 2008

Is there a function in scipy.special equivalent to special.nbdtrc that also works when the second coefficient is smaller than one? >>> special.nbdtrc(2, 4.5, 0.42) 0.79201423936000004 >>> special.nbdtrc(2, 0.4 , 0.42) 1.#QNAN the reason: according to http://scipy.org/scipy/scipy/ticket/583 nbinom should also work for all parameters smaller than one, but the survival function uses special.nbdtrc class nbinom_gen(rv_discrete): def _sf(self, x, n, pr): k = floor(x) return special.nbdtrc(k,n,pr) >>> print stats.nbinom.sf(0,4.4, 0.42) 0.96888304 >>> print stats.nbinom.sf(0,0.4, 0.42) 1.#QNAN I can work around the problem, since I can easily calculate the sf without specifying a specific formula for it. Josef

1 0

memoize temporary values and scipy special questions
by josef.pktd＠gmail.com Nov. 16, 2008

Nov. 16, 2008

2 observations and 2 questions observation 1: gausshyper distribution uses the generic way of generating random variables, which is pretty slow for larger samples. the _pdf recalculates the same temporary value each time again. Cinv is the normalization constant for the pdf, which depends on the parametersn(a, b, c, z) but not on the sample point x class gausshyper_gen(rv_continuous): def _pdf(self, x, a, b, c, z): Cinv = gam(a)*gam(b)/gam(a+b)*special.hyp2f1(c,a,a+b,-z) return 1.0/Cinv * x**(a-1.0) * (1.0-x)**(b-1.0) / (1.0+z*x)**c Note: gam = special.gamma observation 2: special.hyp2f1 seems to have problems over some parameter range question 1: When calculating the pdf for different x but same parameters, then the normalization constant could be temporarily stored somewhere as in memoizing or as a lazy attribute. Is there a precedence and pattern for this in numpy/scipy, or do I try whatever might work in the individual case, e.g. create a cache for repeated calculation of the same values? question 2: In the above case, the normalization constant could also be obtained through direct numerical integration (over x in interval [0,1]). What are the relative merits of using scipy.special versus numerical integration, in terms of speed, accuracy, and range for which correct answers are produced? These are just some ideas I had, while I playing around with gausshyper to see why it doesn't work very well (it is slow and convergence of fit method to true parameters is not very good). I haven't profiled anything and rvs is additionally slow because the generic way of generating random variables when only the pdf is given is pretty roundabout: pdf->cdf->ppf->rvs If somebody can provide some information on this, then I don't have to figure it out myself. Thanks, Josef

1 0

about http://www.scipy.org/LoktaVolterraTutorial
by Johann Cohen-Tanugi Nov. 15, 2008

Nov. 15, 2008

Hello, Icame across this web page browsing my stored emails from the mailing list, and I wondered how I could have found it other than by receiving the email that pointed to it.... Should it be moved inside the cookbook? Are there many other examples like that that are hard to find in the web site? thanks, Johann

2 1

docs scipy.stats has incorrect introduction - instances of discrete distributions
by josef.pktd＠gmail.com Nov. 15, 2008

Nov. 15, 2008

The docs for stats distributions does not mention discrete distributions, i.e change and add something like this: Each included distribution is an instance of the class rv_continous or of rv_discrete: .. autosummary:: :toctree: generated/ rv_discrete rv_discrete.pmf rv_discrete.cdf rv_discrete.sf rv_discrete.ppf rv_discrete.isf rv_discrete.stats additional methods, moment, entropy, nnlf, fit, could also be added (but I haven't finished testing all of those) Josef

2 1