Coverage for /pythoncovmergedfiles/medio/medio/usr/local/lib/python3.8/site-packages/numpy/ma/extras.py: 16%
Shortcuts on this page
r m x toggle line displays
j k next/prev highlighted chunk
0 (zero) top of page
1 (one) first highlighted chunk
1"""
2Masked arrays add-ons.
4A collection of utilities for `numpy.ma`.
6:author: Pierre Gerard-Marchant
7:contact: pierregm_at_uga_dot_edu
8:version: $Id: extras.py 3473 2007-10-29 15:18:13Z jarrod.millman $
10"""
11__all__ = [
12 'apply_along_axis', 'apply_over_axes', 'atleast_1d', 'atleast_2d',
13 'atleast_3d', 'average', 'clump_masked', 'clump_unmasked', 'column_stack',
14 'compress_cols', 'compress_nd', 'compress_rowcols', 'compress_rows',
15 'count_masked', 'corrcoef', 'cov', 'diagflat', 'dot', 'dstack', 'ediff1d',
16 'flatnotmasked_contiguous', 'flatnotmasked_edges', 'hsplit', 'hstack',
17 'isin', 'in1d', 'intersect1d', 'mask_cols', 'mask_rowcols', 'mask_rows',
18 'masked_all', 'masked_all_like', 'median', 'mr_', 'ndenumerate',
19 'notmasked_contiguous', 'notmasked_edges', 'polyfit', 'row_stack',
20 'setdiff1d', 'setxor1d', 'stack', 'unique', 'union1d', 'vander', 'vstack',
21 ]
23import itertools
24import warnings
26from . import core as ma
27from .core import (
28 MaskedArray, MAError, add, array, asarray, concatenate, filled, count,
29 getmask, getmaskarray, make_mask_descr, masked, masked_array, mask_or,
30 nomask, ones, sort, zeros, getdata, get_masked_subclass, dot,
31 mask_rowcols
32 )
34import numpy as np
35from numpy import ndarray, array as nxarray
36from numpy.core.multiarray import normalize_axis_index
37from numpy.core.numeric import normalize_axis_tuple
38from numpy.lib.function_base import _ureduce
39from numpy.lib.index_tricks import AxisConcatenator
42def issequence(seq):
43 """
44 Is seq a sequence (ndarray, list or tuple)?
46 """
47 return isinstance(seq, (ndarray, tuple, list))
50def count_masked(arr, axis=None):
51 """
52 Count the number of masked elements along the given axis.
54 Parameters
55 ----------
56 arr : array_like
57 An array with (possibly) masked elements.
58 axis : int, optional
59 Axis along which to count. If None (default), a flattened
60 version of the array is used.
62 Returns
63 -------
64 count : int, ndarray
65 The total number of masked elements (axis=None) or the number
66 of masked elements along each slice of the given axis.
68 See Also
69 --------
70 MaskedArray.count : Count non-masked elements.
72 Examples
73 --------
74 >>> import numpy.ma as ma
75 >>> a = np.arange(9).reshape((3,3))
76 >>> a = ma.array(a)
77 >>> a[1, 0] = ma.masked
78 >>> a[1, 2] = ma.masked
79 >>> a[2, 1] = ma.masked
80 >>> a
81 masked_array(
82 data=[[0, 1, 2],
83 [--, 4, --],
84 [6, --, 8]],
85 mask=[[False, False, False],
86 [ True, False, True],
87 [False, True, False]],
88 fill_value=999999)
89 >>> ma.count_masked(a)
90 3
92 When the `axis` keyword is used an array is returned.
94 >>> ma.count_masked(a, axis=0)
95 array([1, 1, 1])
96 >>> ma.count_masked(a, axis=1)
97 array([0, 2, 1])
99 """
100 m = getmaskarray(arr)
101 return m.sum(axis)
104def masked_all(shape, dtype=float):
105 """
106 Empty masked array with all elements masked.
108 Return an empty masked array of the given shape and dtype, where all the
109 data are masked.
111 Parameters
112 ----------
113 shape : int or tuple of ints
114 Shape of the required MaskedArray, e.g., ``(2, 3)`` or ``2``.
115 dtype : dtype, optional
116 Data type of the output.
118 Returns
119 -------
120 a : MaskedArray
121 A masked array with all data masked.
123 See Also
124 --------
125 masked_all_like : Empty masked array modelled on an existing array.
127 Examples
128 --------
129 >>> import numpy.ma as ma
130 >>> ma.masked_all((3, 3))
131 masked_array(
132 data=[[--, --, --],
133 [--, --, --],
134 [--, --, --]],
135 mask=[[ True, True, True],
136 [ True, True, True],
137 [ True, True, True]],
138 fill_value=1e+20,
139 dtype=float64)
141 The `dtype` parameter defines the underlying data type.
143 >>> a = ma.masked_all((3, 3))
144 >>> a.dtype
145 dtype('float64')
146 >>> a = ma.masked_all((3, 3), dtype=np.int32)
147 >>> a.dtype
148 dtype('int32')
150 """
151 a = masked_array(np.empty(shape, dtype),
152 mask=np.ones(shape, make_mask_descr(dtype)))
153 return a
156def masked_all_like(arr):
157 """
158 Empty masked array with the properties of an existing array.
160 Return an empty masked array of the same shape and dtype as
161 the array `arr`, where all the data are masked.
163 Parameters
164 ----------
165 arr : ndarray
166 An array describing the shape and dtype of the required MaskedArray.
168 Returns
169 -------
170 a : MaskedArray
171 A masked array with all data masked.
173 Raises
174 ------
175 AttributeError
176 If `arr` doesn't have a shape attribute (i.e. not an ndarray)
178 See Also
179 --------
180 masked_all : Empty masked array with all elements masked.
182 Examples
183 --------
184 >>> import numpy.ma as ma
185 >>> arr = np.zeros((2, 3), dtype=np.float32)
186 >>> arr
187 array([[0., 0., 0.],
188 [0., 0., 0.]], dtype=float32)
189 >>> ma.masked_all_like(arr)
190 masked_array(
191 data=[[--, --, --],
192 [--, --, --]],
193 mask=[[ True, True, True],
194 [ True, True, True]],
195 fill_value=1e+20,
196 dtype=float32)
198 The dtype of the masked array matches the dtype of `arr`.
200 >>> arr.dtype
201 dtype('float32')
202 >>> ma.masked_all_like(arr).dtype
203 dtype('float32')
205 """
206 a = np.empty_like(arr).view(MaskedArray)
207 a._mask = np.ones(a.shape, dtype=make_mask_descr(a.dtype))
208 return a
211#####--------------------------------------------------------------------------
212#---- --- Standard functions ---
213#####--------------------------------------------------------------------------
214class _fromnxfunction:
215 """
216 Defines a wrapper to adapt NumPy functions to masked arrays.
219 An instance of `_fromnxfunction` can be called with the same parameters
220 as the wrapped NumPy function. The docstring of `newfunc` is adapted from
221 the wrapped function as well, see `getdoc`.
223 This class should not be used directly. Instead, one of its extensions that
224 provides support for a specific type of input should be used.
226 Parameters
227 ----------
228 funcname : str
229 The name of the function to be adapted. The function should be
230 in the NumPy namespace (i.e. ``np.funcname``).
232 """
234 def __init__(self, funcname):
235 self.__name__ = funcname
236 self.__doc__ = self.getdoc()
238 def getdoc(self):
239 """
240 Retrieve the docstring and signature from the function.
242 The ``__doc__`` attribute of the function is used as the docstring for
243 the new masked array version of the function. A note on application
244 of the function to the mask is appended.
246 Parameters
247 ----------
248 None
250 """
251 npfunc = getattr(np, self.__name__, None)
252 doc = getattr(npfunc, '__doc__', None)
253 if doc:
254 sig = self.__name__ + ma.get_object_signature(npfunc)
255 doc = ma.doc_note(doc, "The function is applied to both the _data "
256 "and the _mask, if any.")
257 return '\n\n'.join((sig, doc))
258 return
260 def __call__(self, *args, **params):
261 pass
264class _fromnxfunction_single(_fromnxfunction):
265 """
266 A version of `_fromnxfunction` that is called with a single array
267 argument followed by auxiliary args that are passed verbatim for
268 both the data and mask calls.
269 """
270 def __call__(self, x, *args, **params):
271 func = getattr(np, self.__name__)
272 if isinstance(x, ndarray):
273 _d = func(x.__array__(), *args, **params)
274 _m = func(getmaskarray(x), *args, **params)
275 return masked_array(_d, mask=_m)
276 else:
277 _d = func(np.asarray(x), *args, **params)
278 _m = func(getmaskarray(x), *args, **params)
279 return masked_array(_d, mask=_m)
282class _fromnxfunction_seq(_fromnxfunction):
283 """
284 A version of `_fromnxfunction` that is called with a single sequence
285 of arrays followed by auxiliary args that are passed verbatim for
286 both the data and mask calls.
287 """
288 def __call__(self, x, *args, **params):
289 func = getattr(np, self.__name__)
290 _d = func(tuple([np.asarray(a) for a in x]), *args, **params)
291 _m = func(tuple([getmaskarray(a) for a in x]), *args, **params)
292 return masked_array(_d, mask=_m)
295class _fromnxfunction_args(_fromnxfunction):
296 """
297 A version of `_fromnxfunction` that is called with multiple array
298 arguments. The first non-array-like input marks the beginning of the
299 arguments that are passed verbatim for both the data and mask calls.
300 Array arguments are processed independently and the results are
301 returned in a list. If only one array is found, the return value is
302 just the processed array instead of a list.
303 """
304 def __call__(self, *args, **params):
305 func = getattr(np, self.__name__)
306 arrays = []
307 args = list(args)
308 while len(args) > 0 and issequence(args[0]):
309 arrays.append(args.pop(0))
310 res = []
311 for x in arrays:
312 _d = func(np.asarray(x), *args, **params)
313 _m = func(getmaskarray(x), *args, **params)
314 res.append(masked_array(_d, mask=_m))
315 if len(arrays) == 1:
316 return res[0]
317 return res
320class _fromnxfunction_allargs(_fromnxfunction):
321 """
322 A version of `_fromnxfunction` that is called with multiple array
323 arguments. Similar to `_fromnxfunction_args` except that all args
324 are converted to arrays even if they are not so already. This makes
325 it possible to process scalars as 1-D arrays. Only keyword arguments
326 are passed through verbatim for the data and mask calls. Arrays
327 arguments are processed independently and the results are returned
328 in a list. If only one arg is present, the return value is just the
329 processed array instead of a list.
330 """
331 def __call__(self, *args, **params):
332 func = getattr(np, self.__name__)
333 res = []
334 for x in args:
335 _d = func(np.asarray(x), **params)
336 _m = func(getmaskarray(x), **params)
337 res.append(masked_array(_d, mask=_m))
338 if len(args) == 1:
339 return res[0]
340 return res
343atleast_1d = _fromnxfunction_allargs('atleast_1d')
344atleast_2d = _fromnxfunction_allargs('atleast_2d')
345atleast_3d = _fromnxfunction_allargs('atleast_3d')
347vstack = row_stack = _fromnxfunction_seq('vstack')
348hstack = _fromnxfunction_seq('hstack')
349column_stack = _fromnxfunction_seq('column_stack')
350dstack = _fromnxfunction_seq('dstack')
351stack = _fromnxfunction_seq('stack')
353hsplit = _fromnxfunction_single('hsplit')
355diagflat = _fromnxfunction_single('diagflat')
358#####--------------------------------------------------------------------------
359#----
360#####--------------------------------------------------------------------------
361def flatten_inplace(seq):
362 """Flatten a sequence in place."""
363 k = 0
364 while (k != len(seq)):
365 while hasattr(seq[k], '__iter__'):
366 seq[k:(k + 1)] = seq[k]
367 k += 1
368 return seq
371def apply_along_axis(func1d, axis, arr, *args, **kwargs):
372 """
373 (This docstring should be overwritten)
374 """
375 arr = array(arr, copy=False, subok=True)
376 nd = arr.ndim
377 axis = normalize_axis_index(axis, nd)
378 ind = [0] * (nd - 1)
379 i = np.zeros(nd, 'O')
380 indlist = list(range(nd))
381 indlist.remove(axis)
382 i[axis] = slice(None, None)
383 outshape = np.asarray(arr.shape).take(indlist)
384 i.put(indlist, ind)
385 res = func1d(arr[tuple(i.tolist())], *args, **kwargs)
386 # if res is a number, then we have a smaller output array
387 asscalar = np.isscalar(res)
388 if not asscalar:
389 try:
390 len(res)
391 except TypeError:
392 asscalar = True
393 # Note: we shouldn't set the dtype of the output from the first result
394 # so we force the type to object, and build a list of dtypes. We'll
395 # just take the largest, to avoid some downcasting
396 dtypes = []
397 if asscalar:
398 dtypes.append(np.asarray(res).dtype)
399 outarr = zeros(outshape, object)
400 outarr[tuple(ind)] = res
401 Ntot = np.product(outshape)
402 k = 1
403 while k < Ntot:
404 # increment the index
405 ind[-1] += 1
406 n = -1
407 while (ind[n] >= outshape[n]) and (n > (1 - nd)):
408 ind[n - 1] += 1
409 ind[n] = 0
410 n -= 1
411 i.put(indlist, ind)
412 res = func1d(arr[tuple(i.tolist())], *args, **kwargs)
413 outarr[tuple(ind)] = res
414 dtypes.append(asarray(res).dtype)
415 k += 1
416 else:
417 res = array(res, copy=False, subok=True)
418 j = i.copy()
419 j[axis] = ([slice(None, None)] * res.ndim)
420 j.put(indlist, ind)
421 Ntot = np.product(outshape)
422 holdshape = outshape
423 outshape = list(arr.shape)
424 outshape[axis] = res.shape
425 dtypes.append(asarray(res).dtype)
426 outshape = flatten_inplace(outshape)
427 outarr = zeros(outshape, object)
428 outarr[tuple(flatten_inplace(j.tolist()))] = res
429 k = 1
430 while k < Ntot:
431 # increment the index
432 ind[-1] += 1
433 n = -1
434 while (ind[n] >= holdshape[n]) and (n > (1 - nd)):
435 ind[n - 1] += 1
436 ind[n] = 0
437 n -= 1
438 i.put(indlist, ind)
439 j.put(indlist, ind)
440 res = func1d(arr[tuple(i.tolist())], *args, **kwargs)
441 outarr[tuple(flatten_inplace(j.tolist()))] = res
442 dtypes.append(asarray(res).dtype)
443 k += 1
444 max_dtypes = np.dtype(np.asarray(dtypes).max())
445 if not hasattr(arr, '_mask'):
446 result = np.asarray(outarr, dtype=max_dtypes)
447 else:
448 result = asarray(outarr, dtype=max_dtypes)
449 result.fill_value = ma.default_fill_value(result)
450 return result
451apply_along_axis.__doc__ = np.apply_along_axis.__doc__
454def apply_over_axes(func, a, axes):
455 """
456 (This docstring will be overwritten)
457 """
458 val = asarray(a)
459 N = a.ndim
460 if array(axes).ndim == 0:
461 axes = (axes,)
462 for axis in axes:
463 if axis < 0:
464 axis = N + axis
465 args = (val, axis)
466 res = func(*args)
467 if res.ndim == val.ndim:
468 val = res
469 else:
470 res = ma.expand_dims(res, axis)
471 if res.ndim == val.ndim:
472 val = res
473 else:
474 raise ValueError("function is not returning "
475 "an array of the correct shape")
476 return val
479if apply_over_axes.__doc__ is not None:
480 apply_over_axes.__doc__ = np.apply_over_axes.__doc__[
481 :np.apply_over_axes.__doc__.find('Notes')].rstrip() + \
482 """
484 Examples
485 --------
486 >>> a = np.ma.arange(24).reshape(2,3,4)
487 >>> a[:,0,1] = np.ma.masked
488 >>> a[:,1,:] = np.ma.masked
489 >>> a
490 masked_array(
491 data=[[[0, --, 2, 3],
492 [--, --, --, --],
493 [8, 9, 10, 11]],
494 [[12, --, 14, 15],
495 [--, --, --, --],
496 [20, 21, 22, 23]]],
497 mask=[[[False, True, False, False],
498 [ True, True, True, True],
499 [False, False, False, False]],
500 [[False, True, False, False],
501 [ True, True, True, True],
502 [False, False, False, False]]],
503 fill_value=999999)
504 >>> np.ma.apply_over_axes(np.ma.sum, a, [0,2])
505 masked_array(
506 data=[[[46],
507 [--],
508 [124]]],
509 mask=[[[False],
510 [ True],
511 [False]]],
512 fill_value=999999)
514 Tuple axis arguments to ufuncs are equivalent:
516 >>> np.ma.sum(a, axis=(0,2)).reshape((1,-1,1))
517 masked_array(
518 data=[[[46],
519 [--],
520 [124]]],
521 mask=[[[False],
522 [ True],
523 [False]]],
524 fill_value=999999)
525 """
528def average(a, axis=None, weights=None, returned=False, *,
529 keepdims=np._NoValue):
530 """
531 Return the weighted average of array over the given axis.
533 Parameters
534 ----------
535 a : array_like
536 Data to be averaged.
537 Masked entries are not taken into account in the computation.
538 axis : int, optional
539 Axis along which to average `a`. If None, averaging is done over
540 the flattened array.
541 weights : array_like, optional
542 The importance that each element has in the computation of the average.
543 The weights array can either be 1-D (in which case its length must be
544 the size of `a` along the given axis) or of the same shape as `a`.
545 If ``weights=None``, then all data in `a` are assumed to have a
546 weight equal to one. The 1-D calculation is::
548 avg = sum(a * weights) / sum(weights)
550 The only constraint on `weights` is that `sum(weights)` must not be 0.
551 returned : bool, optional
552 Flag indicating whether a tuple ``(result, sum of weights)``
553 should be returned as output (True), or just the result (False).
554 Default is False.
555 keepdims : bool, optional
556 If this is set to True, the axes which are reduced are left
557 in the result as dimensions with size one. With this option,
558 the result will broadcast correctly against the original `a`.
559 *Note:* `keepdims` will not work with instances of `numpy.matrix`
560 or other classes whose methods do not support `keepdims`.
562 .. versionadded:: 1.23.0
564 Returns
565 -------
566 average, [sum_of_weights] : (tuple of) scalar or MaskedArray
567 The average along the specified axis. When returned is `True`,
568 return a tuple with the average as the first element and the sum
569 of the weights as the second element. The return type is `np.float64`
570 if `a` is of integer type and floats smaller than `float64`, or the
571 input data-type, otherwise. If returned, `sum_of_weights` is always
572 `float64`.
574 Examples
575 --------
576 >>> a = np.ma.array([1., 2., 3., 4.], mask=[False, False, True, True])
577 >>> np.ma.average(a, weights=[3, 1, 0, 0])
578 1.25
580 >>> x = np.ma.arange(6.).reshape(3, 2)
581 >>> x
582 masked_array(
583 data=[[0., 1.],
584 [2., 3.],
585 [4., 5.]],
586 mask=False,
587 fill_value=1e+20)
588 >>> avg, sumweights = np.ma.average(x, axis=0, weights=[1, 2, 3],
589 ... returned=True)
590 >>> avg
591 masked_array(data=[2.6666666666666665, 3.6666666666666665],
592 mask=[False, False],
593 fill_value=1e+20)
595 With ``keepdims=True``, the following result has shape (3, 1).
597 >>> np.ma.average(x, axis=1, keepdims=True)
598 masked_array(
599 data=[[0.5],
600 [2.5],
601 [4.5]],
602 mask=False,
603 fill_value=1e+20)
604 """
605 a = asarray(a)
606 m = getmask(a)
608 # inspired by 'average' in numpy/lib/function_base.py
610 if keepdims is np._NoValue:
611 # Don't pass on the keepdims argument if one wasn't given.
612 keepdims_kw = {}
613 else:
614 keepdims_kw = {'keepdims': keepdims}
616 if weights is None:
617 avg = a.mean(axis, **keepdims_kw)
618 scl = avg.dtype.type(a.count(axis))
619 else:
620 wgt = asarray(weights)
622 if issubclass(a.dtype.type, (np.integer, np.bool_)):
623 result_dtype = np.result_type(a.dtype, wgt.dtype, 'f8')
624 else:
625 result_dtype = np.result_type(a.dtype, wgt.dtype)
627 # Sanity checks
628 if a.shape != wgt.shape:
629 if axis is None:
630 raise TypeError(
631 "Axis must be specified when shapes of a and weights "
632 "differ.")
633 if wgt.ndim != 1:
634 raise TypeError(
635 "1D weights expected when shapes of a and weights differ.")
636 if wgt.shape[0] != a.shape[axis]:
637 raise ValueError(
638 "Length of weights not compatible with specified axis.")
640 # setup wgt to broadcast along axis
641 wgt = np.broadcast_to(wgt, (a.ndim-1)*(1,) + wgt.shape, subok=True)
642 wgt = wgt.swapaxes(-1, axis)
644 if m is not nomask:
645 wgt = wgt*(~a.mask)
646 wgt.mask |= a.mask
648 scl = wgt.sum(axis=axis, dtype=result_dtype, **keepdims_kw)
649 avg = np.multiply(a, wgt,
650 dtype=result_dtype).sum(axis, **keepdims_kw) / scl
652 if returned:
653 if scl.shape != avg.shape:
654 scl = np.broadcast_to(scl, avg.shape).copy()
655 return avg, scl
656 else:
657 return avg
660def median(a, axis=None, out=None, overwrite_input=False, keepdims=False):
661 """
662 Compute the median along the specified axis.
664 Returns the median of the array elements.
666 Parameters
667 ----------
668 a : array_like
669 Input array or object that can be converted to an array.
670 axis : int, optional
671 Axis along which the medians are computed. The default (None) is
672 to compute the median along a flattened version of the array.
673 out : ndarray, optional
674 Alternative output array in which to place the result. It must
675 have the same shape and buffer length as the expected output
676 but the type will be cast if necessary.
677 overwrite_input : bool, optional
678 If True, then allow use of memory of input array (a) for
679 calculations. The input array will be modified by the call to
680 median. This will save memory when you do not need to preserve
681 the contents of the input array. Treat the input as undefined,
682 but it will probably be fully or partially sorted. Default is
683 False. Note that, if `overwrite_input` is True, and the input
684 is not already an `ndarray`, an error will be raised.
685 keepdims : bool, optional
686 If this is set to True, the axes which are reduced are left
687 in the result as dimensions with size one. With this option,
688 the result will broadcast correctly against the input array.
690 .. versionadded:: 1.10.0
692 Returns
693 -------
694 median : ndarray
695 A new array holding the result is returned unless out is
696 specified, in which case a reference to out is returned.
697 Return data-type is `float64` for integers and floats smaller than
698 `float64`, or the input data-type, otherwise.
700 See Also
701 --------
702 mean
704 Notes
705 -----
706 Given a vector ``V`` with ``N`` non masked values, the median of ``V``
707 is the middle value of a sorted copy of ``V`` (``Vs``) - i.e.
708 ``Vs[(N-1)/2]``, when ``N`` is odd, or ``{Vs[N/2 - 1] + Vs[N/2]}/2``
709 when ``N`` is even.
711 Examples
712 --------
713 >>> x = np.ma.array(np.arange(8), mask=[0]*4 + [1]*4)
714 >>> np.ma.median(x)
715 1.5
717 >>> x = np.ma.array(np.arange(10).reshape(2, 5), mask=[0]*6 + [1]*4)
718 >>> np.ma.median(x)
719 2.5
720 >>> np.ma.median(x, axis=-1, overwrite_input=True)
721 masked_array(data=[2.0, 5.0],
722 mask=[False, False],
723 fill_value=1e+20)
725 """
726 if not hasattr(a, 'mask'):
727 m = np.median(getdata(a, subok=True), axis=axis,
728 out=out, overwrite_input=overwrite_input,
729 keepdims=keepdims)
730 if isinstance(m, np.ndarray) and 1 <= m.ndim:
731 return masked_array(m, copy=False)
732 else:
733 return m
735 return _ureduce(a, func=_median, keepdims=keepdims, axis=axis, out=out,
736 overwrite_input=overwrite_input)
739def _median(a, axis=None, out=None, overwrite_input=False):
740 # when an unmasked NaN is present return it, so we need to sort the NaN
741 # values behind the mask
742 if np.issubdtype(a.dtype, np.inexact):
743 fill_value = np.inf
744 else:
745 fill_value = None
746 if overwrite_input:
747 if axis is None:
748 asorted = a.ravel()
749 asorted.sort(fill_value=fill_value)
750 else:
751 a.sort(axis=axis, fill_value=fill_value)
752 asorted = a
753 else:
754 asorted = sort(a, axis=axis, fill_value=fill_value)
756 if axis is None:
757 axis = 0
758 else:
759 axis = normalize_axis_index(axis, asorted.ndim)
761 if asorted.shape[axis] == 0:
762 # for empty axis integer indices fail so use slicing to get same result
763 # as median (which is mean of empty slice = nan)
764 indexer = [slice(None)] * asorted.ndim
765 indexer[axis] = slice(0, 0)
766 indexer = tuple(indexer)
767 return np.ma.mean(asorted[indexer], axis=axis, out=out)
769 if asorted.ndim == 1:
770 idx, odd = divmod(count(asorted), 2)
771 mid = asorted[idx + odd - 1:idx + 1]
772 if np.issubdtype(asorted.dtype, np.inexact) and asorted.size > 0:
773 # avoid inf / x = masked
774 s = mid.sum(out=out)
775 if not odd:
776 s = np.true_divide(s, 2., casting='safe', out=out)
777 s = np.lib.utils._median_nancheck(asorted, s, axis)
778 else:
779 s = mid.mean(out=out)
781 # if result is masked either the input contained enough
782 # minimum_fill_value so that it would be the median or all values
783 # masked
784 if np.ma.is_masked(s) and not np.all(asorted.mask):
785 return np.ma.minimum_fill_value(asorted)
786 return s
788 counts = count(asorted, axis=axis, keepdims=True)
789 h = counts // 2
791 # duplicate high if odd number of elements so mean does nothing
792 odd = counts % 2 == 1
793 l = np.where(odd, h, h-1)
795 lh = np.concatenate([l,h], axis=axis)
797 # get low and high median
798 low_high = np.take_along_axis(asorted, lh, axis=axis)
800 def replace_masked(s):
801 # Replace masked entries with minimum_full_value unless it all values
802 # are masked. This is required as the sort order of values equal or
803 # larger than the fill value is undefined and a valid value placed
804 # elsewhere, e.g. [4, --, inf].
805 if np.ma.is_masked(s):
806 rep = (~np.all(asorted.mask, axis=axis, keepdims=True)) & s.mask
807 s.data[rep] = np.ma.minimum_fill_value(asorted)
808 s.mask[rep] = False
810 replace_masked(low_high)
812 if np.issubdtype(asorted.dtype, np.inexact):
813 # avoid inf / x = masked
814 s = np.ma.sum(low_high, axis=axis, out=out)
815 np.true_divide(s.data, 2., casting='unsafe', out=s.data)
817 s = np.lib.utils._median_nancheck(asorted, s, axis)
818 else:
819 s = np.ma.mean(low_high, axis=axis, out=out)
821 return s
824def compress_nd(x, axis=None):
825 """Suppress slices from multiple dimensions which contain masked values.
827 Parameters
828 ----------
829 x : array_like, MaskedArray
830 The array to operate on. If not a MaskedArray instance (or if no array
831 elements are masked), `x` is interpreted as a MaskedArray with `mask`
832 set to `nomask`.
833 axis : tuple of ints or int, optional
834 Which dimensions to suppress slices from can be configured with this
835 parameter.
836 - If axis is a tuple of ints, those are the axes to suppress slices from.
837 - If axis is an int, then that is the only axis to suppress slices from.
838 - If axis is None, all axis are selected.
840 Returns
841 -------
842 compress_array : ndarray
843 The compressed array.
844 """
845 x = asarray(x)
846 m = getmask(x)
847 # Set axis to tuple of ints
848 if axis is None:
849 axis = tuple(range(x.ndim))
850 else:
851 axis = normalize_axis_tuple(axis, x.ndim)
853 # Nothing is masked: return x
854 if m is nomask or not m.any():
855 return x._data
856 # All is masked: return empty
857 if m.all():
858 return nxarray([])
859 # Filter elements through boolean indexing
860 data = x._data
861 for ax in axis:
862 axes = tuple(list(range(ax)) + list(range(ax + 1, x.ndim)))
863 data = data[(slice(None),)*ax + (~m.any(axis=axes),)]
864 return data
867def compress_rowcols(x, axis=None):
868 """
869 Suppress the rows and/or columns of a 2-D array that contain
870 masked values.
872 The suppression behavior is selected with the `axis` parameter.
874 - If axis is None, both rows and columns are suppressed.
875 - If axis is 0, only rows are suppressed.
876 - If axis is 1 or -1, only columns are suppressed.
878 Parameters
879 ----------
880 x : array_like, MaskedArray
881 The array to operate on. If not a MaskedArray instance (or if no array
882 elements are masked), `x` is interpreted as a MaskedArray with
883 `mask` set to `nomask`. Must be a 2D array.
884 axis : int, optional
885 Axis along which to perform the operation. Default is None.
887 Returns
888 -------
889 compressed_array : ndarray
890 The compressed array.
892 Examples
893 --------
894 >>> x = np.ma.array(np.arange(9).reshape(3, 3), mask=[[1, 0, 0],
895 ... [1, 0, 0],
896 ... [0, 0, 0]])
897 >>> x
898 masked_array(
899 data=[[--, 1, 2],
900 [--, 4, 5],
901 [6, 7, 8]],
902 mask=[[ True, False, False],
903 [ True, False, False],
904 [False, False, False]],
905 fill_value=999999)
907 >>> np.ma.compress_rowcols(x)
908 array([[7, 8]])
909 >>> np.ma.compress_rowcols(x, 0)
910 array([[6, 7, 8]])
911 >>> np.ma.compress_rowcols(x, 1)
912 array([[1, 2],
913 [4, 5],
914 [7, 8]])
916 """
917 if asarray(x).ndim != 2:
918 raise NotImplementedError("compress_rowcols works for 2D arrays only.")
919 return compress_nd(x, axis=axis)
922def compress_rows(a):
923 """
924 Suppress whole rows of a 2-D array that contain masked values.
926 This is equivalent to ``np.ma.compress_rowcols(a, 0)``, see
927 `compress_rowcols` for details.
929 See Also
930 --------
931 compress_rowcols
933 """
934 a = asarray(a)
935 if a.ndim != 2:
936 raise NotImplementedError("compress_rows works for 2D arrays only.")
937 return compress_rowcols(a, 0)
940def compress_cols(a):
941 """
942 Suppress whole columns of a 2-D array that contain masked values.
944 This is equivalent to ``np.ma.compress_rowcols(a, 1)``, see
945 `compress_rowcols` for details.
947 See Also
948 --------
949 compress_rowcols
951 """
952 a = asarray(a)
953 if a.ndim != 2:
954 raise NotImplementedError("compress_cols works for 2D arrays only.")
955 return compress_rowcols(a, 1)
958def mask_rows(a, axis=np._NoValue):
959 """
960 Mask rows of a 2D array that contain masked values.
962 This function is a shortcut to ``mask_rowcols`` with `axis` equal to 0.
964 See Also
965 --------
966 mask_rowcols : Mask rows and/or columns of a 2D array.
967 masked_where : Mask where a condition is met.
969 Examples
970 --------
971 >>> import numpy.ma as ma
972 >>> a = np.zeros((3, 3), dtype=int)
973 >>> a[1, 1] = 1
974 >>> a
975 array([[0, 0, 0],
976 [0, 1, 0],
977 [0, 0, 0]])
978 >>> a = ma.masked_equal(a, 1)
979 >>> a
980 masked_array(
981 data=[[0, 0, 0],
982 [0, --, 0],
983 [0, 0, 0]],
984 mask=[[False, False, False],
985 [False, True, False],
986 [False, False, False]],
987 fill_value=1)
989 >>> ma.mask_rows(a)
990 masked_array(
991 data=[[0, 0, 0],
992 [--, --, --],
993 [0, 0, 0]],
994 mask=[[False, False, False],
995 [ True, True, True],
996 [False, False, False]],
997 fill_value=1)
999 """
1000 if axis is not np._NoValue:
1001 # remove the axis argument when this deprecation expires
1002 # NumPy 1.18.0, 2019-11-28
1003 warnings.warn(
1004 "The axis argument has always been ignored, in future passing it "
1005 "will raise TypeError", DeprecationWarning, stacklevel=2)
1006 return mask_rowcols(a, 0)
1009def mask_cols(a, axis=np._NoValue):
1010 """
1011 Mask columns of a 2D array that contain masked values.
1013 This function is a shortcut to ``mask_rowcols`` with `axis` equal to 1.
1015 See Also
1016 --------
1017 mask_rowcols : Mask rows and/or columns of a 2D array.
1018 masked_where : Mask where a condition is met.
1020 Examples
1021 --------
1022 >>> import numpy.ma as ma
1023 >>> a = np.zeros((3, 3), dtype=int)
1024 >>> a[1, 1] = 1
1025 >>> a
1026 array([[0, 0, 0],
1027 [0, 1, 0],
1028 [0, 0, 0]])
1029 >>> a = ma.masked_equal(a, 1)
1030 >>> a
1031 masked_array(
1032 data=[[0, 0, 0],
1033 [0, --, 0],
1034 [0, 0, 0]],
1035 mask=[[False, False, False],
1036 [False, True, False],
1037 [False, False, False]],
1038 fill_value=1)
1039 >>> ma.mask_cols(a)
1040 masked_array(
1041 data=[[0, --, 0],
1042 [0, --, 0],
1043 [0, --, 0]],
1044 mask=[[False, True, False],
1045 [False, True, False],
1046 [False, True, False]],
1047 fill_value=1)
1049 """
1050 if axis is not np._NoValue:
1051 # remove the axis argument when this deprecation expires
1052 # NumPy 1.18.0, 2019-11-28
1053 warnings.warn(
1054 "The axis argument has always been ignored, in future passing it "
1055 "will raise TypeError", DeprecationWarning, stacklevel=2)
1056 return mask_rowcols(a, 1)
1059#####--------------------------------------------------------------------------
1060#---- --- arraysetops ---
1061#####--------------------------------------------------------------------------
1063def ediff1d(arr, to_end=None, to_begin=None):
1064 """
1065 Compute the differences between consecutive elements of an array.
1067 This function is the equivalent of `numpy.ediff1d` that takes masked
1068 values into account, see `numpy.ediff1d` for details.
1070 See Also
1071 --------
1072 numpy.ediff1d : Equivalent function for ndarrays.
1074 """
1075 arr = ma.asanyarray(arr).flat
1076 ed = arr[1:] - arr[:-1]
1077 arrays = [ed]
1078 #
1079 if to_begin is not None:
1080 arrays.insert(0, to_begin)
1081 if to_end is not None:
1082 arrays.append(to_end)
1083 #
1084 if len(arrays) != 1:
1085 # We'll save ourselves a copy of a potentially large array in the common
1086 # case where neither to_begin or to_end was given.
1087 ed = hstack(arrays)
1088 #
1089 return ed
1092def unique(ar1, return_index=False, return_inverse=False):
1093 """
1094 Finds the unique elements of an array.
1096 Masked values are considered the same element (masked). The output array
1097 is always a masked array. See `numpy.unique` for more details.
1099 See Also
1100 --------
1101 numpy.unique : Equivalent function for ndarrays.
1103 Examples
1104 --------
1105 >>> import numpy.ma as ma
1106 >>> a = [1, 2, 1000, 2, 3]
1107 >>> mask = [0, 0, 1, 0, 0]
1108 >>> masked_a = ma.masked_array(a, mask)
1109 >>> masked_a
1110 masked_array(data=[1, 2, --, 2, 3],
1111 mask=[False, False, True, False, False],
1112 fill_value=999999)
1113 >>> ma.unique(masked_a)
1114 masked_array(data=[1, 2, 3, --],
1115 mask=[False, False, False, True],
1116 fill_value=999999)
1117 >>> ma.unique(masked_a, return_index=True)
1118 (masked_array(data=[1, 2, 3, --],
1119 mask=[False, False, False, True],
1120 fill_value=999999), array([0, 1, 4, 2]))
1121 >>> ma.unique(masked_a, return_inverse=True)
1122 (masked_array(data=[1, 2, 3, --],
1123 mask=[False, False, False, True],
1124 fill_value=999999), array([0, 1, 3, 1, 2]))
1125 >>> ma.unique(masked_a, return_index=True, return_inverse=True)
1126 (masked_array(data=[1, 2, 3, --],
1127 mask=[False, False, False, True],
1128 fill_value=999999), array([0, 1, 4, 2]), array([0, 1, 3, 1, 2]))
1129 """
1130 output = np.unique(ar1,
1131 return_index=return_index,
1132 return_inverse=return_inverse)
1133 if isinstance(output, tuple):
1134 output = list(output)
1135 output[0] = output[0].view(MaskedArray)
1136 output = tuple(output)
1137 else:
1138 output = output.view(MaskedArray)
1139 return output
1142def intersect1d(ar1, ar2, assume_unique=False):
1143 """
1144 Returns the unique elements common to both arrays.
1146 Masked values are considered equal one to the other.
1147 The output is always a masked array.
1149 See `numpy.intersect1d` for more details.
1151 See Also
1152 --------
1153 numpy.intersect1d : Equivalent function for ndarrays.
1155 Examples
1156 --------
1157 >>> x = np.ma.array([1, 3, 3, 3], mask=[0, 0, 0, 1])
1158 >>> y = np.ma.array([3, 1, 1, 1], mask=[0, 0, 0, 1])
1159 >>> np.ma.intersect1d(x, y)
1160 masked_array(data=[1, 3, --],
1161 mask=[False, False, True],
1162 fill_value=999999)
1164 """
1165 if assume_unique:
1166 aux = ma.concatenate((ar1, ar2))
1167 else:
1168 # Might be faster than unique( intersect1d( ar1, ar2 ) )?
1169 aux = ma.concatenate((unique(ar1), unique(ar2)))
1170 aux.sort()
1171 return aux[:-1][aux[1:] == aux[:-1]]
1174def setxor1d(ar1, ar2, assume_unique=False):
1175 """
1176 Set exclusive-or of 1-D arrays with unique elements.
1178 The output is always a masked array. See `numpy.setxor1d` for more details.
1180 See Also
1181 --------
1182 numpy.setxor1d : Equivalent function for ndarrays.
1184 """
1185 if not assume_unique:
1186 ar1 = unique(ar1)
1187 ar2 = unique(ar2)
1189 aux = ma.concatenate((ar1, ar2))
1190 if aux.size == 0:
1191 return aux
1192 aux.sort()
1193 auxf = aux.filled()
1194# flag = ediff1d( aux, to_end = 1, to_begin = 1 ) == 0
1195 flag = ma.concatenate(([True], (auxf[1:] != auxf[:-1]), [True]))
1196# flag2 = ediff1d( flag ) == 0
1197 flag2 = (flag[1:] == flag[:-1])
1198 return aux[flag2]
1201def in1d(ar1, ar2, assume_unique=False, invert=False):
1202 """
1203 Test whether each element of an array is also present in a second
1204 array.
1206 The output is always a masked array. See `numpy.in1d` for more details.
1208 We recommend using :func:`isin` instead of `in1d` for new code.
1210 See Also
1211 --------
1212 isin : Version of this function that preserves the shape of ar1.
1213 numpy.in1d : Equivalent function for ndarrays.
1215 Notes
1216 -----
1217 .. versionadded:: 1.4.0
1219 """
1220 if not assume_unique:
1221 ar1, rev_idx = unique(ar1, return_inverse=True)
1222 ar2 = unique(ar2)
1224 ar = ma.concatenate((ar1, ar2))
1225 # We need this to be a stable sort, so always use 'mergesort'
1226 # here. The values from the first array should always come before
1227 # the values from the second array.
1228 order = ar.argsort(kind='mergesort')
1229 sar = ar[order]
1230 if invert:
1231 bool_ar = (sar[1:] != sar[:-1])
1232 else:
1233 bool_ar = (sar[1:] == sar[:-1])
1234 flag = ma.concatenate((bool_ar, [invert]))
1235 indx = order.argsort(kind='mergesort')[:len(ar1)]
1237 if assume_unique:
1238 return flag[indx]
1239 else:
1240 return flag[indx][rev_idx]
1243def isin(element, test_elements, assume_unique=False, invert=False):
1244 """
1245 Calculates `element in test_elements`, broadcasting over
1246 `element` only.
1248 The output is always a masked array of the same shape as `element`.
1249 See `numpy.isin` for more details.
1251 See Also
1252 --------
1253 in1d : Flattened version of this function.
1254 numpy.isin : Equivalent function for ndarrays.
1256 Notes
1257 -----
1258 .. versionadded:: 1.13.0
1260 """
1261 element = ma.asarray(element)
1262 return in1d(element, test_elements, assume_unique=assume_unique,
1263 invert=invert).reshape(element.shape)
1266def union1d(ar1, ar2):
1267 """
1268 Union of two arrays.
1270 The output is always a masked array. See `numpy.union1d` for more details.
1272 See Also
1273 --------
1274 numpy.union1d : Equivalent function for ndarrays.
1276 """
1277 return unique(ma.concatenate((ar1, ar2), axis=None))
1280def setdiff1d(ar1, ar2, assume_unique=False):
1281 """
1282 Set difference of 1D arrays with unique elements.
1284 The output is always a masked array. See `numpy.setdiff1d` for more
1285 details.
1287 See Also
1288 --------
1289 numpy.setdiff1d : Equivalent function for ndarrays.
1291 Examples
1292 --------
1293 >>> x = np.ma.array([1, 2, 3, 4], mask=[0, 1, 0, 1])
1294 >>> np.ma.setdiff1d(x, [1, 2])
1295 masked_array(data=[3, --],
1296 mask=[False, True],
1297 fill_value=999999)
1299 """
1300 if assume_unique:
1301 ar1 = ma.asarray(ar1).ravel()
1302 else:
1303 ar1 = unique(ar1)
1304 ar2 = unique(ar2)
1305 return ar1[in1d(ar1, ar2, assume_unique=True, invert=True)]
1308###############################################################################
1309# Covariance #
1310###############################################################################
1313def _covhelper(x, y=None, rowvar=True, allow_masked=True):
1314 """
1315 Private function for the computation of covariance and correlation
1316 coefficients.
1318 """
1319 x = ma.array(x, ndmin=2, copy=True, dtype=float)
1320 xmask = ma.getmaskarray(x)
1321 # Quick exit if we can't process masked data
1322 if not allow_masked and xmask.any():
1323 raise ValueError("Cannot process masked data.")
1324 #
1325 if x.shape[0] == 1:
1326 rowvar = True
1327 # Make sure that rowvar is either 0 or 1
1328 rowvar = int(bool(rowvar))
1329 axis = 1 - rowvar
1330 if rowvar:
1331 tup = (slice(None), None)
1332 else:
1333 tup = (None, slice(None))
1334 #
1335 if y is None:
1336 xnotmask = np.logical_not(xmask).astype(int)
1337 else:
1338 y = array(y, copy=False, ndmin=2, dtype=float)
1339 ymask = ma.getmaskarray(y)
1340 if not allow_masked and ymask.any():
1341 raise ValueError("Cannot process masked data.")
1342 if xmask.any() or ymask.any():
1343 if y.shape == x.shape:
1344 # Define some common mask
1345 common_mask = np.logical_or(xmask, ymask)
1346 if common_mask is not nomask:
1347 xmask = x._mask = y._mask = ymask = common_mask
1348 x._sharedmask = False
1349 y._sharedmask = False
1350 x = ma.concatenate((x, y), axis)
1351 xnotmask = np.logical_not(np.concatenate((xmask, ymask), axis)).astype(int)
1352 x -= x.mean(axis=rowvar)[tup]
1353 return (x, xnotmask, rowvar)
1356def cov(x, y=None, rowvar=True, bias=False, allow_masked=True, ddof=None):
1357 """
1358 Estimate the covariance matrix.
1360 Except for the handling of missing data this function does the same as
1361 `numpy.cov`. For more details and examples, see `numpy.cov`.
1363 By default, masked values are recognized as such. If `x` and `y` have the
1364 same shape, a common mask is allocated: if ``x[i,j]`` is masked, then
1365 ``y[i,j]`` will also be masked.
1366 Setting `allow_masked` to False will raise an exception if values are
1367 missing in either of the input arrays.
1369 Parameters
1370 ----------
1371 x : array_like
1372 A 1-D or 2-D array containing multiple variables and observations.
1373 Each row of `x` represents a variable, and each column a single
1374 observation of all those variables. Also see `rowvar` below.
1375 y : array_like, optional
1376 An additional set of variables and observations. `y` has the same
1377 shape as `x`.
1378 rowvar : bool, optional
1379 If `rowvar` is True (default), then each row represents a
1380 variable, with observations in the columns. Otherwise, the relationship
1381 is transposed: each column represents a variable, while the rows
1382 contain observations.
1383 bias : bool, optional
1384 Default normalization (False) is by ``(N-1)``, where ``N`` is the
1385 number of observations given (unbiased estimate). If `bias` is True,
1386 then normalization is by ``N``. This keyword can be overridden by
1387 the keyword ``ddof`` in numpy versions >= 1.5.
1388 allow_masked : bool, optional
1389 If True, masked values are propagated pair-wise: if a value is masked
1390 in `x`, the corresponding value is masked in `y`.
1391 If False, raises a `ValueError` exception when some values are missing.
1392 ddof : {None, int}, optional
1393 If not ``None`` normalization is by ``(N - ddof)``, where ``N`` is
1394 the number of observations; this overrides the value implied by
1395 ``bias``. The default value is ``None``.
1397 .. versionadded:: 1.5
1399 Raises
1400 ------
1401 ValueError
1402 Raised if some values are missing and `allow_masked` is False.
1404 See Also
1405 --------
1406 numpy.cov
1408 """
1409 # Check inputs
1410 if ddof is not None and ddof != int(ddof):
1411 raise ValueError("ddof must be an integer")
1412 # Set up ddof
1413 if ddof is None:
1414 if bias:
1415 ddof = 0
1416 else:
1417 ddof = 1
1419 (x, xnotmask, rowvar) = _covhelper(x, y, rowvar, allow_masked)
1420 if not rowvar:
1421 fact = np.dot(xnotmask.T, xnotmask) * 1. - ddof
1422 result = (dot(x.T, x.conj(), strict=False) / fact).squeeze()
1423 else:
1424 fact = np.dot(xnotmask, xnotmask.T) * 1. - ddof
1425 result = (dot(x, x.T.conj(), strict=False) / fact).squeeze()
1426 return result
1429def corrcoef(x, y=None, rowvar=True, bias=np._NoValue, allow_masked=True,
1430 ddof=np._NoValue):
1431 """
1432 Return Pearson product-moment correlation coefficients.
1434 Except for the handling of missing data this function does the same as
1435 `numpy.corrcoef`. For more details and examples, see `numpy.corrcoef`.
1437 Parameters
1438 ----------
1439 x : array_like
1440 A 1-D or 2-D array containing multiple variables and observations.
1441 Each row of `x` represents a variable, and each column a single
1442 observation of all those variables. Also see `rowvar` below.
1443 y : array_like, optional
1444 An additional set of variables and observations. `y` has the same
1445 shape as `x`.
1446 rowvar : bool, optional
1447 If `rowvar` is True (default), then each row represents a
1448 variable, with observations in the columns. Otherwise, the relationship
1449 is transposed: each column represents a variable, while the rows
1450 contain observations.
1451 bias : _NoValue, optional
1452 Has no effect, do not use.
1454 .. deprecated:: 1.10.0
1455 allow_masked : bool, optional
1456 If True, masked values are propagated pair-wise: if a value is masked
1457 in `x`, the corresponding value is masked in `y`.
1458 If False, raises an exception. Because `bias` is deprecated, this
1459 argument needs to be treated as keyword only to avoid a warning.
1460 ddof : _NoValue, optional
1461 Has no effect, do not use.
1463 .. deprecated:: 1.10.0
1465 See Also
1466 --------
1467 numpy.corrcoef : Equivalent function in top-level NumPy module.
1468 cov : Estimate the covariance matrix.
1470 Notes
1471 -----
1472 This function accepts but discards arguments `bias` and `ddof`. This is
1473 for backwards compatibility with previous versions of this function. These
1474 arguments had no effect on the return values of the function and can be
1475 safely ignored in this and previous versions of numpy.
1476 """
1477 msg = 'bias and ddof have no effect and are deprecated'
1478 if bias is not np._NoValue or ddof is not np._NoValue:
1479 # 2015-03-15, 1.10
1480 warnings.warn(msg, DeprecationWarning, stacklevel=2)
1481 # Get the data
1482 (x, xnotmask, rowvar) = _covhelper(x, y, rowvar, allow_masked)
1483 # Compute the covariance matrix
1484 if not rowvar:
1485 fact = np.dot(xnotmask.T, xnotmask) * 1.
1486 c = (dot(x.T, x.conj(), strict=False) / fact).squeeze()
1487 else:
1488 fact = np.dot(xnotmask, xnotmask.T) * 1.
1489 c = (dot(x, x.T.conj(), strict=False) / fact).squeeze()
1490 # Check whether we have a scalar
1491 try:
1492 diag = ma.diagonal(c)
1493 except ValueError:
1494 return 1
1495 #
1496 if xnotmask.all():
1497 _denom = ma.sqrt(ma.multiply.outer(diag, diag))
1498 else:
1499 _denom = diagflat(diag)
1500 _denom._sharedmask = False # We know return is always a copy
1501 n = x.shape[1 - rowvar]
1502 if rowvar:
1503 for i in range(n - 1):
1504 for j in range(i + 1, n):
1505 _x = mask_cols(vstack((x[i], x[j]))).var(axis=1)
1506 _denom[i, j] = _denom[j, i] = ma.sqrt(ma.multiply.reduce(_x))
1507 else:
1508 for i in range(n - 1):
1509 for j in range(i + 1, n):
1510 _x = mask_cols(
1511 vstack((x[:, i], x[:, j]))).var(axis=1)
1512 _denom[i, j] = _denom[j, i] = ma.sqrt(ma.multiply.reduce(_x))
1513 return c / _denom
1515#####--------------------------------------------------------------------------
1516#---- --- Concatenation helpers ---
1517#####--------------------------------------------------------------------------
1519class MAxisConcatenator(AxisConcatenator):
1520 """
1521 Translate slice objects to concatenation along an axis.
1523 For documentation on usage, see `mr_class`.
1525 See Also
1526 --------
1527 mr_class
1529 """
1530 concatenate = staticmethod(concatenate)
1532 @classmethod
1533 def makemat(cls, arr):
1534 # There used to be a view as np.matrix here, but we may eventually
1535 # deprecate that class. In preparation, we use the unmasked version
1536 # to construct the matrix (with copy=False for backwards compatibility
1537 # with the .view)
1538 data = super().makemat(arr.data, copy=False)
1539 return array(data, mask=arr.mask)
1541 def __getitem__(self, key):
1542 # matrix builder syntax, like 'a, b; c, d'
1543 if isinstance(key, str):
1544 raise MAError("Unavailable for masked array.")
1546 return super().__getitem__(key)
1549class mr_class(MAxisConcatenator):
1550 """
1551 Translate slice objects to concatenation along the first axis.
1553 This is the masked array version of `lib.index_tricks.RClass`.
1555 See Also
1556 --------
1557 lib.index_tricks.RClass
1559 Examples
1560 --------
1561 >>> np.ma.mr_[np.ma.array([1,2,3]), 0, 0, np.ma.array([4,5,6])]
1562 masked_array(data=[1, 2, 3, ..., 4, 5, 6],
1563 mask=False,
1564 fill_value=999999)
1566 """
1567 def __init__(self):
1568 MAxisConcatenator.__init__(self, 0)
1570mr_ = mr_class()
1573#####--------------------------------------------------------------------------
1574#---- Find unmasked data ---
1575#####--------------------------------------------------------------------------
1577def ndenumerate(a, compressed=True):
1578 """
1579 Multidimensional index iterator.
1581 Return an iterator yielding pairs of array coordinates and values,
1582 skipping elements that are masked. With `compressed=False`,
1583 `ma.masked` is yielded as the value of masked elements. This
1584 behavior differs from that of `numpy.ndenumerate`, which yields the
1585 value of the underlying data array.
1587 Notes
1588 -----
1589 .. versionadded:: 1.23.0
1591 Parameters
1592 ----------
1593 a : array_like
1594 An array with (possibly) masked elements.
1595 compressed : bool, optional
1596 If True (default), masked elements are skipped.
1598 See Also
1599 --------
1600 numpy.ndenumerate : Equivalent function ignoring any mask.
1602 Examples
1603 --------
1604 >>> a = np.ma.arange(9).reshape((3, 3))
1605 >>> a[1, 0] = np.ma.masked
1606 >>> a[1, 2] = np.ma.masked
1607 >>> a[2, 1] = np.ma.masked
1608 >>> a
1609 masked_array(
1610 data=[[0, 1, 2],
1611 [--, 4, --],
1612 [6, --, 8]],
1613 mask=[[False, False, False],
1614 [ True, False, True],
1615 [False, True, False]],
1616 fill_value=999999)
1617 >>> for index, x in np.ma.ndenumerate(a):
1618 ... print(index, x)
1619 (0, 0) 0
1620 (0, 1) 1
1621 (0, 2) 2
1622 (1, 1) 4
1623 (2, 0) 6
1624 (2, 2) 8
1626 >>> for index, x in np.ma.ndenumerate(a, compressed=False):
1627 ... print(index, x)
1628 (0, 0) 0
1629 (0, 1) 1
1630 (0, 2) 2
1631 (1, 0) --
1632 (1, 1) 4
1633 (1, 2) --
1634 (2, 0) 6
1635 (2, 1) --
1636 (2, 2) 8
1637 """
1638 for it, mask in zip(np.ndenumerate(a), getmaskarray(a).flat):
1639 if not mask:
1640 yield it
1641 elif not compressed:
1642 yield it[0], masked
1645def flatnotmasked_edges(a):
1646 """
1647 Find the indices of the first and last unmasked values.
1649 Expects a 1-D `MaskedArray`, returns None if all values are masked.
1651 Parameters
1652 ----------
1653 a : array_like
1654 Input 1-D `MaskedArray`
1656 Returns
1657 -------
1658 edges : ndarray or None
1659 The indices of first and last non-masked value in the array.
1660 Returns None if all values are masked.
1662 See Also
1663 --------
1664 flatnotmasked_contiguous, notmasked_contiguous, notmasked_edges
1665 clump_masked, clump_unmasked
1667 Notes
1668 -----
1669 Only accepts 1-D arrays.
1671 Examples
1672 --------
1673 >>> a = np.ma.arange(10)
1674 >>> np.ma.flatnotmasked_edges(a)
1675 array([0, 9])
1677 >>> mask = (a < 3) | (a > 8) | (a == 5)
1678 >>> a[mask] = np.ma.masked
1679 >>> np.array(a[~a.mask])
1680 array([3, 4, 6, 7, 8])
1682 >>> np.ma.flatnotmasked_edges(a)
1683 array([3, 8])
1685 >>> a[:] = np.ma.masked
1686 >>> print(np.ma.flatnotmasked_edges(a))
1687 None
1689 """
1690 m = getmask(a)
1691 if m is nomask or not np.any(m):
1692 return np.array([0, a.size - 1])
1693 unmasked = np.flatnonzero(~m)
1694 if len(unmasked) > 0:
1695 return unmasked[[0, -1]]
1696 else:
1697 return None
1700def notmasked_edges(a, axis=None):
1701 """
1702 Find the indices of the first and last unmasked values along an axis.
1704 If all values are masked, return None. Otherwise, return a list
1705 of two tuples, corresponding to the indices of the first and last
1706 unmasked values respectively.
1708 Parameters
1709 ----------
1710 a : array_like
1711 The input array.
1712 axis : int, optional
1713 Axis along which to perform the operation.
1714 If None (default), applies to a flattened version of the array.
1716 Returns
1717 -------
1718 edges : ndarray or list
1719 An array of start and end indexes if there are any masked data in
1720 the array. If there are no masked data in the array, `edges` is a
1721 list of the first and last index.
1723 See Also
1724 --------
1725 flatnotmasked_contiguous, flatnotmasked_edges, notmasked_contiguous
1726 clump_masked, clump_unmasked
1728 Examples
1729 --------
1730 >>> a = np.arange(9).reshape((3, 3))
1731 >>> m = np.zeros_like(a)
1732 >>> m[1:, 1:] = 1
1734 >>> am = np.ma.array(a, mask=m)
1735 >>> np.array(am[~am.mask])
1736 array([0, 1, 2, 3, 6])
1738 >>> np.ma.notmasked_edges(am)
1739 array([0, 6])
1741 """
1742 a = asarray(a)
1743 if axis is None or a.ndim == 1:
1744 return flatnotmasked_edges(a)
1745 m = getmaskarray(a)
1746 idx = array(np.indices(a.shape), mask=np.asarray([m] * a.ndim))
1747 return [tuple([idx[i].min(axis).compressed() for i in range(a.ndim)]),
1748 tuple([idx[i].max(axis).compressed() for i in range(a.ndim)]), ]
1751def flatnotmasked_contiguous(a):
1752 """
1753 Find contiguous unmasked data in a masked array.
1755 Parameters
1756 ----------
1757 a : array_like
1758 The input array.
1760 Returns
1761 -------
1762 slice_list : list
1763 A sorted sequence of `slice` objects (start index, end index).
1765 .. versionchanged:: 1.15.0
1766 Now returns an empty list instead of None for a fully masked array
1768 See Also
1769 --------
1770 flatnotmasked_edges, notmasked_contiguous, notmasked_edges
1771 clump_masked, clump_unmasked
1773 Notes
1774 -----
1775 Only accepts 2-D arrays at most.
1777 Examples
1778 --------
1779 >>> a = np.ma.arange(10)
1780 >>> np.ma.flatnotmasked_contiguous(a)
1781 [slice(0, 10, None)]
1783 >>> mask = (a < 3) | (a > 8) | (a == 5)
1784 >>> a[mask] = np.ma.masked
1785 >>> np.array(a[~a.mask])
1786 array([3, 4, 6, 7, 8])
1788 >>> np.ma.flatnotmasked_contiguous(a)
1789 [slice(3, 5, None), slice(6, 9, None)]
1790 >>> a[:] = np.ma.masked
1791 >>> np.ma.flatnotmasked_contiguous(a)
1792 []
1794 """
1795 m = getmask(a)
1796 if m is nomask:
1797 return [slice(0, a.size)]
1798 i = 0
1799 result = []
1800 for (k, g) in itertools.groupby(m.ravel()):
1801 n = len(list(g))
1802 if not k:
1803 result.append(slice(i, i + n))
1804 i += n
1805 return result
1808def notmasked_contiguous(a, axis=None):
1809 """
1810 Find contiguous unmasked data in a masked array along the given axis.
1812 Parameters
1813 ----------
1814 a : array_like
1815 The input array.
1816 axis : int, optional
1817 Axis along which to perform the operation.
1818 If None (default), applies to a flattened version of the array, and this
1819 is the same as `flatnotmasked_contiguous`.
1821 Returns
1822 -------
1823 endpoints : list
1824 A list of slices (start and end indexes) of unmasked indexes
1825 in the array.
1827 If the input is 2d and axis is specified, the result is a list of lists.
1829 See Also
1830 --------
1831 flatnotmasked_edges, flatnotmasked_contiguous, notmasked_edges
1832 clump_masked, clump_unmasked
1834 Notes
1835 -----
1836 Only accepts 2-D arrays at most.
1838 Examples
1839 --------
1840 >>> a = np.arange(12).reshape((3, 4))
1841 >>> mask = np.zeros_like(a)
1842 >>> mask[1:, :-1] = 1; mask[0, 1] = 1; mask[-1, 0] = 0
1843 >>> ma = np.ma.array(a, mask=mask)
1844 >>> ma
1845 masked_array(
1846 data=[[0, --, 2, 3],
1847 [--, --, --, 7],
1848 [8, --, --, 11]],
1849 mask=[[False, True, False, False],
1850 [ True, True, True, False],
1851 [False, True, True, False]],
1852 fill_value=999999)
1853 >>> np.array(ma[~ma.mask])
1854 array([ 0, 2, 3, 7, 8, 11])
1856 >>> np.ma.notmasked_contiguous(ma)
1857 [slice(0, 1, None), slice(2, 4, None), slice(7, 9, None), slice(11, 12, None)]
1859 >>> np.ma.notmasked_contiguous(ma, axis=0)
1860 [[slice(0, 1, None), slice(2, 3, None)], [], [slice(0, 1, None)], [slice(0, 3, None)]]
1862 >>> np.ma.notmasked_contiguous(ma, axis=1)
1863 [[slice(0, 1, None), slice(2, 4, None)], [slice(3, 4, None)], [slice(0, 1, None), slice(3, 4, None)]]
1865 """
1866 a = asarray(a)
1867 nd = a.ndim
1868 if nd > 2:
1869 raise NotImplementedError("Currently limited to at most 2D array.")
1870 if axis is None or nd == 1:
1871 return flatnotmasked_contiguous(a)
1872 #
1873 result = []
1874 #
1875 other = (axis + 1) % 2
1876 idx = [0, 0]
1877 idx[axis] = slice(None, None)
1878 #
1879 for i in range(a.shape[other]):
1880 idx[other] = i
1881 result.append(flatnotmasked_contiguous(a[tuple(idx)]))
1882 return result
1885def _ezclump(mask):
1886 """
1887 Finds the clumps (groups of data with the same values) for a 1D bool array.
1889 Returns a series of slices.
1890 """
1891 if mask.ndim > 1:
1892 mask = mask.ravel()
1893 idx = (mask[1:] ^ mask[:-1]).nonzero()
1894 idx = idx[0] + 1
1896 if mask[0]:
1897 if len(idx) == 0:
1898 return [slice(0, mask.size)]
1900 r = [slice(0, idx[0])]
1901 r.extend((slice(left, right)
1902 for left, right in zip(idx[1:-1:2], idx[2::2])))
1903 else:
1904 if len(idx) == 0:
1905 return []
1907 r = [slice(left, right) for left, right in zip(idx[:-1:2], idx[1::2])]
1909 if mask[-1]:
1910 r.append(slice(idx[-1], mask.size))
1911 return r
1914def clump_unmasked(a):
1915 """
1916 Return list of slices corresponding to the unmasked clumps of a 1-D array.
1917 (A "clump" is defined as a contiguous region of the array).
1919 Parameters
1920 ----------
1921 a : ndarray
1922 A one-dimensional masked array.
1924 Returns
1925 -------
1926 slices : list of slice
1927 The list of slices, one for each continuous region of unmasked
1928 elements in `a`.
1930 Notes
1931 -----
1932 .. versionadded:: 1.4.0
1934 See Also
1935 --------
1936 flatnotmasked_edges, flatnotmasked_contiguous, notmasked_edges
1937 notmasked_contiguous, clump_masked
1939 Examples
1940 --------
1941 >>> a = np.ma.masked_array(np.arange(10))
1942 >>> a[[0, 1, 2, 6, 8, 9]] = np.ma.masked
1943 >>> np.ma.clump_unmasked(a)
1944 [slice(3, 6, None), slice(7, 8, None)]
1946 """
1947 mask = getattr(a, '_mask', nomask)
1948 if mask is nomask:
1949 return [slice(0, a.size)]
1950 return _ezclump(~mask)
1953def clump_masked(a):
1954 """
1955 Returns a list of slices corresponding to the masked clumps of a 1-D array.
1956 (A "clump" is defined as a contiguous region of the array).
1958 Parameters
1959 ----------
1960 a : ndarray
1961 A one-dimensional masked array.
1963 Returns
1964 -------
1965 slices : list of slice
1966 The list of slices, one for each continuous region of masked elements
1967 in `a`.
1969 Notes
1970 -----
1971 .. versionadded:: 1.4.0
1973 See Also
1974 --------
1975 flatnotmasked_edges, flatnotmasked_contiguous, notmasked_edges
1976 notmasked_contiguous, clump_unmasked
1978 Examples
1979 --------
1980 >>> a = np.ma.masked_array(np.arange(10))
1981 >>> a[[0, 1, 2, 6, 8, 9]] = np.ma.masked
1982 >>> np.ma.clump_masked(a)
1983 [slice(0, 3, None), slice(6, 7, None), slice(8, 10, None)]
1985 """
1986 mask = ma.getmask(a)
1987 if mask is nomask:
1988 return []
1989 return _ezclump(mask)
1992###############################################################################
1993# Polynomial fit #
1994###############################################################################
1997def vander(x, n=None):
1998 """
1999 Masked values in the input array result in rows of zeros.
2001 """
2002 _vander = np.vander(x, n)
2003 m = getmask(x)
2004 if m is not nomask:
2005 _vander[m] = 0
2006 return _vander
2008vander.__doc__ = ma.doc_note(np.vander.__doc__, vander.__doc__)
2011def polyfit(x, y, deg, rcond=None, full=False, w=None, cov=False):
2012 """
2013 Any masked values in x is propagated in y, and vice-versa.
2015 """
2016 x = asarray(x)
2017 y = asarray(y)
2019 m = getmask(x)
2020 if y.ndim == 1:
2021 m = mask_or(m, getmask(y))
2022 elif y.ndim == 2:
2023 my = getmask(mask_rows(y))
2024 if my is not nomask:
2025 m = mask_or(m, my[:, 0])
2026 else:
2027 raise TypeError("Expected a 1D or 2D array for y!")
2029 if w is not None:
2030 w = asarray(w)
2031 if w.ndim != 1:
2032 raise TypeError("expected a 1-d array for weights")
2033 if w.shape[0] != y.shape[0]:
2034 raise TypeError("expected w and y to have the same length")
2035 m = mask_or(m, getmask(w))
2037 if m is not nomask:
2038 not_m = ~m
2039 if w is not None:
2040 w = w[not_m]
2041 return np.polyfit(x[not_m], y[not_m], deg, rcond, full, w, cov)
2042 else:
2043 return np.polyfit(x, y, deg, rcond, full, w, cov)
2045polyfit.__doc__ = ma.doc_note(np.polyfit.__doc__, polyfit.__doc__)