Repository: scipy/scipy
File: scipy/stats/_stats_py.py (Python)

Code before:
    `tsem` uses unbiased sample standard deviation, i.e. it uses a
    correction factor ``n / (n - 1)``.

    Examples
    --------
    >>> import numpy as np
    >>> from scipy import stats
    >>> x = np.arange(20)
    >>> stats.tsem(x)
    1.3228756555322954
    >>> stats.tsem(x, (3,17))
    1.1547005383792515

    """
    xp = array_namespace(a)
    a, _ = _put_val_to_limits(a, limits, inclusive, xp=xp)

    with warnings.catch_warnings():
        warnings.simplefilter("ignore", SmallSampleWarning)
        # Currently, this behaves like nan_policy='omit' for alternative array
        # backends, but nan_policy='propagate' will be handled for other backends
        # by the axis_nan_policy decorator shortly.
        sd = _xp_var(a, correction=ddof, axis=axis, nan_policy='omit', xp=xp)**0.5

    not_nan = xp.astype(~xp.isnan(a), a.dtype)
    n_obs = xp.sum(not_nan, axis=axis, dtype=sd.dtype)
    return sd / n_obs**0.5


#####################################
#              MOMENTS              #
#####################################


def _moment_outputs(kwds, default_order=1):
    order = np.atleast_1d(kwds.get('order', default_order))
    message = "`order` must be a scalar or a non-empty 1D array."
    if order.size == 0 or order.ndim > 1:
        raise ValueError(message)
    return len(order)


def _moment_result_object(*args):
    if len(args) == 1:
        return args[0]
    xp = array_namespace(*args)
    return xp.stack(args)


# When `order` is array-like with size > 1, moment produces an *array*
# rather than a tuple, but the zeroth dimension is to be treated like

Review comment (suggestion):
See code comment. For better or for worse, this is supposed to iterate over the zeroth axis, and we have to do that manually for `array_api_strict`.

Code after:
#####################################
#              MOMENTS              #
#####################################


def _moment_outputs(kwds, default_order=1):
    order = np.atleast_1d(kwds.get('order', default_order))
    message = "`order` must be a scalar or a non-empty 1D array."
    if order.size == 0 or order.ndim > 1:
        raise ValueError(message)
    return len(order)


def _moment_result_object(*args):
    if len(args) == 1:
        return args[0]
    xp = array_namespace(*args)
    return xp.stack(args)


# When `order` is array-like with size > 1, moment produces an *array*
# rather than a tuple, but the zeroth dimension is to be treated like
# separate outputs. It is important to make the distinction between
# separate outputs when adding the reduced axes back (`keepdims=True`).
def _moment_tuple(x, n_out):
    return tuple(x[i, ...] for i in range(x.shape[0])) if n_out > 1 else (x,)


# `moment` fits into the `_axis_nan_policy` pattern, but it is a bit unusual
# because the number of outputs is variable. Specifically,
# `result_to_tuple=lambda x: (x,)` may be surprising for a function that
# can produce more than one output, but it is intended here.
# When `moment is called to produce the output:
# - `result_to_tuple` packs the returned array into a single-element tuple,
# - `_moment_result_object` extracts and returns that single element.
# However, when the input array is empty, `moment` is never called. Instead,
# - `_check_empty_inputs` is used to produce an empty array with the
#   appropriate dimensions.
# - A list comprehension creates the appropriate number of copies of this
#   array, depending on `n_outputs`.
# - This list - which may have multiple elements - is passed into
#   `_moment_result_object`.
# - If there is a single output, `_moment_result_object` extracts and returns
#   the single output from the list.
# - If there are multiple outputs, and therefore multiple elements in the list,
#   `_moment_result_object` converts the list of arrays to a single array and
#   returns it.
# Currently, this leads to a slight inconsistency: when the input array is
# empty, there is no distinction between the `moment` function being called
# with parameter `order=1` and `order=[1]`; the latter *should* produce
# the same as the former but with a singleton zeroth dimension.