diff options
Diffstat (limited to 'numpy/lib/function_base.py')
| -rw-r--r-- | numpy/lib/function_base.py | 192 |
1 files changed, 120 insertions, 72 deletions
diff --git a/numpy/lib/function_base.py b/numpy/lib/function_base.py index 651343a2c..783d45c2f 100644 --- a/numpy/lib/function_base.py +++ b/numpy/lib/function_base.py @@ -8,7 +8,7 @@ import numpy as np import numpy.core.numeric as _nx from numpy.core import transpose from numpy.core.numeric import ( - ones, zeros, arange, concatenate, array, asarray, asanyarray, empty, + ones, zeros_like, arange, concatenate, array, asarray, asanyarray, empty, ndarray, around, floor, ceil, take, dot, where, intp, integer, isscalar, absolute ) @@ -431,10 +431,13 @@ def asarray_chkfinite(a, dtype=None, order=None): of lists and ndarrays. Success requires no NaNs or Infs. dtype : data-type, optional By default, the data-type is inferred from the input data. - order : {'C', 'F'}, optional - Whether to use row-major (C-style) or - column-major (Fortran-style) memory representation. - Defaults to 'C'. + order : {'C', 'F', 'A', 'K'}, optional + Memory layout. 'A' and 'K' depend on the order of input array a. + 'C' row-major (C-style), + 'F' column-major (Fortran-style) memory representation. + 'A' (any) means 'F' if `a` is Fortran contiguous, 'C' otherwise + 'K' (keep) preserve input order + Defaults to 'C'. Returns ------- @@ -590,7 +593,7 @@ def piecewise(x, condlist, funclist, *args, **kw): not isinstance(condlist[0], (list, ndarray)) and x.ndim != 0): condlist = [condlist] - condlist = array(condlist, dtype=bool) + condlist = asarray(condlist, dtype=bool) n = len(condlist) if n == n2 - 1: # compute the "otherwise" condition. @@ -603,7 +606,7 @@ def piecewise(x, condlist, funclist, *args, **kw): .format(n, n, n+1) ) - y = zeros(x.shape, x.dtype) + y = zeros_like(x) for cond, func in zip(condlist, funclist): if not isinstance(func, collections.abc.Callable): y[cond] = func @@ -668,11 +671,22 @@ def select(condlist, choicelist, default=0): raise ValueError("select with an empty condition list is not possible") choicelist = [np.asarray(choice) for choice in choicelist] - choicelist.append(np.asarray(default)) + + try: + intermediate_dtype = np.result_type(*choicelist) + except TypeError as e: + msg = f'Choicelist elements do not have a common dtype: {e}' + raise TypeError(msg) from None + default_array = np.asarray(default) + choicelist.append(default_array) # need to get the result type before broadcasting for correct scalar # behaviour - dtype = np.result_type(*choicelist) + try: + dtype = np.result_type(intermediate_dtype, default_array) + except TypeError as e: + msg = f'Choicelists and default value do not have a common dtype: {e}' + raise TypeError(msg) from None # Convert conditions to arrays and broadcast conditions and choices # as the shape is needed for the result. Doing it separately optimizes @@ -843,7 +857,7 @@ def gradient(f, *varargs, axis=None, edge_order=1): Returns ------- gradient : ndarray or list of ndarray - A set of ndarrays (or a single ndarray if there is only one dimension) + A list of ndarrays (or a single ndarray if there is only one dimension) corresponding to the derivatives of f with respect to each dimension. Each derivative has the same shape as f. @@ -1287,7 +1301,7 @@ def _interp_dispatcher(x, xp, fp, left=None, right=None, period=None): @array_function_dispatch(_interp_dispatcher) def interp(x, xp, fp, left=None, right=None, period=None): """ - One-dimensional linear interpolation. + One-dimensional linear interpolation for monotonically increasing sample points. Returns the one-dimensional piecewise linear interpolant to a function with given discrete data points (`xp`, `fp`), evaluated at `x`. @@ -1334,8 +1348,8 @@ def interp(x, xp, fp, left=None, right=None, period=None): -------- scipy.interpolate - Notes - ----- + Warnings + -------- The x-coordinate sequence is expected to be increasing, but this is not explicitly enforced. However, if the sequence `xp` is non-increasing, interpolation results are meaningless. @@ -1447,7 +1461,7 @@ def angle(z, deg=False): The counterclockwise angle from the positive real axis on the complex plane in the range ``(-pi, pi]``, with dtype as numpy.float64. - ..versionchanged:: 1.16.0 + .. versionchanged:: 1.16.0 This function works on subclasses of ndarray like `ma.array`. See Also @@ -1661,6 +1675,7 @@ def trim_zeros(filt, trim='fb'): [1, 2] """ + first = 0 trim = trim.upper() if 'F' in trim: @@ -1970,8 +1985,8 @@ class vectorize: .. versionadded:: 1.7.0 cache : bool, optional - If `True`, then cache the first function call that determines the number - of outputs if `otypes` is not provided. + If `True`, then cache the first function call that determines the number + of outputs if `otypes` is not provided. .. versionadded:: 1.7.0 @@ -2226,15 +2241,14 @@ class vectorize: ufunc, otypes = self._get_ufunc_and_otypes(func=func, args=args) # Convert args to object arrays first - inputs = [array(a, copy=False, subok=True, dtype=object) - for a in args] + inputs = [asanyarray(a, dtype=object) for a in args] outputs = ufunc(*inputs) if ufunc.nout == 1: - res = array(outputs, copy=False, subok=True, dtype=otypes[0]) + res = asanyarray(outputs, dtype=otypes[0]) else: - res = tuple([array(x, copy=False, subok=True, dtype=t) + res = tuple([asanyarray(x, dtype=t) for x, t in zip(outputs, otypes)]) return res @@ -2303,13 +2317,13 @@ class vectorize: def _cov_dispatcher(m, y=None, rowvar=None, bias=None, ddof=None, - fweights=None, aweights=None): + fweights=None, aweights=None, *, dtype=None): return (m, y, fweights, aweights) @array_function_dispatch(_cov_dispatcher) def cov(m, y=None, rowvar=True, bias=False, ddof=None, fweights=None, - aweights=None): + aweights=None, *, dtype=None): """ Estimate a covariance matrix, given data and weights. @@ -2360,6 +2374,11 @@ def cov(m, y=None, rowvar=True, bias=False, ddof=None, fweights=None, weights can be used to assign probabilities to observation vectors. .. versionadded:: 1.10 + dtype : data-type, optional + Data-type of the result. By default, the return data-type will have + at least `numpy.float64` precision. + + .. versionadded:: 1.20 Returns ------- @@ -2435,13 +2454,16 @@ def cov(m, y=None, rowvar=True, bias=False, ddof=None, fweights=None, if m.ndim > 2: raise ValueError("m has more than 2 dimensions") - if y is None: - dtype = np.result_type(m, np.float64) - else: + if y is not None: y = np.asarray(y) if y.ndim > 2: raise ValueError("y has more than 2 dimensions") - dtype = np.result_type(m, y, np.float64) + + if dtype is None: + if y is None: + dtype = np.result_type(m, np.float64) + else: + dtype = np.result_type(m, y, np.float64) X = array(m, ndmin=2, dtype=dtype) if not rowvar and X.shape[0] != 1: @@ -2521,12 +2543,14 @@ def cov(m, y=None, rowvar=True, bias=False, ddof=None, fweights=None, return c.squeeze() -def _corrcoef_dispatcher(x, y=None, rowvar=None, bias=None, ddof=None): +def _corrcoef_dispatcher(x, y=None, rowvar=None, bias=None, ddof=None, *, + dtype=None): return (x, y) @array_function_dispatch(_corrcoef_dispatcher) -def corrcoef(x, y=None, rowvar=True, bias=np._NoValue, ddof=np._NoValue): +def corrcoef(x, y=None, rowvar=True, bias=np._NoValue, ddof=np._NoValue, *, + dtype=None): """ Return Pearson product-moment correlation coefficients. @@ -2560,6 +2584,11 @@ def corrcoef(x, y=None, rowvar=True, bias=np._NoValue, ddof=np._NoValue): Has no effect, do not use. .. deprecated:: 1.10.0 + dtype : data-type, optional + Data-type of the result. By default, the return data-type will have + at least `numpy.float64` precision. + + .. versionadded:: 1.20 Returns ------- @@ -2582,11 +2611,11 @@ def corrcoef(x, y=None, rowvar=True, bias=np._NoValue, ddof=np._NoValue): for backwards compatibility with previous versions of this function. These arguments had no effect on the return values of the function and can be safely ignored in this and previous versions of numpy. - + Examples - -------- + -------- In this example we generate two random arrays, ``xarr`` and ``yarr``, and - compute the row-wise and column-wise Pearson correlation coefficients, + compute the row-wise and column-wise Pearson correlation coefficients, ``R``. Since ``rowvar`` is true by default, we first find the row-wise Pearson correlation coefficients between the variables of ``xarr``. @@ -2602,11 +2631,11 @@ def corrcoef(x, y=None, rowvar=True, bias=np._NoValue, ddof=np._NoValue): array([[ 1. , 0.99256089, -0.68080986], [ 0.99256089, 1. , -0.76492172], [-0.68080986, -0.76492172, 1. ]]) - - If we add another set of variables and observations ``yarr``, we can + + If we add another set of variables and observations ``yarr``, we can compute the row-wise Pearson correlation coefficients between the variables in ``xarr`` and ``yarr``. - + >>> yarr = rng.random((3, 3)) >>> yarr array([[0.45038594, 0.37079802, 0.92676499], @@ -2628,7 +2657,7 @@ def corrcoef(x, y=None, rowvar=True, bias=np._NoValue, ddof=np._NoValue): 1. ]]) Finally if we use the option ``rowvar=False``, the columns are now - being treated as the variables and we will find the column-wise Pearson + being treated as the variables and we will find the column-wise Pearson correlation coefficients between variables in ``xarr`` and ``yarr``. >>> R3 = np.corrcoef(xarr, yarr, rowvar=False) @@ -2651,7 +2680,7 @@ def corrcoef(x, y=None, rowvar=True, bias=np._NoValue, ddof=np._NoValue): # 2015-03-15, 1.10 warnings.warn('bias and ddof have no effect and are deprecated', DeprecationWarning, stacklevel=3) - c = cov(x, y, rowvar) + c = cov(x, y, rowvar, dtype=dtype) try: d = diag(c) except ValueError: @@ -2768,8 +2797,8 @@ def blackman(M): return array([]) if M == 1: return ones(1, float) - n = arange(0, M) - return 0.42 - 0.5*cos(2.0*pi*n/(M-1)) + 0.08*cos(4.0*pi*n/(M-1)) + n = arange(1-M, M, 2) + return 0.42 + 0.5*cos(pi*n/(M-1)) + 0.08*cos(2.0*pi*n/(M-1)) @set_module('numpy') @@ -2877,8 +2906,8 @@ def bartlett(M): return array([]) if M == 1: return ones(1, float) - n = arange(0, M) - return where(less_equal(n, (M-1)/2.0), 2.0*n/(M-1), 2.0 - 2.0*n/(M-1)) + n = arange(1-M, M, 2) + return where(less_equal(n, 0), 1 + n/(M-1), 1 - n/(M-1)) @set_module('numpy') @@ -2981,8 +3010,8 @@ def hanning(M): return array([]) if M == 1: return ones(1, float) - n = arange(0, M) - return 0.5 - 0.5*cos(2.0*pi*n/(M-1)) + n = arange(1-M, M, 2) + return 0.5 + 0.5*cos(pi*n/(M-1)) @set_module('numpy') @@ -3081,8 +3110,8 @@ def hamming(M): return array([]) if M == 1: return ones(1, float) - n = arange(0, M) - return 0.54 - 0.46*cos(2.0*pi*n/(M-1)) + n = arange(1-M, M, 2) + return 0.54 + 0.46*cos(pi*n/(M-1)) ## Code from cephes for i0 @@ -3178,25 +3207,18 @@ def i0(x): """ Modified Bessel function of the first kind, order 0. - Usually denoted :math:`I_0`. This function does broadcast, but will *not* - "up-cast" int dtype arguments unless accompanied by at least one float or - complex dtype argument (see Raises below). + Usually denoted :math:`I_0`. Parameters ---------- - x : array_like, dtype float or complex + x : array_like of float Argument of the Bessel function. Returns ------- - out : ndarray, shape = x.shape, dtype = x.dtype + out : ndarray, shape = x.shape, dtype = float The modified Bessel function evaluated at each of the elements of `x`. - Raises - ------ - TypeError: array cannot be safely cast to required type - If argument consists exclusively of int dtypes. - See Also -------- scipy.special.i0, scipy.special.iv, scipy.special.ive @@ -3226,12 +3248,16 @@ def i0(x): Examples -------- >>> np.i0(0.) - array(1.0) # may vary - >>> np.i0([0., 1. + 2j]) - array([ 1.00000000+0.j , 0.18785373+0.64616944j]) # may vary + array(1.0) + >>> np.i0([0, 1, 2, 3]) + array([1. , 1.26606588, 2.2795853 , 4.88079259]) """ x = np.asanyarray(x) + if x.dtype.kind == 'c': + raise TypeError("i0 not supported for complex values") + if x.dtype.kind != 'f': + x = x.astype(float) x = np.abs(x) return piecewise(x, [x <= 8.0], [_i0_1, _i0_2]) @@ -3970,11 +3996,10 @@ def _quantile_is_valid(q): # avoid expensive reductions, relevant for arrays with < O(1000) elements if q.ndim == 1 and q.size < 10: for i in range(q.size): - if q[i] < 0.0 or q[i] > 1.0: + if not (0.0 <= q[i] <= 1.0): return False else: - # faster than any() - if np.count_nonzero(q < 0.0) or np.count_nonzero(q > 1.0): + if not (np.all(0 <= q) and np.all(q <= 1)): return False return True @@ -4103,11 +4128,18 @@ def _trapz_dispatcher(y, x=None, dx=None, axis=None): @array_function_dispatch(_trapz_dispatcher) def trapz(y, x=None, dx=1.0, axis=-1): - """ + r""" Integrate along the given axis using the composite trapezoidal rule. - Integrate `y` (`x`) along given axis. - + If `x` is provided, the integration happens in sequence along its + elements - they are not sorted. + + Integrate `y` (`x`) along each 1d slice on the given axis, compute + :math:`\int y(x) dx`. + When `x` is specified, this integrates along the parametric curve, + computing :math:`\int_t y(t) dt = + \int_t y(t) \left.\frac{dx}{dt}\right|_{x=x(t)} dt`. + Parameters ---------- y : array_like @@ -4123,9 +4155,12 @@ def trapz(y, x=None, dx=1.0, axis=-1): Returns ------- - trapz : float - Definite integral as approximated by trapezoidal rule. - + trapz : float or ndarray + Definite integral of 'y' = n-dimensional array as approximated along + a single axis by the trapezoidal rule. If 'y' is a 1-dimensional array, + then the result is a float. If 'n' is greater than 1, then the result + is an 'n-1' dimensional array. + See Also -------- sum, cumsum @@ -4154,6 +4189,20 @@ def trapz(y, x=None, dx=1.0, axis=-1): 8.0 >>> np.trapz([1,2,3], dx=2) 8.0 + + Using a decreasing `x` corresponds to integrating in reverse: + + >>> np.trapz([1,2,3], x=[8,6,4]) + -8.0 + + More generally `x` is used to integrate along a parametric curve. + This finds the area of a circle, noting we repeat the sample which closes + the curve: + + >>> theta = np.linspace(0, 2 * np.pi, num=1000, endpoint=True) + >>> np.trapz(np.cos(theta), x=np.sin(theta)) + 3.141571941375841 + >>> a = np.arange(6).reshape(2, 3) >>> a array([[0, 1, 2], @@ -4162,7 +4211,6 @@ def trapz(y, x=None, dx=1.0, axis=-1): array([1.5, 2.5, 3.5]) >>> np.trapz(a, axis=1) array([2., 8.]) - """ y = asanyarray(y) if x is None: @@ -4267,10 +4315,9 @@ def meshgrid(*xi, copy=True, sparse=False, indexing='xy'): See Also -------- - index_tricks.mgrid : Construct a multi-dimensional "meshgrid" - using indexing notation. - index_tricks.ogrid : Construct an open multi-dimensional "meshgrid" - using indexing notation. + mgrid : Construct a multi-dimensional "meshgrid" using indexing notation. + ogrid : Construct an open multi-dimensional "meshgrid" using indexing + notation. Examples -------- @@ -4298,7 +4345,8 @@ def meshgrid(*xi, copy=True, sparse=False, indexing='xy'): >>> y = np.arange(-5, 5, 0.1) >>> xx, yy = np.meshgrid(x, y, sparse=True) >>> z = np.sin(xx**2 + yy**2) / (xx**2 + yy**2) - >>> h = plt.contourf(x,y,z) + >>> h = plt.contourf(x, y, z) + >>> plt.axis('scaled') >>> plt.show() """ |