BUG: stats: wrong shape from median_absolute_deviation for arrays with size 0 and an axis.

[WarrenWeckesser]

Currently, median_absolute_deviation returns the scalar nan if the size of the input array is 0, even if an axis argument is given. This is not correct. The function must respect the reduction-like behavior for the shape of the arrays. For example, if x has shape (m, n, p) and axis=1, the shape of the result must be (m, p), even if any of m, n, or p are 0.

For example,

In [21]: z = np.ones((3, 0))                                                                                            

Here, just like n.sum(z, axis=0), the result should be an array with shape (0,), but in the master branch, we incorrectly get the scalar nan:

In [22]: median_absolute_deviation(z, axis=0)  # from the master branch, not correct!
Out[22]: nan

With axis=1, the result should be an array with shape (3,). The values in the array should be the same as the value produced by median_absolute_deviation([]), which is nan, so the actual result should be array([nan, nan, nan]). But in the master branch, the behavior is again incorrect--instead of getting the array with shape (3,), we get the scalar nan:

In [23]: median_absolute_deviation(z, axis=1)  # from the master branch, not correct!
Out[23]: nan

The above code was function with the current master branch:

In [28]: import scipy

In [29]: scipy.__version__
Out[29]: '1.5.0.dev0+69f44c0'

[WarrenWeckesser]

gh-11735 wasn't a complete fix, so I'm reopening the issue.

[WarrenWeckesser]

This was fixed in #11431. median_absolute_deviation is now deprecated, but it does give the correct results for arrays with size 0.