Re: Fixing gh-17749 requires a small API change to `stats.moment`

Hi team,

To summarize, the issue (https://github.com/scipy/scipy/issues/17749)
requested that `stats.moment` be enhanced to allow the computation of
non-central sample moments.

The proposed solution is to add an optional keyword argument `center` (or
`mean` in the original email).

   - For central moments (current behavior), use the default value
   (`None`). The sample mean will be calculated from the provided sample (as
   it is now).
   - For efficiency, pass the sample mean (if it is already available).
   - For non-central moments, pass `0`.

Since this requires an addition to the API, we'd like to know whether the
proposed enhancement sounds good and, if so, what the new parameter should
be called.

Matt

On Sat, Mar 11, 2023 at 1:54 PM Doron Behar <[email protected]> wrote:

> Hey all,
>
> I suggested a rather simple fix for [gh-17749][1], and it's a sort of
> non-breaking API change. I was replied in a comment that such a change
> should go through this mailing list. Here's the Gist of the fix, without
> unit tests and without docs updates. I will appreciate comments/reviews:
>
> [1]: https://github.com/scipy/scipy/issues/17749
>
> ```patch
> commit f16d2282bd34726960d511c1ed277ba2754d285e
> Author: Doron Behar <[email protected]>
> Date:   Sat Feb 4 14:50:02 2023 +0200
>
>     stats: Allow computing non central moments
>
> diff --git a/result b/result
> new file mode 120000
> index 000000000..5cff41b90
> --- /dev/null
> +++ b/result
> @@ -0,0 +1 @@
>
> +/nix/store/38g1b2ai2k35kq7jw0ffm5la90zaqa0h-python3.10-scipy-1.11.0-dirty-flake
> \ No newline at end of file
> diff --git a/scipy/stats/_mstats_basic.py b/scipy/stats/_mstats_basic.py
> index 2d10561d2..7676c7798 100644
> --- a/scipy/stats/_mstats_basic.py
> +++ b/scipy/stats/_mstats_basic.py
> @@ -2551,7 +2551,7 @@ def winsorize(a, limits=None, inclusive=(True,
> True), inplace=False,
>                                     upinc, contains_nan, nan_policy)
>
>
> -def moment(a, moment=1, axis=0):
> +def moment(a, moment=1, mean=None, axis=0):
>      """
>      Calculates the nth moment about the mean for a sample.
>
> @@ -2561,13 +2561,15 @@ def moment(a, moment=1, axis=0):
>         data
>      moment : int, optional
>         order of central moment that is returned
> +    mean : float, optional
> +       If True, then the moments are computed centralized. Default is
> True.
>      axis : int or None, optional
>         Axis along which the central moment is computed. Default is 0.
>         If None, compute over the whole array `a`.
>
>      Returns
>      -------
> -    n-th central moment : ndarray or float
> +    n-th (central) moment : ndarray or float
>         The appropriate moment along the given axis or over all values if
> axis
>         is None. The denominator for the moment calculation is the number
> of
>         observations, no degrees of freedom correction is done.
> @@ -2592,8 +2594,11 @@ def moment(a, moment=1, axis=0):
>
>      # for array_like moment input, return a value for each.
>      if not np.isscalar(moment):
> -        mean = a.mean(axis, keepdims=True)
> -        mmnt = [_moment(a, i, axis, mean=mean) for i in moment]
> +        if mean is None:
> +            _mean = a.mean(axis, keepdims=True)
> +        else:
> +            _mean = mean
> +        mmnt = [_moment(a, i, axis, mean=_mean) for i in moment]
>          return ma.array(mmnt)
>      else:
>          return _moment(a, moment, axis)
> diff --git a/scipy/stats/_stats_py.py b/scipy/stats/_stats_py.py
> index 32445bedf..70d098baa 100644
> --- a/scipy/stats/_stats_py.py
> +++ b/scipy/stats/_stats_py.py
> @@ -1046,14 +1046,14 @@ def moment(a, moment=1, axis=0,
> nan_policy='propagate'):
>
>      Notes
>      -----
> -    The k-th central moment of a data sample is:
> +    The k-th moment of a data sample is:
>
>      .. math::
>
> -        m_k = \frac{1}{n} \sum_{i = 1}^n (x_i - \bar{x})^k
> +        m_k = \frac{1}{n} \sum_{i = 1}^n (x_i)^k
>
> -    Where n is the number of samples and x-bar is the mean. This function
> uses
> -    exponentiation by squares [1]_ for efficiency.
> +    Where n is the number of samples. This function uses exponentiation by
> +    squares [1]_ for efficiency.
>
>      Note that, if `a` is an empty array (``a.size == 0``), array `moment`
> with
>      one element (`moment.size == 1`) is treated the same as scalar
> `moment`
> ```
> _______________________________________________
> SciPy-Dev mailing list -- [email protected]
> To unsubscribe send an email to [email protected]
> https://mail.python.org/mailman3/lists/scipy-dev.python.org/
> Member address: [email protected]
>


-- 
Matt Haberland
Assistant Professor
BioResource and Agricultural Engineering
08A-3K, Cal Poly

_______________________________________________
SciPy-Dev mailing list -- [email protected]
To unsubscribe send an email to [email protected]
https://mail.python.org/mailman3/lists/scipy-dev.python.org/
Member address: [email protected]