related -- what is the canonical way to trigger build scipy (out of git) without any optimization flags for gcc? tried FFLAGS='-O0 -g' CXXFLAGS='-O0 -g' CFLAGS='-O0 -g' python-dbg setup.py build_ext --inplace --debug but still getting e.g. Fortran fix compiler: /usr/bin/gfortran -Wall -ffixed-form -fno-second-underscore -Wall -fno-second-underscore -O0 -g -O3 -funroll-loops where it gets overriden On Fri, 30 Mar 2012, josef.pktd@gmail.com wrote:

On Fri, Mar 30, 2012 at 11:43 AM, Yaroslav Halchenko <lists@onerussian.com> wrote:

...

I have reported this issue some time ago on Debian http://bugs.debian.org/cgi-bin/bugreport.cgi?bug=653948 and then forgot about it until now that I ran into it again.

...

With recent scipy (including git master 0fbfdbc)

...

scipy.stats.ksone.fit

...

seems to stall (never return) on big-endian boxes and return (1.0, nan, nan) on x86.  It seems that it is  working correctly with scipy 0.7.2 (as it is now in Debian stable).

...

So -- is that anything known?

No, never seen before.

I didn't think anyone would ever fit ksone, It doesn't even have a pdf defined in the source. It's mainly included for Kolmogorov-Smirnov test.

I don't know how to interpret gdp. Is it clear that it is cephes_smirnov? Do you know the values that trigger it? If not, you could put a print in the ksone_gen._cdf to see where it gets stuck. (it might print a lot)

I think it's a bug for scipy.stats that ksone doesn't define the support boundary .b, it looks like it should be (0,1) (default a=0 looks ok.) There will be lots of nans during fit(). I don't think the generic fit is smart enough to figure out non-nan or non-inf starting values for your dataset. (There is a ticket to avoid inf in starting values but it hasn't been included yet.)

Josef

...

snippet to replicate:

...

import numpy as np import scipy.stats as ss

...

d = np.array([-0.18879233,  0.15734249,  0.18695107,  0.27908787, -0.248649,              -0.2171497 ,  0.12233512,  0.15126419,  0.03119282,  0.4365294 ,              0.08930393, -0.23509903,  0.28231224, -0.09974875, -0.25196048,              0.11102028,  0.1427649 ,  0.10176452,  0.18754054,  0.25826724,              0.05988819,  0.0531668 ,  0.21906056,  0.32106729,  0.2117662 ,              0.10886442,  0.09375789,  0.24583286, -0.22968366, -0.07842391,              -0.31195432, -0.21271196,  0.1114243 , -0.13293002,  0.01331725,              -0.04330977, -0.09485776, -0.28434547,  0.22245721, -0.18518199,              -0.10943985, -0.35243174,  0.06897665, -0.03553363, -0.0701746 ,              -0.06037974,  0.37670779, -0.21684405])

...

print "Fitting now" print ss.ksone.fit(d)

...

-- =------------------------------------------------------------------= Keep in touch                                     www.onerussian.com Yaroslav Halchenko                 www.ohloh.net/accounts/yarikoptic _______________________________________________ SciPy-Dev mailing list SciPy-Dev@scipy.org http://mail.scipy.org/mailman/listinfo/scipy-dev

---

SciPy-Dev mailing list SciPy-Dev@scipy.org http://mail.scipy.org/mailman/listinfo/scipy-dev

-- =------------------------------------------------------------------= Keep in touch www.onerussian.com Yaroslav Halchenko www.ohloh.net/accounts/yarikoptic

yikes -- thanks! -- I simply had no clue (ashamed since I should have known this one) for clarity -- FOPT is numpy/scipy specific and doesn't come from distutils interestingly looks the code handling of OPT in distuils 2.7: if compiler.compiler_type == "unix": (cc, cxx, opt, cflags, opt, extra_cflags, basecflags, ccshared, ldshared, so_ext) = \ get_config_vars('CC', 'CXX', 'OPT', 'CFLAGS', 'OPT', 'EXTRA_CFLAGS', 'BASECFLAGS', 'CCSHARED', 'LDSHARED', 'SO') I guess OPT is queried twice just to be sure ;) Thanks again! On Fri, 30 Mar 2012, Pauli Virtanen wrote:

30.03.2012 18:48, Yaroslav Halchenko kirjoitti:

...

related -- what is the canonical way to trigger build scipy (out of git) without any optimization flags for gcc?

FOPT='-O0 -g3' OPT='-O0 -g3' python setup.py build_ext --inplace

This is inherited from Python's distutils, and so should work for any Python package. -- =------------------------------------------------------------------= Keep in touch www.onerussian.com Yaroslav Halchenko www.ohloh.net/accounts/yarikoptic

would you by any chance also know how to make those statements verbose? (i.e. to list actual invocations to the compiler)? ... gcc: scipy/special/cephes/kolmogorov.c ... On Fri, 30 Mar 2012, Pauli Virtanen wrote:

30.03.2012 18:48, Yaroslav Halchenko kirjoitti:

...

related -- what is the canonical way to trigger build scipy (out of git) without any optimization flags for gcc?

FOPT='-O0 -g3' OPT='-O0 -g3' python setup.py build_ext --inplace

This is inherited from Python's distutils, and so should work for any Python package. -- =------------------------------------------------------------------= Keep in touch www.onerussian.com Yaroslav Halchenko www.ohloh.net/accounts/yarikoptic

ok -- here is the reason: (gdb) print e $27 = nan(0x8000000000000) (gdb) print ((double) n * (1.0 - e)) $26 = nan(0x100000001) (gdb) print (floor ((double) n * (1.0 - e))) $25 = 2146435073 and then it goes into the loop with (gdb) print nn $32 = 2147483647 so it might eventually return (didn't wait long enough) but that would take a while ;) what confuses ignorant me is why floor on sparc * why it doesn't handle nan correctly as manpage says: If x is integral, +0, -0, NaN, or an infinity, x itself is returned. * why it doesn't return double for double as manpage (looking at x86 box but that should remain valid I guess) says: double floor(double x); float floorf(float x); long double floorl(long double x); but returns float : (gdb) print sizeof(floor ((double) n * (1.0 - e))) $1 = 4 * and I guess I need to learn about different types of nans (gdb) print e $27 = nan(0x8000000000000) (gdb) print sizeof(e) $28 = 8 (gdb) print sizeof(NPY_NAN) $29 = 8 (gdb) print NPY_NAN $31 = nan(0x100000001) On Fri, 30 Mar 2012, Yaroslav Halchenko wrote:

quite positive (may be not straight in cephes_smirnov but nowhere near python ;) ) -- I will rebuild ones more -- previously gdb was confused with source file was changed during build... hopefully this time it would be cleaner -- it takes a bit though on that sparc ;)

On Fri, 30 Mar 2012, Pauli Virtanen wrote:

...

...

With recent scipy (including git master 0fbfdbc)

...

...

scipy.stats.ksone.fit

...

...

seems to stall (never return) on big-endian boxes and return (1.0, nan, nan) on x86. It seems that it is working correctly with scipy 0.7.2 (as it is now in Debian stable).

...

Are you sure it really hangs in cephes_smirnov? I'd run the code with 'python -m pdb xxx.py' and step through the execution to be sure that the problem is not on the Python level.

...

The only change made to cephes_smirnov since 0.7.2 is to make it return a NAN in undefined cases, so it seems unlikely the problem is there. To check if the problem is in cephes_smirnov, set a breakpoint to cephes_smirnov in GDB and follow the execution. -- =------------------------------------------------------------------= Keep in touch www.onerussian.com Yaroslav Halchenko www.ohloh.net/accounts/yarikoptic

well -- imho it should have not even got to that point if e is NaN. Just started rebuilding with following patch: - if (n <= 0 || e < 0.0 || e > 1.0) + # This comparison should assure returning NaN whenever + # e is NaN itself. In original || form it would proceed + if !(n > 0 && e >= 0.0 && e <= 1.0) return (NPY_NAN); On eri, 30 Mar 2012, Pnuli Virtanen wrote:

Hi,

30.03.2012 22:37, Yaroslav Halchenko kirjoitti:

...

ok -- here is the reason:

...

(gdb) print e $27 = nan(0x8000000000000) (gdb) print ((double) n * (1.0 - e)) $26 = nan(0x100000001) (gdb) print (floor ((double) n * (1.0 - e))) $25 = 2146435073

Thanks a lot:

kolmogorov.c:41 nn = (int) (floor ((double) n * (1.0 - e)));

This is just wrong -- if `e` happens to be NAN, the result of the integer cast is unspecified (as per C99 standard).

Probably there should be a blanket input finiteness check in most of the routines. I'd be willing to bet that this is not the only bug of this type in there. Detecting these automatically would require some control flow analysis, so I guess the only option is to go through each function manually :/

Pauli

_______________________________________________ SciPy-Dev mailing list SciPy-Dev@scipy.org http://mail.scipy.org/mailman/listinfo/scipy-dev

-- =------------------------------------------------------------------= Keep in touch www.onerussian.com Yaroslav Halchenko www.ohloh.net/accounts/yarikoptic

On Fri, Mar 30, 2012 at 9:50 PM, Yaroslav Halchenko <lists@onerussian.com> wrote:

yeap -- that

-  if (n <= 0 || e < 0.0 || e > 1.0) +  /* This comparison should assure returning NaN whenever +     e is NaN itself.  In original || form it would proceed */ +  if (!(n > 0 && e >= 0.0 && e <= 1.0))

resolved the stalling issue and now I am getting the same (1.0, nan, nan) as on x86 ... sent pull request

https://github.com/scipy/scipy/pull/187

patch is attached here as well

so next I guess is to make it return sensible values for the .fit as it did before? ;)

sensible? or starting values? Fitting now [ 1. 0. 1.]

...

...

import scipy scipy.__version__ '0.7.2'

np.__version__ '1.4.1'

Josef

On Fri, 30 Mar 2012, Yaroslav Halchenko wrote:

...

well -- imho it should have not even got to that point if e is NaN.  Just started rebuilding with following patch:

...

-  if (n <= 0 || e < 0.0 || e > 1.0) +  # This comparison should assure returning NaN whenever +  # e is NaN itself.  In original || form it would proceed +  if !(n > 0 && e >= 0.0 && e <= 1.0)      return (NPY_NAN);

-- =------------------------------------------------------------------= Keep in touch                                     www.onerussian.com Yaroslav Halchenko                 www.ohloh.net/accounts/yarikoptic _______________________________________________ SciPy-Dev mailing list SciPy-Dev@scipy.org http://mail.scipy.org/mailman/listinfo/scipy-dev

On Fri, 30 Mar 2012, josef.pktd@gmail.com wrote:

...

...

so next I guess is to make it return sensible values for the .fit as it did before? ;) sensible? or starting values?

if starting values are the most sensible -- then yeap -- them ;) if I ask to 'fit' something, getting some fit is better than getting no fit (as NaNs in output suggest)

Sorry if Debian is getting some noise from my side today. I have problems paying attention to reply versus reply-all.

well -- it was my fault anyways trying to kill two birds at once ;) Debian BTS would survive that just fine, no worries -- even might appreciate having more context for the report ;) -- =------------------------------------------------------------------= Keep in touch www.onerussian.com Yaroslav Halchenko www.ohloh.net/accounts/yarikoptic

Probably you are right Josef -- especially since I am only distantly familiar with KS test -- but lets keep the dialog open a bit longer ;) :

But what's the point in fitting ksone?

for me it was just that it has .fit() ;) You might recall (I believe I appeared on the list long ago with similar whining and that is how we got introduced to each other) our evil/silly function in PyMVPA match_distributions which simply tries to choose the best matching distribution given the data -- that is the reason how ksone got involved

...

if starting values are the most sensible -- then yeap -- them ;) if I ask to 'fit' something, getting some fit is better than getting no fit (as NaNs in output suggest)

getting the starting values back doesn't mean that you have "some" fit.

If my brief playing with it today is correct, then the starting values don't make sense, for example you have points outside of the support of the distribution with estimated parameters (if you have negative values in the sample)

NaN would be better, then at least you know it doesn't make sense.

1. to me the big question became: what ARE the logical values here? followed docstring/example on http://docs.scipy.org/doc/scipy/reference/generated/scipy.stats.ksone.html -- got NaNs then given that In [44]: ksone.a, ksone.b Out[44]: (0.0, inf) I still failed to get any sensible fit() for positive values or even for its own creation, e.g. ss.ksone.fit(ss.ksone(5).rvs(size=100)) results in bulk of warnings and then (1.0, nan, nan). Looking in detail -- rvs is happily generating NaNs (especially for small n's). b. Also the range of sensible values of the parameter n isn't specified anywhere for KS test newbies like me, which I guess adds the confusion:

support of the sample would help. I have no idea about good starting values for the shape parameter (n is sample size for kstest)

aga -- so the 'demo' value of 0.9 indeed makes no sense ;) Might be worth adjusting somehow? 2. BTW -- trying to familiarize myself with the distribution plotted its pdf, e.g.: x = np.linspace(0, 3, 1000); plt.plot(x, ksone(10).pdf(x)) and it looks weirdish: http://www.onerussian.com/tmp/ksone-ns.png in that it is not smooth and my algebra-forgotten eyes do not see obvious points with no 2nd derivative of cdf given on http://en.wikipedia.org/wiki/Kolmogorov_Smirnov Also why ksone.b is inf -- shouldn't it be 1? -- =------------------------------------------------------------------= Keep in touch www.onerussian.com Yaroslav Halchenko www.ohloh.net/accounts/yarikoptic

On Sat, Mar 31, 2012 at 12:02 PM, <josef.pktd@gmail.com> wrote:

On Sat, Mar 31, 2012 at 11:15 AM, Yaroslav Halchenko <lists@onerussian.com> wrote:

...

Probably you are right Josef -- especially since I am only distantly familiar with KS test -- but lets keep the dialog open a bit longer ;) :

...

But what's the point in fitting ksone?

for me it was just that it has .fit() ;)    You might recall (I believe I appeared on the list long ago with similar whining and that is how we got introduced to each other) our evil/silly function in PyMVPA match_distributions which simply tries to choose the best matching distribution given the data -- that is the reason how ksone got involved

I remember and if I remember correctly, then I recommended using a blacklist of distributions to avoid.

The last time I looked at the source of pymvpa, you used all distribution in the fit and then reported the best fitting ones. At the bottom of this ranking there should be some distributions that will (almost) never be a good match because fit doesn't work for them. The only time you see how bad they are is in extreme cases like going off to neverland.

...

...

...

if starting values are the most sensible -- then yeap -- them ;) if I ask to 'fit' something, getting some fit is better than getting no fit (as NaNs in output suggest)

...

getting the starting values back doesn't mean that you have "some" fit.

...

If my brief playing with it today is correct, then the starting values don't make sense, for example you have points outside of the support of the distribution with estimated parameters (if you have negative values in the sample)

...

NaN would be better, then at least you know it doesn't make sense.

1. to me the big question became: what ARE the logical values here?

if you look at my second message above, you see some examples, where fit returns numbers. I didn't check how good they are.

...

followed docstring/example on http://docs.scipy.org/doc/scipy/reference/generated/scipy.stats.ksone.html -- got NaNs

then given that

In [44]: ksone.a, ksone.b Out[44]: (0.0, inf)

I still failed to get any sensible fit() for positive values or even for its own creation, e.g.

ss.ksone.fit(ss.ksone(5).rvs(size=100))

...

...

...

rv = stats.ksone(50).rvs(size=1000) plt.hist(rv, bins=30, normed=True, cumulative=True) x = np.linspace(0, rv.max(), 1000); plt.plot(x, stats.ksone(50).cdf(x)) plt.show()

...

...

...

stats.ksone.fit(rv, 100, loc=-0.01, scale=1)

(181.94347728444751, -3.8554246919087482e-05, 1.9277121337713585)

...

...

...

stats.ksone.fit(rv, 10, loc=-0.01, scale=1) (13.999896396912176, -0.010783712808254388, 0.57818285700694405)

...

results in bulk of warnings and then (1.0, nan, nan).

Looking in detail -- rvs is happily generating NaNs (especially for small n's).

b. Also the range of sensible values of the parameter n isn't specified anywhere for KS test newbies like me, which I guess adds the confusion:

...

support of the sample would help. I have no idea about good starting values for the shape parameter (n is sample size for kstest)

aga -- so the 'demo' value of 0.9 indeed makes no sense ;)  Might be worth adjusting somehow?

2.

BTW -- trying to familiarize myself with the distribution plotted its pdf, e.g.:

x = np.linspace(0, 3, 1000); plt.plot(x, ksone(10).pdf(x))

and it looks weirdish: http://www.onerussian.com/tmp/ksone-ns.png in that it is not smooth and my algebra-forgotten eyes do not see obvious points with no 2nd derivative of cdf given on http://en.wikipedia.org/wiki/Kolmogorov_Smirnov

IIRC (no time to check again right now): ksone is, I think, a small sample distribution, kstwobign is the distribution of the max/sup of a Brownian Bridge, which is the asymptotic distribution for Kolmogorov-Smirnov

I needed to check the source: the c source says smirnov is the distribution for one-sided test, and I had forgotten that I had added the one-sided option to kstest. algorithm also used by R: " The formula of Birnbaum & Tingey (1951) is used for the one-sample one-sided case." http://www.jstor.org/stable/2236929 ks_2samp is still missing the one-sided options Josef anderson darling is most of the time more powerful than KS https://github.com/aarchiba/kuiper no license

as distribution we are mainly interested in cdf and ppf (both look reasonably good in  a plot), and mainly in the right tail ksone looks like a piecewise approximation, where they didn't care much about the lower part.

(I'm a bit rushed right now so there might be parts missing in my reply)

Josef

...

Also why ksone.b is inf -- shouldn't it be 1?

-- =------------------------------------------------------------------= Keep in touch                                     www.onerussian.com Yaroslav Halchenko                 www.ohloh.net/accounts/yarikoptic _______________________________________________ SciPy-Dev mailing list SciPy-Dev@scipy.org http://mail.scipy.org/mailman/listinfo/scipy-dev