# Maths and statistics


<!-- WARNING: THIS FILE WAS AUTOGENERATED! DO NOT EDIT! -->

``` python
import rustymath as rm
```

## Elementary functions

Functions accept scalars or arrays. Integer inputs to float functions
compute at an appropriate float width. Float arrays retain their width:

``` python
rm.exp([0., 1, 2]), rm.log([1., 10, 100]), rm.sin([0., 0.5, 1])
```

    (Array([1.0, 2.7182818, 7.3890561], dtype=float64),
     Array([0.0, 2.3025851, 4.6051702], dtype=float64),
     Array([0.0, 0.47942554, 0.84147098], dtype=float64))

Functions that NumPy lacks keep their own names. `root(n, y)` gives the
`n`th root of `y`:

``` python
rm.root(3, [8., 27., 1000.])
```

    Array([2.0, 3.0, 10.0], dtype=float64)

`round(x, decimals)` rounds values to a decimal place, with ties to
even. `rint` rounds to whole values. Neither changes how many digits the
array’s `repr` displays. `degrees` and `radians` convert angular units.

``` python
rm.round([1.2345, -2.71828], 2), rm.degrees([0., 1.5707963267948966])
```

    (Array([1.23, -2.72], dtype=float64), Array([0.0, 90.0], dtype=float64))

The Rust names follow Rust’s standard library instead. Python’s `log`,
`power`, `arcsin` and `arctan2` are Rust’s `ln`, `pow`, `asin` and
`atan2`.

## Complex numbers

A complex array gives complex results, at each function’s principal
value. For `exp`, `log`, `power`, `sin`, `cos`, `sinh` and `cosh`,
rustymath computes a complex array through its vectorized real loops:

``` python
rm.sqrt([-4+0j, 3+4j]), rm.log([-1+0j])
```

    (Array([2j, (2+1j)], dtype=complex128), Array([3.1415927j], dtype=complex128))

## Special functions

The special functions take SciPy’s names, such as `gamma`, `gammaln`,
`erf`, `erfinv`, `betainc` and `lambertw`:

``` python
rm.gamma([0.5, 5, 10]), rm.erfinv(0.5)
```

    (Array([1.7724539, 24.0, 362880.0], dtype=float64), 0.4769362762044699)

## Reductions and statistics

With no `axis`, reductions summarize all items. Pass an axis position,
name or tuple of axes to summarize those axes. Here each row holds one
series. Clipping bounds the original values in one pass:

``` python
measurements = rm.array([[1., 2, 3], [4, 5, 6]])
rm.mean(measurements, axis=1), rm.std(measurements, axis=1), rm.clip(measurements, 2., 5.)
```

    (Array([2.0, 5.0], dtype=float64),
     Array([0.81649658, 0.81649658], dtype=float64),
     Array([[2.0, 2.0, 3.0], [4.0, 5.0, 5.0]], dtype=float64))

`diff` takes adjacent differences along an axis; `ptp` gives maximum
minus minimum, and `trace` sums the diagonal of each matrix in the last
two axes. Integer differences wrap where they don’t fit, as NumPy’s do.
Boolean differences report a change (XOR).

``` python
rm.diff(measurements), rm.ptp(measurements, axis=1), rm.trace(measurements)
```

    (Array([[1.0, 1.0], [1.0, 1.0]], dtype=float64),
     Array([2.0, 2.0], dtype=float64),
     6.0)

`var` and `std` divide by `n - ddof`, with `ddof=0` by default.
`average` accepts weights under the same leading-axis broadcasting rule
as other functions. `cov` and `corrcoef` treat rows as variables unless
`rowvar=False`. Float16 summaries accumulate and return float32; complex
variance and norms return real values:

``` python
rm.cov(measurements), rm.corrcoef(measurements)
```

    (Array([[1.0, 1.0], [1.0, 1.0]], dtype=float64),
     Array([[1.0, 1.0], [1.0, 1.0]], dtype=float64))

`norm` computes an elementwise p-norm over the selected axes. It does
not switch to a matrix-induced norm for two-dimensional input:

``` python
rm.norm([3., 4.]), rm.isfinite([1., float('inf'), float('nan')])
```

    (5.0, Array([True, False, False], dtype=bool))

Scans keep a running result. `accumulate(f, x)` applies a native binary
function from the left, over the flattened input or one `axis` by
position or name. A scalar `start` precedes each lane without adding an
output item. Here subtraction spends a budget of 10 in order.
`logical_scan` takes any Boolean dyad as a 2 by 2 truth table; this
table computes running exclusive-or:

``` python
rm.accumulate(rm.subtract, [3., 2, 1], start=10), rm.logical_scan([False, True, True, False], [[False, True], [True, False]])
```

    (Array([7.0, 5.0, 4.0], dtype=float64),
     Array([False, True, False, False], dtype=bool))

`sliding(f, x, width, axis=-1)` folds each window of `width`
neighbouring items by a native binary function with an associative
reduction: `add`, `multiply`, `maximum`, `minimum` or a logical one. A
lane of n items gives n - width + 1 folds. Windows of more than 8 items
fold by van Herk-Gil-Werman’s method, at a cost for each item that
doesn’t grow with the width. `sliding_argmax` and `sliding_argmin` give
the position along the axis of each window’s first maximum or minimum:
5, where the 9 is, for each window that holds it.

``` python
v = [3., 1, 4, 1, 5, 9, 2, 6]
rm.sliding(rm.add, [1., 2, 3, 4, 5], 3), rm.sliding(rm.maximum, v, 4), rm.sliding_argmax(v, 4)
```

    (Array([6.0, 9.0, 12.0], dtype=float64),
     Array([4.0, 5.0, 9.0, 9.0, 9.0], dtype=float64),
     Array([2, 4, 5, 5, 5], dtype=int64))

## Correlation and convolution

`correlate(x, kernel, mode='valid')` slides `kernel` along `x` and gives
the sum of products of the window and the kernel at each position.
`convolve` reverses the kernel first. The modes are SciPy’s: `'valid'`
keeps windows inside `x`, `'full'` keeps every overlap, and `'same'`
gives one position for each of `x`’s, centred as NumPy centres them. The
kernel isn’t conjugated.

``` python
rm.correlate([1., 2, 3, 4], [1., 10]), rm.convolve([1., 2, 3, 4], [1., 10], 'full')
```

    (Array([21.0, 32.0, 43.0], dtype=float64),
     Array([1.0, 12.0, 23.0, 34.0, 40.0], dtype=float64))

A kernel can have more axes than the axes it slides along. `slide` sets
how many leading axes of `x` it slides along, and the kernel contracts
`x`’s other axes, as a convolution layer contracts an image’s channels.
The kernel’s remaining axes give the outputs at each position. Here a 3
by 3 kernel moves 2 positions at a time over a 9 by 8 image with 3
channels, and gives 4 outputs at each of 4 by 3 positions:

``` python
rm.correlate(rm.ones((9, 8, 3)), rm.ones((3, 3, 3, 4)), slide=2, stride=2).shape
```

    (4, 3, 4)

`f` and `g` replace the sum and the product, as in `linalg.inner`. The
maximum over each window of sums is a dilation:

``` python
rm.correlate([0., 5, 0, 0, 3], [0., 1, 0], f=rm.maximum, g=rm.add)
```

    Array([6.0, 5.0, 3.0], dtype=float64)

`overlap_add(frames, step=1, axis=-2)` is the transpose of taking
windows. It adds each frame back where its window was read. Items that
several frames cover add up. An inverse short-time Fourier transform
uses it:

``` python
rm.overlap_add([[1., 1, 1], [1., 1, 1], [1., 1, 1]], 2)
```

    Array([1.0, 1.0, 2.0, 1.0, 2.0, 1.0, 1.0], dtype=float64)

## Tolerant comparison

APL compares floats within a relative tolerance. Two floats are equal
when their difference is less than the tolerance times the larger
magnitude. With the usual tolerance, `0.1 + 0.2` equals `0.3`.
`rm.equal`, `not_equal`, `less`, `less_equal`, `greater` and
`greater_equal` compare this way. Their tolerance `ct` defaults to
`rm.CT`, 1e-14. Pass `ct=0` to compare exactly:

``` python
rm.equal(0.1 + 0.2, 0.3), rm.equal(0.1 + 0.2, 0.3, ct=0)
```

    (True, False)

`floor`, `ceil`, `remainder`, `gcd` and `lcm` take `ct` too. A value
within tolerance of a whole number counts as that whole number:

``` python
rm.floor([2.9999999999999996, 2.5])
```

    Array([3.0, 2.0], dtype=float64)

`remainder` takes its arguments in NumPy’s order. `rm.remainder(x, y)`
is the Rust `residue(y, x)`.

## Checking accuracy

``` python
import numpy as np
```

To check the accuracy yourself, compare with NumPy. `np.spacing` gives
the size of one ulp at each value. Dividing the difference by it gives
the distance between the two results in ulps:

``` python
xs = np.linspace(-10, 10, 100_001)
y = np.exp(xs)
(np.abs(rm.exp(xs) - y) / np.spacing(y)).max()
```

    np.float64(2.0)