PEP: Improving the basic statistical functions in Scipy
Hi, I do apologize in advance if this is considered inappropriate but my goal is to advance the stats capabilities in Scipy. I do recognize that there is a large chunk of excellent work so part of this is simply to ensure that there is adequate documentation and tests. The following is my attempt of a PEP to provide some direction on how to improve the basic statistical functions within Scipy. I do have a list of the individual functions and the arguments involved but I decided it was in appropriate to attach it here. Probably the main aspect that I would like feedback is on whether or not there should be a single interface to these basic statistical functions. Thanks Bruce PEP: Improving the basic statistical functions in Scipy Authors: Bruce Southey Created: 26-Feb-2009 Abstract ======== This current PEP is orientated towards addressing the fundamental problems with the basic statistical functions in Scipy. The outcome is to provide Scipy with a consistent, well-tested and documented set of basic statistical functions that are available to different array types. Motivation ======== This PEP addresses the basic statistical functions available in the stats component of Scipy. These functions are defined in the following files: stats.py – Defines many statistical functions and imports statlib. morestats.py – Adds additional statistical functions to stats.py _support.py - Defines the functions used in stats.py but also is a circular because it also imports stats mstats.py – Just imports functions from mstats_basic.py and mstats_extras.py mstats_basic.py – Defines statistical functions for masked arrays mstats_extras.py– Defines additional statistical functions for masked arrays In total there are 178 unique functions defined in these files, some of which are private or internal and some have the same name but are defined slightly differently between standard and masked arrays. A list of theses functions is available. While the functions are defined for standard arrays and masked arrays, not all functions are available for both array types. For example, the majority of functions defined in stats.py for standard arrays are available in mstats_basic.py. But none of the standard arrays functions defined in morestats.py are available for masked arrays. Also none of these are functions are directly supported for other array types available to Scipy such record arrays, record arrays that contain masked data and sparse arrays. Specification ======== 1) Provide the same basic statistical functions with the same arguments for standard and masked arrays. 2) Utilize a single interface. For example, the gmean function (note _chk_asarray is defined differently): stats.py: def gmean(a, axis=0): a, axis = _chk_asarray(a, axis) log_a = np.log(a) return np.exp(log_a.mean(axis=axis)) mstats_basic.py: def gmean(a, axis=0): a, axis = _chk_asarray(a, axis) log_a = ma.log(a) return ma.exp(log_a.mean(axis=axis)) Rather a single function can be defined as: def gmean(a, axis=0): log_a = np.log(a) return np.exp(log_a.mean(axis=axis)) import numpy as np import numpy.ma as ma X=[1,2,3,4,5] a=np.array(X) m=ma.array(X, mask=[0,0,0,0,0]) np.exp((np.log(X).mean())) #2.6051710846973517 np.exp((np.log(a).mean())) #2.6051710846973517 np.exp((np.log(m).mean())) #2.6051710846973517 3) Depreciation and removal of unnecessary functions such as linregress. 4) Cleanup styles issues including: a) White space usage b) Consistent arguments such as 'a' vs 'x' and the usage of *args c) Uniquely identifying functions. i) rootfunc and tempfunc defined two and three times, respectively, in morestats.py but have different arguments. ii) makestr is defined twice in _support.py, once a main function and once as a subfunction of printcc. 5) Ensure info.py is complete and correct. 6) Improve the documentation of basic statistical functions in connection with the Scipy documentation Marathon (http://www.scipy.org/Developer_Zone/DocMarathon2008) 7) Improve the tests of the basic statistical functions: i) All functions should have at least have basic test coverage that indicates whether or not it is functional. ii) Important functions should have tests that include unexpected elements like Nan's, positive and negative infinity and other unexpected inputs. iii) Ideally there should be tests that check the function accuracy. 8) Extension of the functions to other array types available to Scipy such record arrays, record arrays that contain masked data and sparse arrays. Perhaps beyond the scope of this PEP. Backwards Compatibility ======== There is no guarantee that the outcome will maintain complete backwards compatibility because a consistent API is required across different array types. However, any changes to existing APIs must be justified such as ensuring the same keywords between functions for different array types.
On Thu, Feb 26, 2009 at 04:26:39PM -0600, Bruce Southey wrote:
1) Provide the same basic statistical functions with the same arguments for standard and masked arrays.
Sounds great.
4) Cleanup styles issues including: a) White space usage b) Consistent arguments such as 'a' vs 'x' and the usage of *args c) Uniquely identifying functions. i) rootfunc and tempfunc defined two and three times, respectively, in morestats.py but have different arguments. ii) makestr is defined twice in _support.py, once a main function and once as a subfunction of printcc.
Excellent.
5) Ensure info.py is complete and correct.
Great.
6) Improve the documentation of basic statistical functions in connection with the Scipy documentation Marathon (http://www.scipy.org/Developer_Zone/DocMarathon2008)
Go for it!
7) Improve the tests of the basic statistical functions: i) All functions should have at least have basic test coverage that indicates whether or not it is functional. ii) Important functions should have tests that include unexpected elements like Nan's, positive and negative infinity and other unexpected inputs. iii) Ideally there should be tests that check the function accuracy.
Good. Hell, I can't say more. If you can get just half of what you have listed up there, it would be fantastic. I'll try to review your work (which I expect you will be putting up on a code review site, as discussed previously), but I don't promise anything: it can be hard for me to find time to sit down and do something serious on top of my current workload. Cheers, Gaël
Hi Bruce 2009/2/27 Bruce Southey <bsouthey@gmail.com>:
The following is my attempt of a PEP to provide some direction on how to improve the basic statistical functions within Scipy. I do have a list of the individual functions and the arguments involved but I decided it was in appropriate to attach it here.
Thank you for all the thoughtful suggestions. I think some of these issues can already be turned into tickets, i.e. API inconsistencies, missing test coverage, broken docs, etc. It might be useful to do so, so that we know what needs to be done next. Regards Stéfan
I think a discussion for a roadmap for stats will be very useful. Currently my priority is still your point 7 iii) Ideally there should be tests that check the function accuracy. I consider this the main point of almost all my work on stats. And there are still some incorrect parts left. The next part for the current code base, that I think about, was to evaluate function whether they are ok, can be generalized, e.g. dimension, or are trivial and should be removed. Next are changes in the interface and combining or comparing mstats and stats. Here, I don't have a clear opinion yet of how far we can or want to consistently generalize all statistical functions to the different type of arrays. In many cases I looked at, the masked array version looked sufficiently different that I would be reluctant to merge them. One radical alternative would be to depreciate stats.stats and expand mstats, since it is already better designed to handle different array types. But I like the "simple" versions in stats, and I'm curious about any speed difference. But general tools to interface to different array types would be useful and should be carefully designed, e.g. function like ols that have a plain ndarray core, but can access the data from structured arrays and masked arrays. After, the changes to the current statistical function, I was considering areas of statistics that have partial but incomplete coverage. Non-parametric tests are well represented, and I have some extension for tests for discrete distributions. I think ANOVA, which I never used myself, has a very incomplete collection, which, I guess is a historical accident since Gary Strangman had, I think more ANOVA functions that are not included in stats. So instead of having a laundry list of functions, (some of which don't seem to have been used for years), I would prefer at least a conceptional grouping around statistical topics. Regression of course is currently MIA. The next large interface issue, especially for enhancements, is whether to use functions or proper classes. I think for some statistical analysis the current statistical function, once cleaned up, work fine. However, even R returns result classes (or whatever their equivalent is) for every statistical test, while in python we use matlab style functions. This will change when models will be included again. I have a list of functions that have no test coverage, a list (not written down) of functions that have bug suspects or known bugs, and it would be useful to get a wider opinion about which functions and interfaces are important Working on the list of functions on the wiki page maybe simpler for collecting comments than going through the statistical review in trac. Overall, I think there is still a lot of work to do before I start to worry about white space issues. Josef
On Thu, Feb 26, 2009 at 17:47, <josef.pktd@gmail.com> wrote:
After, the changes to the current statistical function, I was considering areas of statistics that have partial but incomplete coverage. Non-parametric tests are well represented, and I have some extension for tests for discrete distributions. I think ANOVA, which I never used myself, has a very incomplete collection, which, I guess is a historical accident since Gary Strangman had, I think more ANOVA functions that are not included in stats.
It's no accident. I removed it. It was a big monster function. Gary had a big disclaimer on it saying that it basically worked for the use case he had at the time, but it was far from a good general implementation. It printed things without being asked. No one knew if it actually worked. Gary removed it from later versions of his code, etc. Some of use decided that it was better for someone with an interest in ANOVA to reimplement it from scratch. -- Robert Kern "I have come to believe that the whole world is an enigma, a harmless enigma that is made terrible by our own mad attempt to interpret it as though it had an underlying truth." -- Umberto Eco
josef.pktd@gmail.com wrote:
I think a discussion for a roadmap for stats will be very useful.
Currently my priority is still your point 7 iii) Ideally there should be tests that check the function accuracy.
I consider this the main point of almost all my work on stats. And there are still some incorrect parts left.
Yes, that is why I added it.
The next part for the current code base, that I think about, was to evaluate function whether they are ok, can be generalized, e.g. dimension, or are trivial and should be removed.
I agree as I do think some are a consequence of the porting process and have never received the appropriate followup over time.
Next are changes in the interface and combining or comparing mstats and stats. Here, I don't have a clear opinion yet of how far we can or want to consistently generalize all statistical functions to the different type of arrays. In many cases I looked at, the masked array version looked sufficiently different that I would be reluctant to merge them. One radical alternative would be to depreciate stats.stats and expand mstats, since it is already better designed to handle different array types. But I like the "simple" versions in stats, and I'm curious about any speed difference.
The main issue that prevents me from going further with this aspect! I do not find it that radical at all to suggest that as I am for just using masked arrays because I do not perceive a speed difference. (Okay I am perhaps unusual in that I work with large datasets and complex models so differences of a few seconds are not that meaningful to me.) It would be less work to convert the missing as there are about 85 functions missing from masked.
But general tools to interface to different array types would be useful and should be carefully designed, e.g. function like ols that have a plain ndarray core, but can access the data from structured arrays and masked arrays.
After, the changes to the current statistical function, I was considering areas of statistics that have partial but incomplete coverage. Non-parametric tests are well represented, and I have some extension for tests for discrete distributions. I think ANOVA, which I never used myself, has a very incomplete collection, which, I guess is a historical accident since Gary Strangman had, I think more ANOVA functions that are not included in stats. So instead of having a laundry list of functions, (some of which don't seem to have been used for years), I would prefer at least a conceptional grouping around statistical topics. Regression of course is currently MIA.
Even after Robert's reply on that, stats.py at least still has linregress (simple regression with one variable) and glm that address these. However, there is a strong case that both of these should also be removed in favor of a better approach. I agree that doing things like general linear models (eg regression and ANOVA assuming normality), generalized linear models and such need a careful design that integrates where possible existing solutions. Even SAS has different procedures and different modules are available for R to do these. But must be a separate discussion.
The next large interface issue, especially for enhancements, is whether to use functions or proper classes. I think for some statistical analysis the current statistical function, once cleaned up, work fine. However, even R returns result classes (or whatever their equivalent is) for every statistical test, while in python we use matlab style functions.
This will change when models will be included again.
Excellent!
I have a list of functions that have no test coverage, a list (not written down) of functions that have bug suspects or known bugs, and it would be useful to get a wider opinion about which functions and interfaces are important Working on the list of functions on the wiki page maybe simpler for collecting comments than going through the statistical review in trac.
I agree that we need to address what functions we really need and what interface is required. From that we can address the required tests and documentation.
Overall, I think there is still a lot of work to do before I start to worry about white space issues.
Yeah, I just figured that we should correct any of these coding styles issues on the way. Thanks for all the comments, Bruce
On Fri, Feb 27, 2009 at 10:03 AM, Bruce Southey <bsouthey@gmail.com> wrote:
josef.pktd@gmail.com wrote:
I think a discussion for a roadmap for stats will be very useful.
Currently my priority is still your point 7 iii) Ideally there should be tests that check the function accuracy.
I consider this the main point of almost all my work on stats. And there are still some incorrect parts left.
Yes, that is why I added it.
The next part for the current code base, that I think about, was to evaluate function whether they are ok, can be generalized, e.g. dimension, or are trivial and should be removed.
I agree as I do think some are a consequence of the porting process and have never received the appropriate followup over time.
Next are changes in the interface and combining or comparing mstats and stats. Here, I don't have a clear opinion yet of how far we can or want to consistently generalize all statistical functions to the different type of arrays. In many cases I looked at, the masked array version looked sufficiently different that I would be reluctant to merge them. One radical alternative would be to depreciate stats.stats and expand mstats, since it is already better designed to handle different array types. But I like the "simple" versions in stats, and I'm curious about any speed difference.
The main issue that prevents me from going further with this aspect!
I do not find it that radical at all to suggest that as I am for just using masked arrays because I do not perceive a speed difference. (Okay I am perhaps unusual in that I work with large datasets and complex models so differences of a few seconds are not that meaningful to me.) It would be less work to convert the missing as there are about 85 functions missing from masked.
I don't know what the current range of use cases for stats is. But for example in matlab, I have some ols estimation in an innerloop where I wouldn't want much overhead. But in this case, it would always be possible to go back to raw linalg.lstsq. The other disadvantage for me is that it is much easier to write functions that work for plain arrays, since I'm not working with masked/missing data. It's ok if the handling of different array types can be done in the interface of the function, but translating some statistical formulas into code or porting it from another language will be more difficult for me if I have to worry about missing values all the time. An example that I looked at recently, is statistical analysis of panel data, with a balanced panel the linear algebra and matrix operations are much easier than with an unbalanced panel. What I would like to do, but didn't have the time yet is to run the tests for stats.stats on stats.mstats. This way even if we would have some duplicate functions, we would have some cross check that they are consistent, and it would be a reminder for bug fixing also the other version.
But general tools to interface to different array types would be useful and should be carefully designed, e.g. function like ols that have a plain ndarray core, but can access the data from structured arrays and masked arrays.
After, the changes to the current statistical function, I was considering areas of statistics that have partial but incomplete coverage. Non-parametric tests are well represented, and I have some extension for tests for discrete distributions. I think ANOVA, which I never used myself, has a very incomplete collection, which, I guess is a historical accident since Gary Strangman had, I think more ANOVA functions that are not included in stats. So instead of having a laundry list of functions, (some of which don't seem to have been used for years), I would prefer at least a conceptional grouping around statistical topics. Regression of course is currently MIA.
Even after Robert's reply on that, stats.py at least still has linregress (simple regression with one variable) and glm that address these. However, there is a strong case that both of these should also be removed in favor of a better approach.
I don't really count linregress as a "serious" statistical function, since the restriction to one explanatory variable has no computational advantage if we have access to linalg. Similarly, I don't know what the purpose of pointbiserial is, if you can use np.corrcoef for the correlation coefficient or stats.pearsonr for the p-values. My impression is that these are historical functions, when there was no easy access to fast computers and full matrix and array packages. stats.glm is a bit of a misnomer it is just a t-test for the regression on one dummy variable, not an estimator. But again I don't see an advantage compared to ols with multivariate regressors and dummy variables.
I agree that doing things like general linear models (eg regression and ANOVA assuming normality), generalized linear models and such need a careful design that integrates where possible existing solutions. Even SAS has different procedures and different modules are available for R to do these. But must be a separate discussion.
The next large interface issue, especially for enhancements, is whether to use functions or proper classes. I think for some statistical analysis the current statistical function, once cleaned up, work fine. However, even R returns result classes (or whatever their equivalent is) for every statistical test, while in python we use matlab style functions.
This will change when models will be included again.
There are still bugs in it, and test coverage is still low. If anyone wants to help in the review, bug hunting or adding test the current version is in nipy at https://code.launchpad.net/~nipy-developers/nipy/trunk-josef-models
Excellent!
I have a list of functions that have no test coverage, a list (not written down) of functions that have bug suspects or known bugs, and it would be useful to get a wider opinion about which functions and interfaces are important Working on the list of functions on the wiki page maybe simpler for collecting comments than going through the statistical review in trac.
I agree that we need to address what functions we really need and what interface is required. From that we can address the required tests and documentation.
Overall, I think there is still a lot of work to do before I start to worry about white space issues.
Yeah, I just figured that we should correct any of these coding styles issues on the way.
I'm slowly getting used to the formatting requirements, and at least during code changes, I try to stick to it.
Thanks for all the comments, Bruce
Josef
One more issue for the design of statistical function is the availability of using weights. I was looking at calculating weighted means and variances and so on, but the current situation doesn't look very good. There is np.average and the new curvefit allows for weights. http://scipy.org/scipy/scipy/ticket/604 has a full set of statistical functions using weights, but I couldn't make up my mind about how this should fit in. Many of the functions are very short wrappers and would increase the number of functions without necessarily a big benefit. But an efficient implementation of statistical functions that allow weights would make the use of dummy variables and the conversion of masked arrays to use the weighted functions easier, (use mask as dummy variable for the weight.) This won't help for all cases where masked arrays are used, but looking at specific functions and coming up with a good general design would be very useful. Josef
josef.pktd@gmail.com wrote: [snip]
What I would like to do, but didn't have the time yet is to run the tests for stats.stats on stats.mstats. This way even if we would have some duplicate functions, we would have some cross check that they are consistent, and it would be a reminder for bug fixing also the other version.
Okay, I do not know how to get timeit to work with numpy/scipy but this is not how I would like it to be. But I managed somehow to (unfairly) compare the geometric means function (gmean) using this code: import timeit stand_t=timeit.Timer('scipy.stats.stats.gmean(X, axis=xs)', 'import numpy, scipy.stats.stats; X=numpy.random.gamma(shape=2, scale=1, size=(1,10)); xs=None').timeit(1000) masked_t=timeit.Timer('scipy.stats.mstats.gmean(X, axis=xs)', 'import numpy, scipy.stats.stats; X=numpy.random.gamma(shape=2, scale=1, size=(1,10)); xs=None').timeit(1000) numpy_t=timeit.Timer('numpy.exp((numpy.log(X).mean()))', 'import numpy, numpy.random; X=numpy.random.gamma(shape=2, scale=1, size=(1,10))').timeit(1000) I use Linux and Python 2.5 but my system is very buzy so perhaps not that fair for benchmarks. numpy.__version__ '1.3.0.dev6338' scipy.__version__ '0.8.0.dev5597' There is a cost of using _chk_asarray in this case which decreases as the array size increases. (I am not sure that _chk_asarray is really needed anyhow.) There is a huge cost for using masked array for small sizes but decreases as the array size increases. For 1 by 10 array, the difference between masked and non masked versions was 0.13 seconds to do it 1000 times with the ratio of masked to non masked = 7.94 For 1 by 10000 array, the difference between masked and non masked versions was 0.07 seconds to do it 1000 times with the ratio of masked to non masked = 2.14 However, briefly looking at some of these functions, I think that numpy/scipy would naturally handle the array type as I know numpy.exp((numpy.log(X).mean())) this works whether X is the usual array or if it is a masked array. If so then there is no reason for different functions unless we need to address masks. Bruce
On Fri, Feb 27, 2009 at 12:42 PM, Bruce Southey <bsouthey@gmail.com> wrote:
josef.pktd@gmail.com wrote: [snip]
What I would like to do, but didn't have the time yet is to run the tests for stats.stats on stats.mstats. This way even if we would have some duplicate functions, we would have some cross check that they are consistent, and it would be a reminder for bug fixing also the other version.
Okay, I do not know how to get timeit to work with numpy/scipy but this is not how I would like it to be. But I managed somehow to (unfairly) compare the geometric means function (gmean) using this code: import timeit stand_t=timeit.Timer('scipy.stats.stats.gmean(X, axis=xs)', 'import numpy, scipy.stats.stats; X=numpy.random.gamma(shape=2, scale=1, size=(1,10)); xs=None').timeit(1000) masked_t=timeit.Timer('scipy.stats.mstats.gmean(X, axis=xs)', 'import numpy, scipy.stats.stats; X=numpy.random.gamma(shape=2, scale=1, size=(1,10)); xs=None').timeit(1000) numpy_t=timeit.Timer('numpy.exp((numpy.log(X).mean()))', 'import numpy, numpy.random; X=numpy.random.gamma(shape=2, scale=1, size=(1,10))').timeit(1000)
I use Linux and Python 2.5 but my system is very buzy so perhaps not that fair for benchmarks. numpy.__version__ '1.3.0.dev6338' scipy.__version__ '0.8.0.dev5597'
There is a cost of using _chk_asarray in this case which decreases as the array size increases. (I am not sure that _chk_asarray is really needed anyhow.) There is a huge cost for using masked array for small sizes but decreases as the array size increases.
For 1 by 10 array, the difference between masked and non masked versions was 0.13 seconds to do it 1000 times with the ratio of masked to non masked = 7.94 For 1 by 10000 array, the difference between masked and non masked versions was 0.07 seconds to do it 1000 times with the ratio of masked to non masked = 2.14
However, briefly looking at some of these functions, I think that numpy/scipy would naturally handle the array type as I know numpy.exp((numpy.log(X).mean())) this works whether X is the usual array or if it is a masked array. If so then there is no reason for different functions unless we need to address masks.
Bruce
I just ran the stats.stats test using mstats instead of stats. I didn't look at the results carefully, but the are some numerical inconsistencies between the two implementation, that need to be checked. I attached the test results to http://scipy.org/scipy/scipy/ticket/845. Your timing numbers don't sound so bad in absolute terms, but if it is inside an optimization loop, eg. for maximum likelihood estimation then an 8-fold slowdown can get painful. The main problem for the basic functions, I think, are those functions that need a loop because the data is not rectangular and cannot use simple broad casting and matrix/array operations. On the other hand, I don't think that the masked array functions have been checked for performance ("premature optimization") , since many of them are still relatively new. Josef
All, I followed the thread without actively participating, but as the author of stats.mstats, I feel compelled to jump in. When I started working on some masked versions of scipy.stats version, numpy.ma wasn't part of numpy per se (if I remember correctly), and the package hadn't been thouroughly checked. Modifying scipy.stats to recognize masked arrays (for example, by changing _chk_array) wasn't really an option at the time, because numpy.ma was still considered as experimental. The easiest was therefore just to duplicate the functions. I checked that the results were consistent with the non- masked versions at the time, checked against R also, so I was fairly confident in the results. However, I didn't strive for exhaustivity: I coded the functions I needed and some of them direct relatives, but never tried to expand some more complex functions. I'm all in favor for merging the masked and non-masked versions: that's cleaner, easier to debug and maintain should there be some changes in signature (or even just doc). There's a few aspects we must keep in mind however: * standard numpy functions usually work well with masked arrays: if the input is MA, the np function should call MA.__array_wrap__ which will transform the result back to a MA. I have the nagging feeling it's not completely fool-proof, however, but the numpy.ma functions should always work. if the input is not a MA, then the output will never be masked, which may be a problem. Consider this example:
x=np.array([0,1,2]) np.log(x).mean() -inf ma.log(x).mean() 0.34657359027997264 np.log(ma.array(x)).mean() 0.34657359027997264
If we don't transform x to a MA, or don't use the numpy.ma function, we just get a NaN/Inf as results. Otherwise, we get a nice float. * Systematically using MA may be problematic: I can picture cases where a standard ndarray is expected as output when a standard ndarray is given as inputt. If we force the conversion, the result will be a MA. Should we convert it back to a ndarray ? Using .filled() ? But then, with which filling_value ? In that case, we may want to consider a "usemask" flag: if usemask=True and the input is a ndarray, then the output will be a MA, otherwise it'll be a ndarray. Using a MA as input would set usemask to True no matter what. * The correlation functions discard masked values pair-wise: we can pre-process the inputs and still use the standard functions, so no problem here. * Some functions (threshold) can work directly w/ MA. * Some functions (the ones based on ranking) should behave differently whether the input has masked values (as missing values must be taken as ties). About optimization and speed test: * There's definitely some room for improvement here: for example, instead of using the count method, we could use the count function to prevent any unnecessary conversion to MA. (I'd need to optimize the count function, but that should be easy...). That'll depend on what we decide for handling MA. * Just running tests w/ the masked versions of the function will always show that they are slower, of course. * Slight differences of the order of 1e-15 should not really matter. All, don't hesitate to contact me on or off-list if you have some specific questions about implementation details.
On Fri, Feb 27, 2009 at 1:54 PM, Pierre GM <pgmdevlist@gmail.com> wrote:
All, I followed the thread without actively participating, but as the author of stats.mstats, I feel compelled to jump in.
When I started working on some masked versions of scipy.stats version, numpy.ma wasn't part of numpy per se (if I remember correctly), and the package hadn't been thouroughly checked. Modifying scipy.stats to recognize masked arrays (for example, by changing _chk_array) wasn't really an option at the time, because numpy.ma was still considered as experimental. The easiest was therefore just to duplicate the functions. I checked that the results were consistent with the non- masked versions at the time, checked against R also, so I was fairly confident in the results. However, I didn't strive for exhaustivity: I coded the functions I needed and some of them direct relatives, but never tried to expand some more complex functions.
I'm all in favor for merging the masked and non-masked versions: that's cleaner, easier to debug and maintain should there be some changes in signature (or even just doc). There's a few aspects we must keep in mind however:
* standard numpy functions usually work well with masked arrays: if the input is MA, the np function should call MA.__array_wrap__ which will transform the result back to a MA. I have the nagging feeling it's not completely fool-proof, however, but the numpy.ma functions should always work. if the input is not a MA, then the output will never be masked, which may be a problem.
Consider this example: >>> x=np.array([0,1,2]) >>> np.log(x).mean() -inf >>> ma.log(x).mean() 0.34657359027997264 >>> np.log(ma.array(x)).mean() 0.34657359027997264
If we don't transform x to a MA, or don't use the numpy.ma function, we just get a NaN/Inf as results. Otherwise, we get a nice float.
* Systematically using MA may be problematic: I can picture cases where a standard ndarray is expected as output when a standard ndarray is given as inputt. If we force the conversion, the result will be a MA. Should we convert it back to a ndarray ? Using .filled() ? But then, with which filling_value ?
In that case, we may want to consider a "usemask" flag: if usemask=True and the input is a ndarray, then the output will be a MA, otherwise it'll be a ndarray. Using a MA as input would set usemask to True no matter what.
* The correlation functions discard masked values pair-wise: we can pre-process the inputs and still use the standard functions, so no problem here.
* Some functions (threshold) can work directly w/ MA.
* Some functions (the ones based on ranking) should behave differently whether the input has masked values (as missing values must be taken as ties).
About optimization and speed test: * There's definitely some room for improvement here: for example, instead of using the count method, we could use the count function to prevent any unnecessary conversion to MA. (I'd need to optimize the count function, but that should be easy...). That'll depend on what we decide for handling MA. * Just running tests w/ the masked versions of the function will always show that they are slower, of course. * Slight differences of the order of 1e-15 should not really matter.
All, don't hesitate to contact me on or off-list if you have some specific questions about implementation details.
I still need to look at several examples, before I get a better feeling of how this will work. I just looked at the implementation of ma.var, ma.cov, ma.exp and a few more, and I think they are very well written and I don't see any way how their performance could be improved. Given that gmean is a very simple function, I was pretty surprised about the difference in timing. Now, I think that the main slowdown is that the mask has to be checked in every operation that calls a ma.* version of a function. As we discussed for the OLS case for larger statistical functions, building the main workload with plain arrays will save a lot of overhead. This works for cases where a single compression or fill is correct for all required numerical operations. If we get the correct setup (interface, conversion) for two kinds of functions, any array to ma core of the function", and "any array to plain core", then it will be easier, at least for me, to follow this pattern when (re)writing functions. One more issue is the treatment of nan and masked values, for example, if a function produces nans because of a zero division, then I would want to treat it differently than a missing value in the data. If it is automatically included in the mask then this distinction is lost. Or is there a different use case for this? In your log example, I wouldn't want to get a nice number back. I want the function to complain. Silently changing the definition of mathematical operations creates a huge potential for errors (that's why I also don't like the silent conversions when casting to int) For example, if this is maximum likelihood estimation, the log likelihood is -inf and not some nice number.
x=np.array([0,1,2]) np.log(x).mean() I think if users want nice numbers, then they should mask them in the first place. Actually, I didn't realize this before, that ma adds additional points to the mask.
But, before we start to rewrite and refactor across the board, I still want to finish cleaning up the existing functions and resolve some of the current inconsistencies. Josef
Given that gmean is a very simple function, I was pretty surprised about the difference in timing. Now, I think that the main slowdown is that the mask has to be checked in every operation that calls a ma.* version of a function.
It's actually a tad more complex: ma.log checks the mask of the input, but also converts the output to a MA when needed, with all the overhead of MA.__array_finalize__.
As we discussed for the OLS case for larger statistical functions, building the main workload with plain arrays will save a lot of overhead. This works for cases where a single compression or fill is correct for all required numerical operations.
That's indeed the way to go: preprocess a MA to transform it into a ndarray (by dropping masked values, or processing them afterwards), perform the operation, revert to MA if needed.
One more issue is the treatment of nan and masked values, for example, if a function produces nans because of a zero division, then I would want to treat it differently than a missing value in the data. If it is automatically included in the mask then this distinction is lost. Or is there a different use case for this?
Nope. If a value get masked by an operation, you won't be able to track it (unless by comparing the mask of the output w/ the mask of the input).
In your log example, I wouldn't want to get a nice number back. I want the function to complain.
Because you work w/ ndarrays. If I work w/ MA, I expect it not to crash but drop the masked values.
But, before we start to rewrite and refactor across the board, I still want to finish cleaning up the existing functions and resolve some of the current inconsistencies.
Well, you may double the workload. One way would be to first agree on how we should refactor/reorganize the functions, then clean the ndarray part of the function. We can always add a NotImplementedError if the input is a MA w/ missing values.
In your log example, I wouldn't want to get a nice number back. I want the function to complain. Silently changing the definition of mathematical operations creates a huge potential for errors (that's why I also don't like the silent conversions when casting to int) For example, if this is maximum likelihood estimation, the log likelihood is -inf and not some nice number.
x=np.array([0,1,2]) np.log(x).mean() I think if users want nice numbers, then they should mask them in the first place.
the more I think, about
np.ma.log([0,1,2]).sum() 0.69314718055994529 np.log([0,1,2]).sum() -inf
the more worried, I get about using ma functions. One example: In the fit method of the distributions with bounded support, if there are observations outside of the bound than the negative log-likelihood is set to inf: cond0 = (x <= self.a) | (x >= self.b) if (any(cond0)): return inf else: N = len(x) return self._nnlf(x, *args) + N*log(scale) In this case, it might still produce the correct result since the check is before the aggregation. However, this is implementation specific. If I had assigned the inf before the summation of the log-likelihood contributions, ma.log would have removed them, and killed the boundary check. So when working with masked array functions, it is necessary to always keep in mind that the math is defined differently, which promises many happy hours of bug hunting. Josef
On Feb 27, 2009, at 3:14 PM, josef.pktd@gmail.com wrote:
One example: In the fit method of the distributions with bounded support, if there are observations outside of the bound than the negative log-likelihood is set to inf:
cond0 = (x <= self.a) | (x >= self.b) if (any(cond0)): return inf else: N = len(x) return self._nnlf(x, *args) + N*log(scale)
In this case, it might still produce the correct result since the check is before the aggregation. However, this is implementation specific. If I had assigned the inf before the summation of the log-likelihood contributions, ma.log would have removed them, and killed the boundary check.
x = ma.array([0,1,2],mask=[0,1,0]) np.log(x) masked_array(data = [-- -- 0.69314718056], mask = [ True True False], fill_value = 1e+20) np.log(x) first work on the data, then call MA.__array_wrap__. This function checks the initial mask, then the context of the function: as it's a domained function, the entries outside the domain are
OK, so you don't want to use the ma functions there. Pb is that you won't be able to use the np versions on MA either transformed into mask. For this kind of problem, the easiest is to decouple: 1. Take a view of the input as a standard ndarray. 2. Process the view 3. Add the mask of the input if needed. With the previous example, that'd be roughly
ma.array(np.log(x.view(ndarray)), mask=ma.getmask(x)) masked_array(data = [-inf -- 0.69314718056], mask = [False True False], fill_value = 1e+20)
You keep the masked entry at index 1, but don't mask the entry at index 0.
So when working with masked array functions, it is necessary to always keep in mind that the math is defined differently, which promises many happy hours of bug hunting.
Indeed. But once again, the masked versions of the function are more for convenience. If you need performance, you have to preprocess the inputs by transforming them into standard ndarrays one way or another. In the case of correlation functions, for example, you can suppress missing values pair-wise (that is, drop the entries of x if the corresponding entries of y are masked, and vice-versa). For basic linear fit, that might be an approach. A second would be to work by intervals, the limits of the intervals being a masked value. For more complex fitting (eg, loess), problems arise. You can bypass them temporarily by raisong a NotImplementedError if the inputs are masked, it'd be up to the user to find a way to fill the inputs.
Indeed. But once again, the masked versions of the function are more for convenience. If you need performance, you have to preprocess the inputs by transforming them into standard ndarrays one way or another. In the case of correlation functions, for example, you can suppress missing values pair-wise (that is, drop the entries of x if the corresponding entries of y are masked, and vice-versa). For basic linear fit, that might be an approach. A second would be to work by intervals, the limits of the intervals being a masked value. For more complex fitting (eg, loess), problems arise. You can bypass them temporarily by raisong a NotImplementedError if the inputs are masked, it'd be up to the user to find a way to fill the inputs.
For most of the current statistical functions, with the exception of different tie handling, I think that we can expand the _chk_asarray to do the necessary preprocessing. I also thought that the return for these functions should be easy, since most of them return statistics and not data arrays. However, looking at some examples, it is not obvious to me which return result and type you would like to have. A good example is `moment`, here is some of the current returns. I think they cover the main return patterns.
x array([ 0., 1., NaN, 2.])
masked arrays with masked nan: ------------------------------------------------ do you need a masked array as return type, since all values are valid? how about for t-statistic and p-values? Do p-values need to be masked arrays?
stats.mstats.moment(np.ma.fix_invalid(np.ma.column_stack([x,x])),3) masked_array(data = [0.0 0.0], mask = [False False], fill_value = 1e+020)
stats.mstats.moment(np.ma.fix_invalid(np.ma.column_stack([x,x])),2) masked_array(data = [0.666666666667 0.666666666667], mask = [False False], fill_value = 1e+020)
stats.mstats.moment(np.ma.fix_invalid(np.ma.column_stack([x,x])),1) #inconsistent return type array([ 0., 0.])
masked array without masked values ----------------------------------------------------- same as above about return type
stats.mstats.moment(np.ma.column_stack([np.arange(4),np.arange(4)]),2) masked_array(data = [ 1.25 1.25], mask = False, fill_value = 1e+020)
masked array with nan that is not masked ------------------------------------------------------------- masked array in, masked array out, nan results converted to mask is this desired?
stats.mstats.moment(np.ma.column_stack([x,x]),3) masked_array(data = [-- --], mask = [ True True], fill_value = 1e+020)
stats.mstats.moment(np.ma.column_stack([x,x]),0) masked_array(data = [-- --], mask = [ True True], fill_value = 1e+020)
stats.mstats.moment(np.ma.column_stack([x,x]),1) array([ 0., 0.])
ndarray with nans ------------------------- converted to masked array, nans are masked. here I want to get ndarray with nans returned
stats.mstats.moment(np.column_stack([x,x]),0) masked_array(data = [-- --], mask = [ True True], fill_value = 1e+020)
stats.mstats.moment(np.column_stack([x,x]),1) array([ 0., 0.])
ndarray without nans ------------------------------ this should return ndarray
stats.mstats.moment(np.column_stack([np.arange(4),np.arange(4)]),2) masked_array(data = [ 1.25 1.25], mask = False, fill_value = 1e+020)
stats.mstats.moment(np.column_stack([np.arange(4),np.arange(4)]),1) array([ 0., 0.])
If this return "API" is specified, then it is possible to work out some examples to see how the merged function works. If you think converting nans that are the result of calculations to masked arrays are important, then we could add a keyword argument that implies a fix_invalid before returning the results. Josef
On Feb 27, 2009, at 4:52 PM, josef.pktd@gmail.com wrote:
For most of the current statistical functions, with the exception of different tie handling, I think that we can expand the _chk_asarray to do the necessary preprocessing.
Mmh. _chk_asarray will always return a MA. Is it what you want? Are you An idea is then to use the 'usemask' parameter I was talking about earlier: * if usemask is False (default), return a ndarray * If usemask is True, return a MA * if the input is a MA (w/ or w/o missing values), set usemask to True, and mask the NaNs/Infs first w/ ma.fix_invalid. That way, we need only one function. If we really need it, we can have duplicate functions in scipy.mstats where usemask is set to True by default. Now, for the actual implementation: * usemask=False and some NaNs: return NaN * usemask=True: use the ma implementation.
stats.mstats.moment(np.ma.fix_invalid(np.ma.column_stack([x,x])), 1) #inconsistent return type array([ 0., 0.])
That's a bug, we should have a MA.
Pierre GM wrote:
On Feb 27, 2009, at 4:52 PM, josef.pktd@gmail.com wrote:
For most of the current statistical functions, with the exception of different tie handling, I think that we can expand the _chk_asarray to do the necessary preprocessing.
Mmh. _chk_asarray will always return a MA. Is it what you want? Are you
An idea is then to use the 'usemask' parameter I was talking about earlier: * if usemask is False (default), return a ndarray * If usemask is True, return a MA * if the input is a MA (w/ or w/o missing values), set usemask to True, and mask the NaNs/Infs first w/ ma.fix_invalid.
This may not be appropriate for scipy, but for my own purposes I included a third option for the similar "masked" kwarg in a simple stats class: http://currents.soest.hawaii.edu/hg/hgwebdir.cgi/pycurrents/file/7b4103d34cc... masked='auto' makes the output masked if and only if the input is masked. Eric
That way, we need only one function. If we really need it, we can have duplicate functions in scipy.mstats where usemask is set to True by default.
Now, for the actual implementation: * usemask=False and some NaNs: return NaN * usemask=True: use the ma implementation.
stats.mstats.moment(np.ma.fix_invalid(np.ma.column_stack([x,x])), 1) #inconsistent return type array([ 0., 0.])
That's a bug, we should have a MA.
_______________________________________________ Scipy-dev mailing list Scipy-dev@scipy.org http://projects.scipy.org/mailman/listinfo/scipy-dev
On Fri, Feb 27, 2009 at 6:01 PM, Eric Firing <efiring@hawaii.edu> wrote:
Pierre GM wrote:
On Feb 27, 2009, at 4:52 PM, josef.pktd@gmail.com wrote:
For most of the current statistical functions, with the exception of different tie handling, I think that we can expand the _chk_asarray to do the necessary preprocessing.
Mmh. _chk_asarray will always return a MA. Is it what you want? Are you
An idea is then to use the 'usemask' parameter I was talking about earlier: * if usemask is False (default), return a ndarray * If usemask is True, return a MA * if the input is a MA (w/ or w/o missing values), set usemask to True, and mask the NaNs/Infs first w/ ma.fix_invalid.
This may not be appropriate for scipy, but for my own purposes I included a third option for the similar "masked" kwarg in a simple stats class:
http://currents.soest.hawaii.edu/hg/hgwebdir.cgi/pycurrents/file/7b4103d34cc...
masked='auto' makes the output masked if and only if the input is masked.
Eric
Yes, your class looks similar to what I have in mind. But I didn't see a license statement to know whether I'm allowed to look. Also, your broadcastable (squeeze) option looks like a very useful idea. Two differences that I think of are to have the main part in ndarrays while your _y is a masked array, and at this stage we won't switch to classes for the basic statistical functions. Additionally, if I rewrite these functions I would like to get also weights in. Josef
On Fri, Feb 27, 2009 at 5:47 PM, Pierre GM <pgmdevlist@gmail.com> wrote:
On Feb 27, 2009, at 4:52 PM, josef.pktd@gmail.com wrote:
For most of the current statistical functions, with the exception of different tie handling, I think that we can expand the _chk_asarray to do the necessary preprocessing.
Mmh. _chk_asarray will always return a MA. Is it what you want? Are you
No, what I meant was, that _chk_asarray is currently called for preprocessing in most functions, so it will be easy to use a replacement function to obtain the preprocessed (e.g. compressed) data, and whatever flags (usemask) we need, in the main body of the function and for the decision about the return type.
An idea is then to use the 'usemask' parameter I was talking about earlier: * if usemask is False (default), return a ndarray * If usemask is True, return a MA * if the input is a MA (w/ or w/o missing values), set usemask to True, and mask the NaNs/Infs first w/ ma.fix_invalid.
That way, we need only one function. If we really need it, we can have duplicate functions in scipy.mstats where usemask is set to True by default.
Now, for the actual implementation: * usemask=False and some NaNs: return NaN * usemask=True: use the ma implementation.
That clarifies the API. I will try to write a prototype, but I spend too much time on scipy this week.
stats.mstats.moment(np.ma.fix_invalid(np.ma.column_stack([x,x])), 1) #inconsistent return type array([ 0., 0.])
That's a bug, we should have a MA.
On Fri, Feb 27, 2009 at 5:13 PM, <josef.pktd@gmail.com> wrote:
On Fri, Feb 27, 2009 at 5:47 PM, Pierre GM <pgmdevlist@gmail.com> wrote:
On Feb 27, 2009, at 4:52 PM, josef.pktd@gmail.com wrote:
For most of the current statistical functions, with the exception of different tie handling, I think that we can expand the _chk_asarray to do the necessary preprocessing.
Mmh. _chk_asarray will always return a MA. Is it what you want? Are you
No, what I meant was, that _chk_asarray is currently called for preprocessing in most functions, so it will be easy to use a replacement function to obtain the preprocessed (e.g. compressed) data, and whatever flags (usemask) we need, in the main body of the function and for the decision about the return type.
I really do not see the requirement for _chk_asarray at all. When a user passes a typical array or masked array then there should be no further processing required. Also _chk_asarray will use ravel() if axis is None but my understanding of many numpy functions operate over a flattened array when there is no axis defined. The only case that needs addressing is when a user supplies an object that can be converted to an array otherwise a error needs to be raised. After conversion to an array no further processing is required and even that conversion in some cases will be done within the existing functions.
An idea is then to use the 'usemask' parameter I was talking about earlier: * if usemask is False (default), return a ndarray * If usemask is True, return a MA * if the input is a MA (w/ or w/o missing values), set usemask to True, and mask the NaNs/Infs first w/ ma.fix_invalid.
That way, we need only one function. If we really need it, we can have duplicate functions in scipy.mstats where usemask is set to True by default.
Now, for the actual implementation: * usemask=False and some NaNs: return NaN * usemask=True: use the ma implementation.
That clarifies the API. I will try to write a prototype, but I spend too much time on scipy this week.
This is a little messy and there has been discussion regarding this elsewhere. In these terms there are two distinct issues: 1) If the array contains non-finite numbers (NaN, positive and negative infinity) then perhaps the user can strip these out first for example R's mean function has the argument 'na.rm = FALSE'. 2) If non-finite elements arise during the function like taking the log of zero then I think that the user must know these have occurred rather than be forced to check the mask - especially if they have already masked values for other reasons like incomplete data. Bruce
On Fri, Feb 27, 2009 at 10:04 PM, Bruce Southey <bsouthey@gmail.com> wrote:
On Fri, Feb 27, 2009 at 5:13 PM, <josef.pktd@gmail.com> wrote:
On Fri, Feb 27, 2009 at 5:47 PM, Pierre GM <pgmdevlist@gmail.com> wrote:
On Feb 27, 2009, at 4:52 PM, josef.pktd@gmail.com wrote:
For most of the current statistical functions, with the exception of different tie handling, I think that we can expand the _chk_asarray to do the necessary preprocessing.
Mmh. _chk_asarray will always return a MA. Is it what you want? Are you
No, what I meant was, that _chk_asarray is currently called for preprocessing in most functions, so it will be easy to use a replacement function to obtain the preprocessed (e.g. compressed) data, and whatever flags (usemask) we need, in the main body of the function and for the decision about the return type.
I really do not see the requirement for _chk_asarray at all. When a user passes a typical array or masked array then there should be no further processing required. Also _chk_asarray will use ravel() if axis is None but my understanding of many numpy functions operate over a flattened array when there is no axis defined.
The only case that needs addressing is when a user supplies an object that can be converted to an array otherwise a error needs to be raised. After conversion to an array no further processing is required and even that conversion in some cases will be done within the existing functions.
The current usage allows to pass lists instead of arrays. This is very convenient for interactive use but might also have other uses, e.g when building a list incrementally. And I thought asarray doesn't have much cost if it is already an array. I didn't look systematically at ravel, but while axis=None works automatically for many numpy functions, for more complex statistical functions more control over the dimension of the input arrays is necessary. Many statistical functions are only designed for 1d or 2d and controlling the dimension at the beginning simplifies the main part of the functions. I had some cases where I was struggling for a while with the dimensions and axis, but in many cases it could be redundant. If we want to handle different array types with the same function then the _chk_asarray call will be replaced by the type specific preprocessing.
An idea is then to use the 'usemask' parameter I was talking about earlier: * if usemask is False (default), return a ndarray * If usemask is True, return a MA * if the input is a MA (w/ or w/o missing values), set usemask to True, and mask the NaNs/Infs first w/ ma.fix_invalid.
That way, we need only one function. If we really need it, we can have duplicate functions in scipy.mstats where usemask is set to True by default.
Now, for the actual implementation: * usemask=False and some NaNs: return NaN * usemask=True: use the ma implementation.
That clarifies the API. I will try to write a prototype, but I spend too much time on scipy this week.
This is a little messy and there has been discussion regarding this elsewhere. In these terms there are two distinct issues: 1) If the array contains non-finite numbers (NaN, positive and negative infinity) then perhaps the user can strip these out first for example R's mean function has the argument 'na.rm = FALSE'.
If the merged functions are able to handle masked arrays and plain ndarrays, then we can also offer the user the option for the treatment of nans, this would make the separate nanmean, ... obsolete. Operation on inf might be too ambiguous and I would think they are the responsibility of the user. And if there is a inf*0 then I want to give them the nan back, and the user can decide what to do. In general inf is a legitimate number and might or should propagate correctly (if the user wants to leave them in) e.g.
stats.norm.cdf(-np.inf) 0.0 stats.norm.cdf(np.inf) 1.0
2) If non-finite elements arise during the function like taking the log of zero then I think that the user must know these have occurred rather than be forced to check the mask - especially if they have already masked values for other reasons like incomplete data.
I agree and I want this behavior for ndarrays, for masked arrays I'm less involved since I'm not using them (yet). I like Erics use of a trivariate (?) choice with "auto" which adds one option for the user: masked='auto' : True|False|'auto' determines the output; if True, output will be a masked array; if False, output will be an ndarray with nan used as a bad flag if necessary; if 'auto', output will match input What the exact definition is for masked arrays in the case "auto", is up to the masked array users. Josef
Hi, I am seeing a few functions that should be made depreciated as these appear to duplicate Numpy or Scipy functions. Do you want these as new or old tickets (for example, samplestd has ticket #81 as part of the Statistics Review)? Would you want a large patch or one for each ticket? These functions are just renamed functions present in scipy.special just with perhaps slightly more informative names: erfc ksprob fprob chisqprob zprob But I do not think we need these as separate functions but there is the issue of depreciation involved if users use these specific functions. There are other like that should be treated as depreciated: samplestd samplevar
import numpy as np import scipy.stats.stats as stats a=np.array([[1,2,3,4,5], [6,7,8,9,10]]) np.std(a,axis=0) array([ 2.5, 2.5, 2.5, 2.5, 2.5]) stats.samplestd(a,axis=0) array([ 2.5, 2.5, 2.5, 2.5, 2.5]) stats.samplestd(a,axis=None) 2.8722813232690143 np.std(a,axis=None) 2.8722813232690143
Also, stats.py has the histogram and histogram2 functions where I agree with the comment in the code about being obsoleted by numpy.histogram. I would think these should be depreciated although the cumfreq and relfreq functions would need to be rewritten, Thanks Bruce
On Mon, Mar 2, 2009 at 2:09 PM, Bruce Southey <bsouthey@gmail.com> wrote:
Hi, I am seeing a few functions that should be made depreciated as these appear to duplicate Numpy or Scipy functions.
Do you want these as new or old tickets (for example, samplestd has ticket #81 as part of the Statistics Review)? Would you want a large patch or one for each ticket?
I agree with all the depreciation, and there might be some more (eg. sem and stderr are essentially the same). For depreciation warnings I would prefer one new ticket with one patch (or easier for me to verify is the changed complete sourcefile of stats.py)
These functions are just renamed functions present in scipy.special just with perhaps slightly more informative names: erfc ksprob fprob chisqprob zprob
Most calls to these functions can be replaced to calls to the distribution, e.g distributions.f.sf, as I did for the t-tests. However, I have seen them used in some external packages, and a release with a depreciation warning might be necessary.
But I do not think we need these as separate functions but there is the issue of depreciation involved if users use these specific functions.
There are other like that should be treated as depreciated: samplestd samplevar
>>> import numpy as np >>> import scipy.stats.stats as stats >>> a=np.array([[1,2,3,4,5], [6,7,8,9,10]]) >>> np.std(a,axis=0) array([ 2.5, 2.5, 2.5, 2.5, 2.5]) >>> stats.samplestd(a,axis=0) array([ 2.5, 2.5, 2.5, 2.5, 2.5]) >>> stats.samplestd(a,axis=None) 2.8722813232690143 >>> np.std(a,axis=None) 2.8722813232690143
Also, stats.py has the histogram and histogram2 functions where I agree with the comment in the code about being obsoleted by numpy.histogram. I would think these should be depreciated although the cumfreq and relfreq functions would need to be rewritten,
I never looked closely at the histogram and histogram2 functions in stats, because I also use the numpy version. So I don't know if they have equivalent functionality. Neither the histogram functions nor cumfreq and relfreq have tests, so before depreciating we should find out what these functions are doing for different cases.
Thanks Bruce
Thank you for checking this Josef
josef.pktd@gmail.com wrote:
On Mon, Mar 2, 2009 at 2:09 PM, Bruce Southey <bsouthey@gmail.com> wrote:
Hi, I am seeing a few functions that should be made depreciated as these appear to duplicate Numpy or Scipy functions.
Do you want these as new or old tickets (for example, samplestd has ticket #81 as part of the Statistics Review)? Would you want a large patch or one for each ticket?
I agree with all the depreciation, and there might be some more (eg. sem and stderr are essentially the same). For depreciation warnings I would prefer one new ticket with one patch (or easier for me to verify is the changed complete sourcefile of stats.py)
These functions are just renamed functions present in scipy.special just with perhaps slightly more informative names: erfc ksprob fprob chisqprob zprob
Most calls to these functions can be replaced to calls to the distribution, e.g distributions.f.sf, as I did for the t-tests. However, I have seen them used in some external packages, and a release with a depreciation warning might be necessary.
I agree that these should be first depreciated. I will try to write these when I get the time.
But I do not think we need these as separate functions but there is the issue of depreciation involved if users use these specific functions.
There are other like that should be treated as depreciated: samplestd samplevar
Okay, I have created two tickets with hopefully suitable patches for: samplevar: 877 samplestd: 878
I did not change the info.py and any tests but these will need to be changed if the patches are applied. Also, if you apply these patches, I think that tickets 80 and 81 can be closed.
Also, stats.py has the histogram and histogram2 functions where I agree with the comment in the code about being obsoleted by numpy.histogram. I would think these should be depreciated although the cumfreq and relfreq functions would need to be rewritten,
I never looked closely at the histogram and histogram2 functions in stats, because I also use the numpy version. So I don't know if they have equivalent functionality.
I have not examined it in detail, histogram is different from numpy in various ways like arguments and implementation. But, after the histogram discussion on the numpy list, I do not consider that it is sufficiently different than the numpy version to justify yet another version. I think histogram2 is just a utility function than a useful function and it not used elsewhere.
Neither the histogram functions nor cumfreq and relfreq have tests, so before depreciating we should find out what these functions are doing for different cases.
When I get to these functions, I will look into providing tests. Also these may need changes depending on what happens with histogram.
Thanks Bruce
Thank you for checking this
Josef
No problems (yet) especially when I will have more 'issues' as I go through these functions. Bruce
participants (7)
-
Bruce Southey -
Eric Firing -
Gael Varoquaux -
josef.pktd＠gmail.com -
Pierre GM -
Robert Kern -
Stéfan van der Walt