summaryrefslogtreecommitdiff
path: root/numpy/lib
diff options
context:
space:
mode:
authorStephan Hoyer <shoyer@climate.com>2015-01-03 18:37:33 +0200
committerStephan Hoyer <shoyer@climate.com>2015-01-05 11:08:16 -0800
commit1b2d2e977d51414121548d4454afdf02d9def3b3 (patch)
tree4cdf43e2cbe8f2ab6df2f95183dfeb23f8e08de3 /numpy/lib
parent3ef77eea0d9c2cd76bc9b89b04a32f1322f842d5 (diff)
downloadnumpy-1b2d2e977d51414121548d4454afdf02d9def3b3.tar.gz
ENH: add np.nanprod
This PR adds an implementation of `nanprod`. The actual function is a two-liner adapted from `nansum`. Most of this PR consists of documentation and tests (for which I took the opportunity to do some consolidation). A method with the same functionality exists in pandas, and I was surprised to discover that it's not in numpy.
Diffstat (limited to 'numpy/lib')
-rw-r--r--numpy/lib/nanfunctions.py77
-rw-r--r--numpy/lib/tests/test_nanfunctions.py229
2 files changed, 157 insertions, 149 deletions
diff --git a/numpy/lib/nanfunctions.py b/numpy/lib/nanfunctions.py
index 7260a35b8..1e3208ac8 100644
--- a/numpy/lib/nanfunctions.py
+++ b/numpy/lib/nanfunctions.py
@@ -9,9 +9,12 @@ Functions
- `nanargmin` -- index of minimum non-NaN value
- `nanargmax` -- index of maximum non-NaN value
- `nansum` -- sum of non-NaN values
+- `nanprod` -- product of non-NaN values
- `nanmean` -- mean of non-NaN values
- `nanvar` -- variance of non-NaN values
- `nanstd` -- standard deviation of non-NaN values
+- `nanmedian` -- median of non-NaN values
+- `nanpercentile` -- qth percentile of non-NaN values
"""
from __future__ import division, absolute_import, print_function
@@ -22,7 +25,7 @@ from numpy.lib.function_base import _ureduce as _ureduce
__all__ = [
'nansum', 'nanmax', 'nanmin', 'nanargmax', 'nanargmin', 'nanmean',
- 'nanmedian', 'nanpercentile', 'nanvar', 'nanstd'
+ 'nanmedian', 'nanpercentile', 'nanvar', 'nanstd', 'nanprod',
]
@@ -510,6 +513,76 @@ def nansum(a, axis=None, dtype=None, out=None, keepdims=0):
return np.sum(a, axis=axis, dtype=dtype, out=out, keepdims=keepdims)
+def nanprod(a, axis=None, dtype=None, out=None, keepdims=0):
+ """
+ Return the product of array elements over a given axis treating Not a
+ Numbers (NaNs) as zero.
+
+ One is returned for slices that are all-NaN or empty.
+
+ .. versionadded:: 1.10.0
+
+ Parameters
+ ----------
+ a : array_like
+ Array containing numbers whose sum is desired. If `a` is not an
+ array, a conversion is attempted.
+ axis : int, optional
+ Axis along which the product is computed. The default is to compute
+ the product of the flattened array.
+ dtype : data-type, optional
+ The type of the returned array and of the accumulator in which the
+ elements are summed. By default, the dtype of `a` is used. An
+ exception is when `a` has an integer type with less precision than
+ the platform (u)intp. In that case, the default will be either
+ (u)int32 or (u)int64 depending on whether the platform is 32 or 64
+ bits. For inexact inputs, dtype must be inexact.
+ out : ndarray, optional
+ Alternate output array in which to place the result. The default
+ is ``None``. If provided, it must have the same shape as the
+ expected output, but the type will be cast if necessary. See
+ `doc.ufuncs` for details. The casting of NaN to integer can yield
+ unexpected results.
+ keepdims : bool, optional
+ If True, the axes which are reduced are left in the result as
+ dimensions with size one. With this option, the result will
+ broadcast correctly against the original `arr`.
+
+ Returns
+ -------
+ y : ndarray or numpy scalar
+
+ See Also
+ --------
+ numpy.prod : Product across array propagating NaNs.
+ isnan : Show which elements are NaN.
+
+ Notes
+ -----
+ Numpy integer arithmetic is modular. If the size of a product exceeds
+ the size of an integer accumulator, its value will wrap around and the
+ result will be incorrect. Specifying ``dtype=double`` can alleviate
+ that problem.
+
+ Examples
+ --------
+ >>> np.nanprod(1)
+ 1
+ >>> np.nanprod([1])
+ 1
+ >>> np.nanprod([1, np.nan])
+ 1.0
+ >>> a = np.array([[1, 2], [3, np.nan]])
+ >>> np.nanprod(a)
+ 6.0
+ >>> np.nanprod(a, axis=0)
+ array([ 3., 2.])
+
+ """
+ a, mask = _replace_nan(a, 1)
+ return np.prod(a, axis=axis, dtype=dtype, out=out, keepdims=keepdims)
+
+
def nanmean(a, axis=None, dtype=None, out=None, keepdims=False):
"""
Compute the arithmetic mean along the specified axis, ignoring NaNs.
@@ -770,6 +843,8 @@ def nanpercentile(a, q, axis=None, out=None, overwrite_input=False,
Returns the qth percentile of the array elements.
+ .. versionadded:: 1.9.0
+
Parameters
----------
a : array_like
diff --git a/numpy/lib/tests/test_nanfunctions.py b/numpy/lib/tests/test_nanfunctions.py
index 35ae86c20..f418504c2 100644
--- a/numpy/lib/tests/test_nanfunctions.py
+++ b/numpy/lib/tests/test_nanfunctions.py
@@ -236,6 +236,11 @@ class TestNanFunctions_IntTypes(TestCase):
for mat in self.integer_arrays():
assert_equal(np.nansum(mat), tgt)
+ def test_nanprod(self):
+ tgt = np.prod(self.mat)
+ for mat in self.integer_arrays():
+ assert_equal(np.nanprod(mat), tgt)
+
def test_nanmean(self):
tgt = np.mean(self.mat)
for mat in self.integer_arrays():
@@ -260,119 +265,7 @@ class TestNanFunctions_IntTypes(TestCase):
assert_equal(np.nanstd(mat, ddof=1), tgt)
-class TestNanFunctions_Sum(TestCase):
-
- def test_mutation(self):
- # Check that passed array is not modified.
- ndat = _ndat.copy()
- np.nansum(ndat)
- assert_equal(ndat, _ndat)
-
- def test_keepdims(self):
- mat = np.eye(3)
- for axis in [None, 0, 1]:
- tgt = np.sum(mat, axis=axis, keepdims=True)
- res = np.nansum(mat, axis=axis, keepdims=True)
- assert_(res.ndim == tgt.ndim)
-
- def test_out(self):
- mat = np.eye(3)
- resout = np.zeros(3)
- tgt = np.sum(mat, axis=1)
- res = np.nansum(mat, axis=1, out=resout)
- assert_almost_equal(res, resout)
- assert_almost_equal(res, tgt)
-
- def test_dtype_from_dtype(self):
- mat = np.eye(3)
- codes = 'efdgFDG'
- for c in codes:
- tgt = np.sum(mat, dtype=np.dtype(c), axis=1).dtype.type
- res = np.nansum(mat, dtype=np.dtype(c), axis=1).dtype.type
- assert_(res is tgt)
- # scalar case
- tgt = np.sum(mat, dtype=np.dtype(c), axis=None).dtype.type
- res = np.nansum(mat, dtype=np.dtype(c), axis=None).dtype.type
- assert_(res is tgt)
-
- def test_dtype_from_char(self):
- mat = np.eye(3)
- codes = 'efdgFDG'
- for c in codes:
- tgt = np.sum(mat, dtype=c, axis=1).dtype.type
- res = np.nansum(mat, dtype=c, axis=1).dtype.type
- assert_(res is tgt)
- # scalar case
- tgt = np.sum(mat, dtype=c, axis=None).dtype.type
- res = np.nansum(mat, dtype=c, axis=None).dtype.type
- assert_(res is tgt)
-
- def test_dtype_from_input(self):
- codes = 'efdgFDG'
- for c in codes:
- mat = np.eye(3, dtype=c)
- tgt = np.sum(mat, axis=1).dtype.type
- res = np.nansum(mat, axis=1).dtype.type
- assert_(res is tgt)
- # scalar case
- tgt = np.sum(mat, axis=None).dtype.type
- res = np.nansum(mat, axis=None).dtype.type
- assert_(res is tgt)
-
- def test_result_values(self):
- tgt = [np.sum(d) for d in _rdat]
- res = np.nansum(_ndat, axis=1)
- assert_almost_equal(res, tgt)
-
- def test_allnans(self):
- # Check for FutureWarning
- with warnings.catch_warnings(record=True) as w:
- warnings.simplefilter('always')
- res = np.nansum([np.nan]*3, axis=None)
- assert_(res == 0, 'result is not 0')
- assert_(len(w) == 0, 'warning raised')
- # Check scalar
- res = np.nansum(np.nan)
- assert_(res == 0, 'result is not 0')
- assert_(len(w) == 0, 'warning raised')
- # Check there is no warning for not all-nan
- np.nansum([0]*3, axis=None)
- assert_(len(w) == 0, 'unwanted warning raised')
-
- def test_empty(self):
- mat = np.zeros((0, 3))
- tgt = [0]*3
- res = np.nansum(mat, axis=0)
- assert_equal(res, tgt)
- tgt = []
- res = np.nansum(mat, axis=1)
- assert_equal(res, tgt)
- tgt = 0
- res = np.nansum(mat, axis=None)
- assert_equal(res, tgt)
-
- def test_scalar(self):
- assert_(np.nansum(0.) == 0.)
-
- def test_matrices(self):
- # Check that it works and that type and
- # shape are preserved
- mat = np.matrix(np.eye(3))
- res = np.nansum(mat, axis=0)
- assert_(isinstance(res, np.matrix))
- assert_(res.shape == (1, 3))
- res = np.nansum(mat, axis=1)
- assert_(isinstance(res, np.matrix))
- assert_(res.shape == (3, 1))
- res = np.nansum(mat)
- assert_(np.isscalar(res))
-
-
-class TestNanFunctions_MeanVarStd(TestCase):
-
- nanfuncs = [np.nanmean, np.nanvar, np.nanstd]
- stdfuncs = [np.mean, np.var, np.std]
-
+class SharedNanFunctionsTestsMixin(object):
def test_mutation(self):
# Check that passed array is not modified.
ndat = _ndat.copy()
@@ -380,17 +273,6 @@ class TestNanFunctions_MeanVarStd(TestCase):
f(ndat)
assert_equal(ndat, _ndat)
- def test_dtype_error(self):
- for f in self.nanfuncs:
- for dtype in [np.bool_, np.int_, np.object]:
- assert_raises(TypeError, f, _ndat, axis=1, dtype=np.int)
-
- def test_out_dtype_error(self):
- for f in self.nanfuncs:
- for dtype in [np.bool_, np.int_, np.object]:
- out = np.empty(_ndat.shape[0], dtype=dtype)
- assert_raises(TypeError, f, _ndat, axis=1, out=out)
-
def test_keepdims(self):
mat = np.eye(3)
for nf, rf in zip(self.nanfuncs, self.stdfuncs):
@@ -447,6 +329,81 @@ class TestNanFunctions_MeanVarStd(TestCase):
res = nf(mat, axis=None).dtype.type
assert_(res is tgt)
+ def test_result_values(self):
+ for nf, rf in zip(self.nanfuncs, self.stdfuncs):
+ tgt = [rf(d) for d in _rdat]
+ res = nf(_ndat, axis=1)
+ assert_almost_equal(res, tgt)
+
+ def test_scalar(self):
+ for f in self.nanfuncs:
+ assert_(f(0.) == 0.)
+
+ def test_matrices(self):
+ # Check that it works and that type and
+ # shape are preserved
+ mat = np.matrix(np.eye(3))
+ for f in self.nanfuncs:
+ res = f(mat, axis=0)
+ assert_(isinstance(res, np.matrix))
+ assert_(res.shape == (1, 3))
+ res = f(mat, axis=1)
+ assert_(isinstance(res, np.matrix))
+ assert_(res.shape == (3, 1))
+ res = f(mat)
+ assert_(np.isscalar(res))
+
+
+class TestNanFunctions_SumProd(TestCase, SharedNanFunctionsTestsMixin):
+
+ nanfuncs = [np.nansum, np.nanprod]
+ stdfuncs = [np.sum, np.prod]
+
+ def test_allnans(self):
+ # Check for FutureWarning
+ with warnings.catch_warnings(record=True) as w:
+ warnings.simplefilter('always')
+ res = np.nansum([np.nan]*3, axis=None)
+ assert_(res == 0, 'result is not 0')
+ assert_(len(w) == 0, 'warning raised')
+ # Check scalar
+ res = np.nansum(np.nan)
+ assert_(res == 0, 'result is not 0')
+ assert_(len(w) == 0, 'warning raised')
+ # Check there is no warning for not all-nan
+ np.nansum([0]*3, axis=None)
+ assert_(len(w) == 0, 'unwanted warning raised')
+
+ def test_empty(self):
+ for f, tgt_value in zip([np.nansum, np.nanprod], [0, 1]):
+ mat = np.zeros((0, 3))
+ tgt = [tgt_value]*3
+ res = f(mat, axis=0)
+ assert_equal(res, tgt)
+ tgt = []
+ res = f(mat, axis=1)
+ assert_equal(res, tgt)
+ tgt = tgt_value
+ res = f(mat, axis=None)
+ assert_equal(res, tgt)
+
+
+class TestNanFunctions_MeanVarStd(TestCase, SharedNanFunctionsTestsMixin):
+
+ nanfuncs = [np.nanmean, np.nanvar, np.nanstd]
+ stdfuncs = [np.mean, np.var, np.std]
+
+ def test_dtype_error(self):
+ for f in self.nanfuncs:
+ for dtype in [np.bool_, np.int_, np.object]:
+ assert_raises(TypeError, f, _ndat, axis=1, dtype=np.int)
+
+ def test_out_dtype_error(self):
+ for f in self.nanfuncs:
+ for dtype in [np.bool_, np.int_, np.object]:
+ out = np.empty(_ndat.shape[0], dtype=dtype)
+ assert_raises(TypeError, f, _ndat, axis=1, out=out)
+
def test_ddof(self):
nanfuncs = [np.nanvar, np.nanstd]
stdfuncs = [np.var, np.std]
@@ -473,12 +430,6 @@ class TestNanFunctions_MeanVarStd(TestCase):
else:
assert_(len(w) == 0)
- def test_result_values(self):
- for nf, rf in zip(self.nanfuncs, self.stdfuncs):
- tgt = [rf(d) for d in _rdat]
- res = nf(_ndat, axis=1)
- assert_almost_equal(res, tgt)
-
def test_allnans(self):
mat = np.array([np.nan]*9).reshape(3, 3)
for f in self.nanfuncs:
@@ -508,24 +459,6 @@ class TestNanFunctions_MeanVarStd(TestCase):
assert_equal(f(mat, axis=axis), np.zeros([]))
assert_(len(w) == 0)
- def test_scalar(self):
- for f in self.nanfuncs:
- assert_(f(0.) == 0.)
-
- def test_matrices(self):
- # Check that it works and that type and
- # shape are preserved
- mat = np.matrix(np.eye(3))
- for f in self.nanfuncs:
- res = f(mat, axis=0)
- assert_(isinstance(res, np.matrix))
- assert_(res.shape == (1, 3))
- res = f(mat, axis=1)
- assert_(isinstance(res, np.matrix))
- assert_(res.shape == (3, 1))
- res = f(mat)
- assert_(np.isscalar(res))
-
class TestNanFunctions_Median(TestCase):