Statistics
**********

Octave has support for various statistical methods.  The emphasis is on basic descriptive statistics, but the Octave Forge statistics package includes probability distributions, statistical tests, random number generation, and much more.

The functions that analyze data all assume that multi-dimensional data is arranged in a matrix where each row is an observation, and each column is a variable.  Thus, the matrix defined by

     a = [ 0.9, 0.7;
           0.1, 0.1;
           0.5, 0.4 ];

contains three observations from a two-dimensional distribution.  While this is the default data arrangement, most functions support different arrangements.

It should be noted that the statistics functions don't test for data containing NaN, NA, or Inf.  These values need to be detected and dealt with explicitly. See isnan, isna, isinf, isfinite.

Descriptive Statistics
======================

One principal goal of descriptive statistics is to represent the essence of a large data set concisely.  Octave provides the mean, median, and mode functions which all summarize a data set with just a single number corresponding to the central tendency of the data.

 -- M =: mean (X)
 -- M =: mean (X, DIM)
 -- M =: mean (X, VECDIM)
 -- M =: mean (X, "all")
 -- M =: mean (..., NANFLAG)
 -- M =: mean (..., OUTTYPE)
     Compute the mean of the elements of X.

     If X is a vector, then mean (X) returns the mean of the elements in X defined as

          mean (X) = SUM_i X(i) / N

     where N is the number of elements in X.

     If X is an array, then mean(X) computes the mean along the first non-singleton dimension of X.

     The optional variable DIM forces mean to operate over the specified dimension, which must be a positive integer-valued number. Specifying any singleton dimension in X, including any dimension exceeding ndims (X), will result in a mean equal to X.

     Specifying the dimensions as  VECDIM, a vector of non-repeating dimensions, will return the mean over the array slice defined by VECDIM.  If VECDIM indexes all dimensions of X, then it is equivalent to the option "all".  Any dimension in VECDIM greater than ndims (X) is ignored.

     Specifying the dimension as "all" will force mean to operate on all elements of X, and is equivalent to mean (X(:)).

     The optional input OUTTYPE specifies the data type that is returned. OUTTYPE can take the following values:

     'default' : Output is of type double, unless the input is
          single in which case the output is of type single.

     'double' : Output is of type double.

     'native' : Output is of the same type as the input as reported
          by (class (X)), unless the input is logical in which case the output is of type double.

     The optional variable NANFLAG specifies whether to include or exclude NaN values from the calculation using any of the previously specified input argument combinations.  The default value for NANFLAG is "includenan" which keeps NaN values in the calculation.  To exclude NaN values set the value of NANFLAG to "omitnan".  The output will still contain NaN values if X consists of all NaN values in the operating dimension.

     median, mode, movmean

 -- M =: median (X)
 -- M =: median (X, DIM)
 -- M =: median (X, VECDIM)
 -- M =: median (X, "all")
 -- M =: median (..., NANFLAG)
 -- M =: median (..., OUTTYPE)
     Compute the median value of the elements of X.

     When the elements of X are sorted, say S = sort (X), the median is defined as

                       |  S(ceil (N/2))          N odd
          median (X) = |
                       | (S(N/2) + S(N/2+1))/2   N even

     If X is an array, then median (X) operates along the first non-singleton dimension of X.

     The optional variable DIM forces median to operate over the specified dimension, which must be a positive integer-valued number. Specifying any singleton dimension in X, including any dimension exceeding ndims (X), will result in a median equal to X.

     Specifying the dimensions as  VECDIM, a vector of non-repeating dimensions, will return the median over the array slice defined by VECDIM.  If VECDIM indexes all dimensions of X, then it is equivalent to the option "all".  Any dimension in VECDIM greater than ndims (X) is ignored.

     Specifying the dimension as "all" will force median to operate on all elements of X, and is equivalent to median (X(:)).

     median (..., OUTTYPE) returns the median with a specified data type, using any of the input arguments in the previous syntaxes. OUTTYPE can take the following values:

     "default"
          Output is of type double, unless the input is single in which case the output is of type single.

     "double"
          Output is of type double.

     "native".
          Output is of the same type as the input (class (X)), unless the input is logical in which case the output is of type double.

     The optional variable NANFLAG specifies whether to include or exclude NaN values from the calculation using any of the previously specified input argument combinations.  The default value for NANFLAG is "includenan" which keeps NaN values in the calculation.  To exclude NaN values set the value of NANFLAG to "omitnan". The output will still contain NaN values if X consists of all NaN values in the operating dimension.

     mean, mode, movmedian

 -- M =: mode (X)
 -- M =: mode (X, DIM)
 -- [M, F, C] =: mode (...)
     Compute the most frequently occurring value in a dataset (mode).

     mode determines the frequency of values along the first non-singleton dimension and returns the value with the highest frequency.  If two, or more, values have the same frequency mode returns the smallest.

     If the optional argument DIM is given, operate along this dimension.

     The return variable F is the number of occurrences of the mode in the dataset.

     The cell array C contains all of the elements with the maximum frequency. mean, median

Using just one number, such as the mean, to represent an entire data set may not give an accurate picture of the data.  One way to characterize the fit is to measure the dispersion of the data.  Octave provides several functions for measuring dispersion.

 -- [S, L] =: bounds (X)
 -- [S, L] =: bounds (X, DIM)
 -- [S, L] =: bounds (..., "nanflag")
     Return the smallest and largest values of the input data X.

     If X is a vector, the bounds are calculated over the elements of X.  If X is a matrix, the bounds are calculated for each column. For a multi-dimensional array, the bounds are calculated over the first non-singleton dimension.

     If the optional argument DIM is given, operate along this dimension.

     The optional argument "nanflag" defaults to "omitnan" which does not include NaN values in the result.  If the argument "includenan" is given, and there is a NaN present, then the result for both smallest (S) and largest (L) elements will be NaN.

     The bounds are a quickly computed measure of the dispersion of a data set, but are less accurate than iqr if there are outlying data points. range, iqr, mad, std

 -- Y =: range (X)
 -- Y =: range (X, DIM)
     Return the range, i.e., the difference between the maximum and the minimum of the input data.

     If X is a vector, the range is calculated over the elements of X.  If X is a matrix, the range is calculated over each column of X.

     If the optional argument DIM is given, operate along this dimension.

     The range is a quickly computed measure of the dispersion of a data set, but is less accurate than iqr if there are outlying data points. bounds, iqr, mad, std

 -- Z =: iqr (X)
 -- Z =: iqr (X, DIM)
 -- Z =: iqr (X, "ALL")
     Return the interquartile range of X, defined as the distance between the 25th and 75th percentile values of X calculated using: quantile (x, [0.25 0.75])

     If X is a vector, iqr (X) will operate on the data in X.

     If X is a matrix, iqr (X) will operate independently on each column in X returning a row vector Z.

     If X is a n-dimensional array, iqr (X) will operate independently on the first non-singleton dimension in X, returning an array Z the same shape as X with the non-singleton dimenion reduced to 1.

     The optional variable DIM can be used to force iqr to operate over the specified dimension.  DIM can either be a scalar dimension or a vector of non-repeating dimensions over which to operate.  In either case DIM must be positive integers.  A vector DIM concatenates all specified dimensions for independent operation by iqr.

     Specifying dimension "ALL" will force iqr to operate on all elements of X, and is equivalent to iqr (X(:)). Similarly, specifying a vector dimension including all non-singleton dimensions of X is equivalent to iqr (X, "ALL").

     If X is a scalar, or only singleton dimensions are specified for DIM, the output will be zeros (size (X)).

     As a measure of dispersion, the interquartile range is less affected by outliers than either range or std.

     bounds, mad, range, std, prctile, quantile

 -- M =: mad (X)
 -- M =: mad (X, OPT)
 -- M =: mad (X, OPT, DIM)
 -- M =: mad (X, OPT, VECDIM)
 -- M =: mad (X, OPT, "all")
     Compute the mean or median absolute deviation (MAD) of the elements of X.

     The mean absolute deviation is defined as

          MAD = mean (abs (X - mean (X)))

     The median absolute deviation is defined as

          MAD = median (abs (X - median (X)))

     If X is a vector, compute mad for each element in X.  If X is an array the calculation is performed over the first non-singleton dimension.

     mad excludes NaN values from calculation similar to using the omitnan option in var, mean, and median.

     The optional argument OPT determines whether mean or median absolute deviation is calculated.  The default is 0 which corresponds to mean absolute deviation; a value of 1 corresponds to median absolute deviation. Passing an empty input [] defaults to mean absolute deviation (OPT = 0).

     The optional argument DIM forces mad to operate along the specified dimension.  Specifying any singleton dimension in X, including any dimension exceeding ndims (X), will result in an output of 0.

     Specifying the dimension as VECDIM, a vector of non-repeating dimensions, will return the mad over the array slice defined by VECDIM.  If VECDIM indexes all dimensions of X, then it is equivalent to the option "all".  Any dimension included in VECDIM greater than ndims (X) is ignored.

     Specifying the dimension as "all" will force mad to operate on all elements of X, and is equivalent to mad (X(:)).

     As a measure of dispersion, mad is less affected by outliers than std. bounds, range, iqr, std, mean, median

 -- Y =: meansq (X)
 -- Y =: meansq (X, DIM)
     Compute the mean square of the elements of the vector X.

     The mean square is defined as

          meansq (X) = 1/N SUM_i X(i)^2

     where N is the length of the X vector.

     If X is a matrix, return a row vector containing the mean square of each column.

     If the optional argument DIM is given, operate along this dimension. var, std, moment

 -- S =: std (X)
 -- S =: std (X, W)
 -- S =: std (X, W, DIM)
 -- S =: std (X, W, VECDIM)
 -- S =: std (X, W, "ALL")
 -- S =: std (..., NANFLAG)
 -- [S, M] =: std (...)
     Compute the standard deviation of the elements of the vector X.

     The standard deviation is defined as

          std (X) = sqrt ((1 / (N-1)) * SUM_i ((X(i) - mean(X))^2))

     where N is the number of elements of X.

     If X is an array, compute the standard deviation along the first non-singleton dimensions of X.

     The optional argument W determines the weighting scheme to use.  Valid values are:

     0 [default]:
          Normalize with N-1 (population standard deviation).  This provides the square root of the best unbiased estimator of the standard deviation.

     1:
          Normalize with N (sample standard deviation).  This provides the square root of the second moment around the mean.

     a vector:
          Compute the weighted standard deviation with non-negative weights. The length of W must equal the size of X in the operating dimension.  NaN values are permitted in W, will be multiplied with the associated values in X, and can be excluded by the NANFLAG option.

     an array:
          Similar to vector weights, but W must be the same size as X.  If the operating dimension is supplied as VECDIM or "all" and W is not a scalar, W must be an same-sized array.

     Note: W must always be specified before specifying any of the following dimension options.  To use the default value for W you may pass an empty input argument [].

     The optional variable DIM forces std to operate over the specified dimension, which must be a positive integer-valued number. Specifying any singleton dimension in X, including any dimension exceeding ndims (X), will result in a standard deviation of 0.

     Specifying the dimensions as  VECDIM, a vector of non-repeating dimensions, will return the standard deviation calculated over the array slice defined by VECDIM.  If VECDIM indexes all dimensions of X, then it is equivalent to the option "all".  Any dimension in VECDIM greater than ndims (X) is ignored.

     Specifying the dimension as "all" will force std to operate on all elements of X, and is equivalent to std (X(:)).

     The optional variable NANFLAG specifies whether to include or exclude NaN values from the calculation using any of the previously specified input argument combinations.  The default value for NANFLAG is "includenan" which keeps NaN values in the calculation.  To exclude NaN values set the value of NANFLAG to "omitnan". The output will still contain NaN values if X consists of all NaN values in the operating dimension.

     The optional second output variable M contains the mean of the elements of X used to calculate the standard deviation.  If V is the weighted standard deviation, then M is also the weighted mean.

     var, bounds, mad, range, iqr, mean, median

In addition to knowing the size of a dispersion it is useful to know the shape of the data set.  For example, are data points massed to the left or right of the mean?  Octave provides several common measures to describe the shape of the data set.  Octave can also calculate moments allowing arbitrary shape measures to be developed.

 -- V =: var (X)
 -- V =: var (X, W)
 -- V =: var (X, W, DIM)
 -- V =: var (X, W, VECDIM)
 -- V =: var (X, W, "all")
 -- V =: var (..., NANFLAG)
 -- [V, M] =: var (...)
     Compute the variance of the elements of the vector X.

     The variance is defined as

          var (X) = (1 / (N-1)) * SUM_i ((X(i) - mean(X))^2)

     where N is the number of elements of X.

     If X is an array, compute the variance along the first non-singleton dimensions of X.

     The optional argument W determines the weighting scheme to use.  Valid values are:

     0 [default]:
          Normalize with N-1 (population variance).  This provides the square root of the best unbiased estimator of the variance.

     1:
          Normalize with N (sample variance).  This provides the square root of the second moment around the mean.

     a vector:
          Compute the weighted variance with non-negative weights.  The length of W must equal the size of X in the operating dimension.  NaN values are permitted in W, will be multiplied with the associated values in X, and can be excluded by the NANFLAG option.

     an array:
          Similar to vector weights, but W must be the same size as X.  If the operating dimension is supplied as VECDIM or "all" and W is not a scalar, W must be an same-sized array.

     Note: W must always be specified before specifying any of the following dimension options.  To use the default value for W you may pass an empty input argument [].

     The optional variable DIM forces var to operate over the specified dimension, which must be a positive integer-valued number. Specifying any singleton dimension in X, including any dimension exceeding ndims (X), will result in a variance of 0.

     Specifying the dimensions as  VECDIM, a vector of non-repeating dimensions, will return the variance calculated over the array slice defined by VECDIM.  If VECDIM indexes all dimensions of X, then it is equivalent to the option "all".  Any dimension in VECDIM greater than ndims (X) is ignored.

     Specifying the dimension as "all" will force var to operate on all elements of X, and is equivalent to var (X(:)).

     The optional variable NANFLAG specifies whether to include or exclude NaN values from the calculation using any of the previously specified input argument combinations.  The default value for NANFLAG is "includenan" which keeps NaN values in the calculation.  To exclude NaN values set the value of NANFLAG to "omitnan". The output will still contain NaN values if X consists of all NaN values in the operating dimension.

     The optional second output variable M contains the mean of the elements of X used to calculate the variance.  If V is the weighted variance, then M is also the weighted mean.

     std, mean, cov, skewness, kurtosis, moment

 -- Y =: skewness (X)
 -- Y =: skewness (X, FLAG)
 -- Y =: skewness (X, FLAG, DIM)
     Compute the sample skewness of the elements of X.

     The sample skewness is defined as

                         mean ((X - mean (X)).^3)
          skewness (X) = ------------------------.
                                std (X).^3

     The optional argument FLAG controls which normalization is used. If FLAG is equal to 1 (default value, used when FLAG is omitted or empty), return the sample skewness as defined above.  If FLAG is equal to 0, return the adjusted skewness coefficient instead:

                            sqrt (N*(N-1))   mean ((X - mean (X)).^3)
          skewness (X, 0) = -------------- * ------------------------.
                                (N - 2)             std (X).^3

     where N is the length of the X vector.

     The adjusted skewness coefficient is obtained by replacing the sample second and third central moments by their bias-corrected versions.

     If X is a matrix, or more generally a multi-dimensional array, return the skewness along the first non-singleton dimension.  If the optional DIM argument is given, operate along this dimension.

     var, kurtosis, moment

 -- Y =: kurtosis (X)
 -- Y =: kurtosis (X, FLAG)
 -- Y =: kurtosis (X, FLAG, DIM)
     Compute the sample kurtosis of the elements of X.

     The sample kurtosis is defined as

               mean ((X - mean (X)).^4)
          k1 = ------------------------
                      std (X).^4

     The optional argument FLAG controls which normalization is used. If FLAG is equal to 1 (default value, used when FLAG is omitted or empty), return the sample kurtosis as defined above.  If FLAG is equal to 0, return the "bias-corrected" kurtosis coefficient instead:

                        N - 1
          k0 = 3 + -------------- * ((N + 1) * k1 - 3 * (N - 1))
                   (N - 2)(N - 3)

     where N is the length of the X vector.

     The bias-corrected kurtosis coefficient is obtained by replacing the sample second and fourth central moments by their unbiased versions.  It is an unbiased estimate of the population kurtosis for normal populations.

     If X is a matrix, or more generally a multi-dimensional array, return the kurtosis along the first non-singleton dimension.  If the optional DIM argument is given, operate along this dimension.

     var, skewness, moment

 -- M =: moment (X, P)
 -- M =: moment (X, P, TYPE)
 -- M =: moment (X, P, DIM)
 -- M =: moment (X, P, TYPE, DIM)
 -- M =: moment (X, P, DIM, TYPE)
     Compute the P-th central moment of the vector X.

     The P-th central moment of X is defined as:

          1/N SUM_i (X(i) - mean(X))^P

     where N is the length of the X vector.

     If X is a matrix, return the row vector containing the P-th central moment of each column.

     If the optional argument DIM is given, operate along this dimension.

     The optional string TYPE specifies the type of moment to be computed. Valid options are:

     "c"
          Central Moment (default).

     "a"
     "ac"
          Absolute Central Moment.  The moment about the mean ignoring sign defined as

               1/N SUM_i (abs (X(i) - mean(X)))^P

     "r"
          Raw Moment.  The moment about zero defined as

               moment (X) = 1/N SUM_i X(i)^P

     "ar"
          Absolute Raw Moment.  The moment about zero ignoring sign defined as

               1/N SUM_i ( abs (X(i)) )^P

     If both TYPE and DIM are given they may appear in any order. var, skewness, kurtosis

 -- Q =: quantile (X)
 -- Q =: quantile (X, P)
 -- Q =: quantile (X, P, DIM)
 -- Q =: quantile (X, P, DIM, METHOD)
     For a sample, X, calculate the quantiles, Q, corresponding to the cumulative probability values in P.  All non-numeric values (NaNs) of X are ignored.

     If X is a matrix, compute the quantiles for each column and return them in a matrix, such that the i-th row of Q contains the P(i)th quantiles of each column of X.

     If P is unspecified, return the quantiles for [0.00 0.25 0.50 0.75 1.00]. The optional argument DIM determines the dimension along which the quantiles are calculated.  If DIM is omitted it defaults to the first non-singleton dimension.

     The methods available to calculate sample quantiles are the nine methods used by R (https://www.r-project.org/).  The default value is METHOD = 5.

     Discontinuous sample quantile methods 1, 2, and 3

      1. Method 1: Inverse of empirical distribution function.

      2. Method 2: Similar to method 1 but with averaging at discontinuities.

      3. Method 3: SAS definition: nearest even order statistic.

     Continuous sample quantile methods 4 through 9, where P(k) is the linear interpolation function respecting each method's representative cdf.

      4. Method 4: P(k) = k / N. That is, linear interpolation of the empirical cdf, where N is the length of P.

      5. Method 5: P(k) = (k - 0.5) / N. That is, a piecewise linear function where the knots are the values midway through the steps of the empirical cdf.

      6. Method 6: P(k) = k / (N + 1).

      7. Method 7: P(k) = (k - 1) / (N - 1).

      8. Method 8: P(k) = (k - 1/3) / (N + 1/3). The resulting quantile estimates are approximately median-unbiased regardless of the distribution of X.

      9. Method 9: P(k) = (k - 3/8) / (N + 1/4). The resulting quantile estimates are approximately unbiased for the expected order statistics if X is normally distributed.

     Hyndman and Fan (1996) recommend method 8.  Maxima, S, and R (versions prior to 2.0.0) use 7 as their default.  Minitab and SPSS use method 6.  MATLAB uses method 5.

     References:

       * Becker, R. A., Chambers, J. M. and Wilks, A. R. (1988) The New S Language.  Wadsworth & Brooks/Cole.

       * Hyndman, R. J. and Fan, Y. (1996) Sample quantiles in statistical packages, American Statistician, 50, 361--365.

       * R: A Language and Environment for Statistical Computing; https://cran.r-project.org/doc/manuals/fullrefman.pdf.

     Examples:

          x = randi (1000, [10, 1]);  # Create empirical data in range 1-1000
          q = quantile (x, [0, 1]);   # Return minimum, maximum of distribution
          q = quantile (x, [0.25 0.5 0.75]); # Return quartiles of distribution

     prctile

 -- Q =: prctile (X)
 -- Q =: prctile (X, P)
 -- Q =: prctile (X, P, DIM)
     For a sample X, compute the quantiles, Q, corresponding to the cumulative probability values, P, in percent.

     If X is a matrix, compute the percentiles for each column and return them in a matrix, such that the i-th row of Q contains the P(i)th percentiles of each column of X.

     If P is unspecified, return the quantiles for [0 25 50 75 100].

     The optional argument DIM determines the dimension along which the percentiles are calculated.  If DIM is omitted it defaults to the first non-singleton dimension.

     Programming Note: All non-numeric values (NaNs) of X are ignored. quantile

A summary view of a data set can be generated quickly with the statistics function.

 -- STATS =: statistics (X)
 -- STATS =: statistics (X, DIM)
     Return a vector with the minimum, first quartile, median, third quartile, maximum, mean, standard deviation, skewness, and kurtosis of the elements of the vector X.

     If X is a matrix, calculate statistics over the first non-singleton dimension.

     If the optional argument DIM is given, operate along this dimension. min, max, median, mean, std, skewness, kurtosis

Statistics on Sliding Windows of Data
=====================================

It is often useful to calculate descriptive statistics over a subsection (i.e., window) of a full dataset.  Octave provides the function movfun which will call an arbitrary function handle with windows of data and accumulate the results.  Many of the most commonly desired functions, such as the moving average over a window of data (movmean), are already provided.

 -- Y =: movfun (FCN, X, WLEN)
 -- Y =: movfun (FCN, X, [NB, NA])
 -- Y =: movfun (..., "PROPERTY", VALUE)

     Apply function FCN to a moving window of length WLEN on data X.

     If WLEN is a scalar, the function FCN is applied to a moving window of length WLEN.  When WLEN is an odd number the window is symmetric and includes (WLEN - 1) / 2 elements on either side of the central element.  For example, when calculating the output at index 5 with a window length of 3, movfun uses data elements [4, 5, 6].  If WLEN is an even number, the window is asymmetric and has WLEN/2 elements to the left of the central element and WLEN/2 - 1 elements to the right of the central element.  For example, when calculating the output at index 5 with a window length of 4, movfun uses data elements [3, 4, 5, 6].

     If WLEN is an array with two elements [NB, NA], the function is applied to a moving window -NB:NA.  This window includes NB number of elements before the current element and NA number of elements after the current element. The current element is always included.  For example, given WLEN = [3, 0], the data used to calculate index 5 is [2, 3, 4, 5].

     During calculations the data input X is reshaped into a 2-dimensional WLEN-by-N matrix and FCN is called on this new matrix. Therefore, FCN must accept an array input argument and apply the computation along dimension 1, i.e., down the columns of the array.

     When applied to an array (possibly multi-dimensional) with N columns, FCN may return a result in either of two formats: Format 1) an array of size 1-by-N-by-DIM3-by-...-by-DIMN.  This is the typical output format from Octave core functions.  Type demo ("movfun", 5) for an example of this use case. Format 2) a row vector of length N * NUMEL_HIGHER_DIMS where NUMEL_HIGHER_DIMS is prod (size (X)(3:end)).  The output of FCN for the i-th input column must be found in the output at indices i:N:(N*NUMEL_HIGHER_DIMS). This format is useful when concatenating functions into arrays, or when using nthargout.  Type demo ("movfun", 6) for an example of this case.

     The calculation can be controlled by specifying PROPERTY/VALUE pairs.  Valid properties are

     "dim"
          Operate along the dimension specified, rather than the default of the first non-singleton dimension.

     "Endpoints"

          This property controls how results are calculated at the boundaries (endpoints) of the window.  Possible values are:

          "shrink"  (default)
               The window is truncated at the beginning and end of the array to exclude elements for which there is no source data.  For example, with a window of length 3, Y(1) = FCN (X(1:2)), and Y(end) = FCN (X(end-1:end)).

          "discard"
               Any Y values that use a window extending beyond the original data array are deleted.  For example, with a 10-element data vector and a window of length 3, the output will contain only 8 elements.  The first element would require calculating the function over indices [0, 1, 2] and is therefore discarded.  The last element would require calculating the function over indices [9, 10, 11] and is therefore discarded.

          "fill"
               Any window elements outside the data array are replaced by NaN.  For example, with a window of length 3, Y(1) = FCN ([NaN, X(1:2)]), and Y(end) = FCN ([X(end-1:end), NaN]). This option usually results in Y having NaN values at the boundaries, although it is influenced by how FCN handles NaN, and also by the property "nancond".

          USER_VALUE
               Any window elements outside the data array are replaced by the specified value USER_VALUE which must be a numeric scalar.  For example, with a window of length 3, Y(1) = FCN ([USER_VALUE, X(1:2)]), and Y(end) = FCN ([X(end-1:end), USER_VALUE]). A common choice for USER_VALUE is 0.

          "same"
               Any window elements outside the data array are replaced by the value of X at the boundary.  For example, with a window of length 3, Y(1) = FCN ([X(1), X(1:2)]), and Y(end) = FCN ([X(end-1:end), X(end)]).

          "periodic"
               The window is wrapped so that any missing data elements are taken from the other side of the data.  For example, with a window of length 3, Y(1) = FCN ([X(end), X(1:2)]), and Y(end) = FCN ([X(end-1:end), X(1)]).

          Note that for some of these choices, the window size at the boundaries is not the same as for the central part, and FCN must work in these cases.

     "nancond"
          Controls whether NaN and NA values should be included (value: "includenan"), or excluded (value: "omitnan"), from the data passed to FCN.  The default is "includenan".  Caution: The "omitnan" option is not yet implemented.

     "outdim"
          A row vector that selects which dimensions of the calculation will appear in the output Y.  This is only useful when FCN returns an N-dimensional array in Format 1.  The default is to return all output dimensions.

     Programming Note: The property "outdim" can be used to save memory when the output of FCN has many dimensions, or when a wrapper to the base function that selects the desired outputs is too costly.  When memory is not an issue, the easiest way to select output dimensions is to first calculate the complete result with movfun and then filter that result with indexing.  If code complexity is not an issue then a wrapper can be created using anonymous functions.  For example, if basefcn is a function returning a K-dimensional row output, and only dimension D is desired, then the following wrapper could be used.

          FCN = @(x) basefcn (x)(:,columns(x) * (D-1) + (1:columns(x)));
          Y = movfun (@fcn, ...);

     movslice, prepad, postpad, permute, reshape

 -- SLCIDX =: movslice (N, WLEN)
 -- [SLCIDX, C, CPRE, CPOST, WIN] =: movslice (...)
     Generate indices to slice a vector of length N into windows of length WLEN.

     The input N must be a positive integer.

     The moving window length input WLEN can either be a scalar not equal to 1 or a 2-element array of integers.  For scalar values, if odd the window is symmetric and includes (WLEN - 1) / 2 elements on either side of the central element.  If WLEN is even the window is asymmetric and has WLEN/2 elements to the left of the central element and WLEN/2 - 1 elements to the right of the central element.  When WLEN is a 2-element array, [NB, NA], the window includes NB elements to the left of the current element and NA elements to the right of the current element.

     The output SLCIDX is an array of indices of the slices that fit fully within the vector, where each column is an individual slice as the window moves from left to right.  The slices have WLEN elements for scalar WLEN, or NB + NA + 1 elements for array valued WLEN.

     Optional output C is an row vector of window center positions where the window stays fully within the vector.

     Optional outputs CPRE and CPOST contain the vector elements at the start and end of the vector, respectively, that result in the window extending beyond the ends of the vector.

     Optional output WIN is a column vector with the same number of rows as SLCIDX that contains the moving window defined as a center relative position stencil.

     movfun

 -- Y =: movmad (X, WLEN)
 -- Y =: movmad (X, [NB, NA])
 -- Y =: movmad (..., DIM)
 -- Y =: movmad (..., "NANCOND")
 -- Y =: movmad (..., PROPERTY, VALUE)
     Calculate the moving mean absolute deviation over a sliding window of length WLEN on data X.

     If WLEN is a scalar, the function mad is applied to a moving window of length WLEN.  When WLEN is an odd number the window is symmetric and includes (WLEN - 1) / 2 elements oneither side of the central element.  For example, when calculating the output at index 5 with a window length of 3, movmad uses data elements [4, 5, 6].  If WLEN is an even number, the window is asymmetric and has WLEN/2 elements to the left of the central element and WLEN/2 - 1 elements to the right of thecentral element.  For example, when calculating the output at index 5 with a window length of 4, movmad uses data elements [3, 4, 5, 6].

     If WLEN is an array with two elements [NB, NA], the function is applied to a moving window -NB:NA.  This window includes NB number of elements before the current element and NA number of elements after the current element. The current element is always included.  For example, given WLEN = [3, 0], the data used to calculate index 5 is [2, 3, 4, 5].

     If the optional argument DIM is given, operate along this dimension.

     The optional string argument "NANCOND" controls whether NaN and NA values should be included ("includenan"), or excluded ("omitnan"), from the data passed to mad.  The default is "includenan".  Caution: the "omitnan" option is not yet implemented.

     The calculation can be controlled by specifying PROPERTY/VALUE pairs.  Valid properties are

     "Endpoints"

          This property controls how results are calculated at the boundaries (endpoints) of the window.  Possible values are:

          "shrink"  (default)
               The window is truncated at the beginning and end of the array to exclude elements for which there is no source data.  For example, with a window of length 3, Y(1) = mad (X(1:2)), and Y(end) = mad (X(end-1:end)).

          "discard"
               Any Y values that use a window extending beyond the original data array are deleted.  For example, with a 10-element data vector and a window of length 3, the output will contain only 8 elements.  The first element would require calculating the function over indices [0, 1, 2] and is therefore discarded.  The last element would require calculating the function over indices [9, 10, 11] and is therefore discarded.

          "fill"
               Any window elements outside the data array are replaced by NaN.  For example, with a window of length 3, Y(1) = mad ([NaN, X(1:2)]), and Y(end) = mad ([X(end-1:end), NaN]). This option usually results in Y having NaN values at the boundaries, although it is influenced by how mad handles NaN, and also by the property "nancond".

          USER_VALUE
               Any window elements outside the data array are replaced by the specified value USER_VALUE which must be a numeric scalar.  For example, with a window of length 3, Y(1) = mad ([USER_VALUE, X(1:2)]), and Y(end) = mad ([X(end-1:end), USER_VALUE]). A common choice for USER_VALUE is 0.

          "same"
               Any window elements outside the data array are replaced by the value of X at the boundary.  For example, with a window of length 3, Y(1) = mad ([X(1), X(1:2)]), and Y(end) = mad ([X(end-1:end), X(end)]).

          "periodic"
               The window is wrapped so that any missing data elements are taken from the other side of the data.  For example, with a window of length 3, Y(1) = mad ([X(end), X(1:2)]), and Y(end) = mad ([X(end-1:end), X(1)]).

     "SamplePoints"
          Caution: This option is not yet implemented.

     Programming Note: This function is a wrapper which calls movfun. For additional options and documentation, see movfun.

     movfun, movslice, movmax, movmean, movmedian, movmin, movprod, movstd, movsum, movvar

 -- Y =: movmax (X, WLEN)
 -- Y =: movmax (X, [NB, NA])
 -- Y =: movmax (..., DIM)
 -- Y =: movmax (..., "NANCOND")
 -- Y =: movmax (..., PROPERTY, VALUE)
     Calculate the moving maximum over a sliding window of length WLEN on data X.

     If WLEN is a scalar, the function max is applied to a moving window of length WLEN.  When WLEN is an odd number the window is symmetric and includes (WLEN - 1) / 2 elements on either side of the central element.  For example, when calculating the output at index 5 with a window length of 3, movmax uses data elements [4, 5, 6].  If WLEN is an even number, the window is asymmetric and has WLEN/2 elements to the left of the central element and WLEN/2 - 1 elements to the right of the central element.  For example, when calculating the output at index 5 with a window length of 4, movmax uses data elements [3, 4, 5, 6].

     If WLEN is an array with two elements [NB, NA], the function is applied to a moving window -NB:NA.  This window includes NB number of elements before the current element and NA number of elements after the current element. The current element is always included.  For example, given WLEN = [3, 0], the data used to calculate index 5 is [2, 3, 4, 5].

     If the optional argument DIM is given, operate along this dimension.

     The optional string argument "NANCOND" controls whether NaN and NA values should be included ("includenan"), or excluded ("omitnan"), from the data passed to max.  The default is "includenan".  Caution: the "omitnan" option is not yet implemented.

     The calculation can be controlled by specifying PROPERTY/VALUE pairs.  Valid properties are

     "Endpoints"

          This property controls how results are calculated at the boundaries (endpoints) of the window.  Possible values are:

          "shrink"  (default)
               The window is truncated at the beginning and end of the array to exclude elements for which there is no source data.  For example, with a window of length 3, Y(1) = max (X(1:2)), and Y(end) = max (X(end-1:end)).

          "discard"
               Any Y values that use a window extending beyond the original data array are deleted.  For example, with a 10-element data vector and a window of length 3, the output will contain only 8 elements.  The first element would require calculating the function over indices [0, 1, 2] and is therefore discarded.  The last element would require calculating the function over indices [9, 10, 11] and is therefore discarded.

          "fill"
               Any window elements outside the data array are replaced by NaN.  For example, with a window of length 3, Y(1) = max ([NaN, X(1:2)]), and Y(end) = max ([X(end-1:end), NaN]). This option usually results in Y having NaN values at the boundaries, although it is influenced by how max handles NaN, and also by the property "nancond".

          USER_VALUE
               Any window elements outside the data array are replaced by the specified value USER_VALUE which must be a numeric scalar.  For example, with a window of length 3, Y(1) = max ([USER_VALUE, X(1:2)]), and Y(end) = max ([X(end-1:end), USER_VALUE]). A common choice for USER_VALUE is 0.

          "same"
               Any window elements outside the data array are replaced by the value of X at the boundary.  For example, with a window of length 3, Y(1) = max ([X(1), X(1:2)]), and Y(end) = max ([X(end-1:end), X(end)]).

          "periodic"
               The window is wrapped so that any missing data elements are taken from the other side of the data.  For example, with a window of length 3, Y(1) = max ([X(end), X(1:2)]), and Y(end) = max ([X(end-1:end), X(1)]).

     "SamplePoints"
          Caution: This option is not yet implemented.

     Programming Note: This function is a wrapper which calls movfun. For additional options and documentation, see movfun.

     movfun, movslice, movmad, movmean, movmedian, movmin, movprod, movstd, movsum, movvar

 -- Y =: movmean (X, WLEN)
 -- Y =: movmean (X, [NB, NA])
 -- Y =: movmean (..., DIM)
 -- Y =: movmean (..., "NANCOND")
 -- Y =: movmean (..., PROPERTY, VALUE)
     Calculate the moving average over a sliding window of length WLEN on data X.

     If WLEN is a scalar, the function mean is applied to a moving window of length WLEN.  When WLEN is an odd number the window is symmetric and includes (WLEN - 1) / 2 elements on either side of the central element.  For example, when calculating the output at index 5 with a window length of 3, movmean uses data elements [4, 5, 6].  If WLEN is an even number, the window is asymmetric and has WLEN/2 elements to the left of the central element and WLEN/2 - 1 elements to the right of the central element.  For example, when calculating the output at index 5 with a window length of 4, movmean uses data elements [3, 4, 5, 6].

     If WLEN is an array with two elements [NB, NA], the function is applied to a moving window -NB:NA.  This window includes NB number of elements before the current element and NA number of elements after the current element. The current element is always included.  For example, given WLEN = [3, 0], the data used to calculate index 5 is [2, 3, 4, 5].

     If the optional argument DIM is given, operate along this dimension.

     The optional string argument "NANCOND" controls whether NaN and NA values should be included ("includenan"), or excluded ("omitnan"), from the data passed to mean.  The default is "includenan".  Caution: the "omitnan" option is not yet implemented.

     The calculation can be controlled by specifying PROPERTY/VALUE pairs.  Valid properties are

     "Endpoints"

          This property controls how results are calculated at the boundaries (endpoints) of the window.  Possible values are:

          "shrink"  (default)
               The window is truncated at the beginning and end of the array to exclude elements for which there is no source data.  For example, with a window of length 3, Y(1) = mean (X(1:2)), and Y(end) = mean (X(end-1:end)).

          "discard"
               Any Y values that use a window extending beyond the original data array are deleted.  For example, with a 10-element data vector and a window of length 3, the output will contain only 8 elements.  The first element would require calculating the function over indices [0, 1, 2] and is therefore discarded.  The last element would require calculating the function over indices [9, 10, 11] and is therefore discarded.

          "fill"
               Any window elements outside the data array are replaced by NaN.  For example, with a window of length 3, Y(1) = mean ([NaN, X(1:2)]), and Y(end) = mean ([X(end-1:end), NaN]). This option usually results in Y having NaN values at the boundaries, although it is influenced by how mean handles NaN, and also by the property "nancond".

          USER_VALUE
               Any window elements outside the data array are replaced by the specified value USER_VALUE which must be a numeric scalar.  For example, with a window of length 3, Y(1) = mean ([USER_VALUE, X(1:2)]), and Y(end) = mean ([X(end-1:end), USER_VALUE]). A common choice for USER_VALUE is 0.

          "same"
               Any window elements outside the data array are replaced by the value of X at the boundary.  For example, with a window of length 3, Y(1) = mean ([X(1), X(1:2)]), and Y(end) = mean ([X(end-1:end), X(end)]).

          "periodic"
               The window is wrapped so that any missing data elements are taken from the other side of the data.  For example, with a window of length 3, Y(1) = mean ([X(end), X(1:2)]), and Y(end) = mean ([X(end-1:end), X(1)]).

     "SamplePoints"
          Caution: This option is not yet implemented.

     Programming Note: This function is a wrapper which calls movfun. For additional options and documentation, see movfun.

     movfun, movslice, movmad, movmax, movmedian, movmin, movprod, movstd, movsum, movvar

 -- Y =: movmedian (X, WLEN)
 -- Y =: movmedian (X, [NB, NA])
 -- Y =: movmedian (..., DIM)
 -- Y =: movmedian (..., "NANCOND")
 -- Y =: movmedian (..., PROPERTY, VALUE)
     Calculate the moving median over a sliding window of length WLEN on data X.

     If WLEN is a scalar, the function movmedian is applied to a moving window of length WLEN.  When WLEN is an odd number the window is symmetric and includes (WLEN - 1) / 2 elements on either side of the central element.  For example, when calculating the output at index 5 with a window length of 3, movmedian uses data elements [4, 5, 6].  If WLEN is an even number, the window is asymmetric and has WLEN/2 elements to the left of the central element and WLEN/2 - 1 elements to the right of the central element.  For example, when calculating the output at index 5 with a window length of 4, movmedian uses data elements [3, 4, 5, 6].

     If WLEN is an array with two elements [NB, NA], the function is applied to a moving window -NB:NA.  This window includes NB number of elements before the current element and NA number of elements after the current element. The current element is always included.  For example, given WLEN = [3, 0], the data used to calculate index 5 is [2, 3, 4, 5].

     If the optional argument DIM is given, operate along this dimension.

     The optional string argument "NANCOND" controls whether NaN and NA values should be included ("includenan"), or excluded ("omitnan"), from the data passed to movmedian.  The default is "includenan".  Caution: the "omitnan" option is not yet implemented.

     The calculation can be controlled by specifying PROPERTY/VALUE pairs.  Valid properties are

     "Endpoints"

          This property controls how results are calculated at the boundaries (endpoints) of the window.  Possible values are:

          "shrink"  (default)
               The window is truncated at the beginning and end of the array to exclude elements for which there is no source data.  For example, with a window of length 3, Y(1) = movmedian (X(1:2)), and Y(end) = movmedian (X(end-1:end)).

          "discard"
               Any Y values that use a window extending beyond the original data array are deleted.  For example, with a 10-element data vector and a window of length 3, the output will contain only 8 elements.  The first element would require calculating the function over indices [0, 1, 2] and is therefore discarded.  The last element would require calculating the function over indices [9, 10, 11] and is therefore discarded.

          "fill"
               Any window elements outside the data array are replaced by NaN.  For example, with a window of length 3, Y(1) = movmedian ([NaN, X(1:2)]), and Y(end) = movmedian ([X(end-1:end), NaN]). This option usually results in Y having NaN values at the boundaries, although it is influenced by how movmedian handles NaN, and also by the property "nancond".

          USER_VALUE
               Any window elements outside the data array are replaced by the specified value USER_VALUE which must be a numeric scalar.  For example, with a window of length 3, Y(1) = movmedian ([USER_VALUE, X(1:2)]), and Y(end) = movmedian ([X(end-1:end), USER_VALUE]). A common choice for USER_VALUE is 0.

          "same"
               Any window elements outside the data array are replaced by the value of X at the boundary.  For example, with a window of length 3, Y(1) = movmedian ([X(1), X(1:2)]), and Y(end) = movmedian ([X(end-1:end), X(end)]).

          "periodic"
               The window is wrapped so that any missing data elements are taken from the other side of the data.  For example, with a window of length 3, Y(1) = movmedian ([X(end), X(1:2)]), and Y(end) = movmedian ([X(end-1:end), X(1)]).

     "SamplePoints"
          Caution: This option is not yet implemented.

     Programming Note: This function is a wrapper which calls movfun. For additional options and documentation, see movfun.

     movfun, movslice, movmad, movmax, movmean, movmin, movprod, movstd, movsum, movvar

 -- Y =: movmin (X, WLEN)
 -- Y =: movmin (X, [NB, NA])
 -- Y =: movmin (..., DIM)
 -- Y =: movmin (..., "NANCOND")
 -- Y =: movmin (..., PROPERTY, VALUE)
     Calculate the moving minimum over a sliding window of length WLEN on data X.

     If WLEN is a scalar, the function min is applied to a moving window of length WLEN.  When WLEN is an odd number the window is symmetric and includes (WLEN - 1) / 2 elements on either side of the central element.  For example, when calculating the output at index 5 with a window length of 3, movmin uses data elements [4, 5, 6].  If WLEN is an even number, the window is asymmetric and has WLEN/2 elements to the left of the central element and WLEN/2 - 1 elements to the right of the central element.  For example, when calculating the output at index 5 with a window length of 4, movmin uses data elements [3, 4, 5, 6].

     If WLEN is an array with two elements [NB, NA], the function is applied to a moving window -NB:NA.  This window includes NB number of elements before the current element and NA number of elements after the current element. The current element is always included.  For example, given WLEN = [3, 0], the data used to calculate index 5 is [2, 3, 4, 5].

     If the optional argument DIM is given, operate along this dimension.

     The optional string argument "NANCOND" controls whether NaN and NA values should be included ("includenan"), or excluded ("omitnan"), from the data passed to min.  The default is "includenan".  Caution: the "omitnan" option is not yet implemented.

     The calculation can be controlled by specifying PROPERTY/VALUE pairs.  Valid properties are

     "Endpoints"

          This property controls how results are calculated at the boundaries (endpoints) of the window.  Possible values are:

          "shrink"  (default)
               The window is truncated at the beginning and end of the array to exclude elements for which there is no source data.  For example, with a window of length 3, Y(1) = min (X(1:2)), and Y(end) = min (X(end-1:end)).

          "discard"
               Any Y values that use a window extending beyond the original data array are deleted.  For example, with a 10-element data vector and a window of length 3, the output will contain only 8 elements.  The first element would require calculating the function over indices [0, 1, 2] and is therefore discarded.  The last element would require calculating the function over indices [9, 10, 11] and is therefore discarded.

          "fill"
               Any window elements outside the data array are replaced by NaN.  For example, with a window of length 3, Y(1) = min ([NaN, X(1:2)]), and Y(end) = min ([X(end-1:end), NaN]). This option usually results in Y having NaN values at the boundaries, although it is influenced by how min handles NaN, and also by the property "nancond".

          USER_VALUE
               Any window elements outside the data array are replaced by the specified value USER_VALUE which must be a numeric scalar.  For example, with a window of length 3, Y(1) = min ([USER_VALUE, X(1:2)]), and Y(end) = min ([X(end-1:end), USER_VALUE]). A common choice for USER_VALUE is 0.

          "same"
               Any window elements outside the data array are replaced by the value of X at the boundary.  For example, with a window of length 3, Y(1) = min ([X(1), X(1:2)]), and Y(end) = min ([X(end-1:end), X(end)]).

          "periodic"
               The window is wrapped so that any missing data elements are taken from the other side of the data.  For example, with a window of length 3, Y(1) = min ([X(end), X(1:2)]), and Y(end) = min ([X(end-1:end), X(1)]).

     "SamplePoints"
          Caution: This option is not yet implemented.

     Programming Note: This function is a wrapper which calls movfun. For additional options and documentation, see movfun.

     movfun, movslice, movmad, movmax, movmean, movmedian, movprod, movstd, movsum, movvar

 -- Y =: movprod (X, WLEN)
 -- Y =: movprod (X, [NB, NA])
 -- Y =: movprod (..., DIM)
 -- Y =: movprod (..., "NANCOND")
 -- Y =: movprod (..., PROPERTY, VALUE)
     Calculate the moving product over a sliding window of length WLEN on data X.

     If WLEN is a scalar, the function movprod is applied to a moving window of length WLEN.  When WLEN is an odd number the window is symmetric and includes (WLEN - 1) / 2 elements on either side of the central element.  For example, when calculating the output at index 5 with a window length of 3, movprod uses data elements [4, 5, 6].  If WLEN is an even number, the window is asymmetric and has WLEN/2 elements to the left of the central element and WLEN/2 - 1 elements to the right of the central element.  For example, when calculating the output at index 5 with a window length of 4, movprod uses data elements [3, 4, 5, 6].

     If WLEN is an array with two elements [NB, NA], the function is applied to a moving window -NB:NA.  This window includes NB number of elements before the current element and NA number of elements after the current element. The current element is always included.  For example, given WLEN = [3, 0], the data used to calculate index 5 is [2, 3, 4, 5].

     If the optional argument DIM is given, operate along this dimension.

     The optional string argument "NANCOND" controls whether NaN and NA values should be included ("includenan"), or excluded ("omitnan"), from the data passed to movprod. The default is "includenan".  Caution: the "omitnan" option is not yet implemented.

     The calculation can be controlled by specifying PROPERTY/VALUE pairs.  Valid properties are

     "Endpoints"

          This property controls how results are calculated at the boundaries (endpoints) of the window.  Possible values are:

          "shrink"  (default)
               The window is truncated at the beginning and end of the array to exclude elements for which there is no source data.  For example, with a window of length 3, Y(1) = movprod (X(1:2)), and Y(end) = movprod (X(end-1:end)).

          "discard"
               Any Y values that use a window extending beyond the original data array are deleted.  For example, with a 10-element data vector and a window of length 3, the output will contain only 8 elements.  The first element would require calculating the function over indices [0, 1, 2] and is therefore discarded.  The last element would require calculating the function over indices [9, 10, 11] and is therefore discarded.

          "fill"
               Any window elements outside the data array are replaced by NaN.  For example, with a window of length 3, Y(1) = movprod ([NaN, X(1:2)]), and Y(end) = movprod ([X(end-1:end), NaN]). This option usually results in Y having NaN values at the boundaries, although it is influenced by how movprod handles NaN, and also by the property "nancond".

          USER_VALUE
               Any window elements outside the data array are replaced by the specified value USER_VALUE which must be a numeric scalar.  For example, with a window of length 3, Y(1) = movprod ([USER_VALUE, X(1:2)]), and Y(end) = movprod ([X(end-1:end), USER_VALUE]). A common choice for USER_VALUE is 0.

          "same"
               Any window elements outside the data array are replaced by the value of X at the boundary.  For example, with a window of length 3, Y(1) = movprod ([X(1), X(1:2)]), and Y(end) = movprod ([X(end-1:end), X(end)]).

          "periodic"
               The window is wrapped so that any missing data elements are taken from the other side of the data.  For example, with a window of length 3, Y(1) = movprod ([X(end), X(1:2)]), and Y(end) = movprod ([X(end-1:end), X(1)]).

     "SamplePoints"
          Caution: This option is not yet implemented.

     Programming Note: This function is a wrapper which calls movfun. For additional options and documentation, see movfun.

     movfun, movslice, movmad, movmax, movmean, movmedian, movmin, movstd, movsum, movvar

 -- Y =: movstd (X, WLEN)
 -- Y =: movstd (X, [NB, NA])
 -- Y =: movstd (..., OPT)
 -- Y =: movstd (..., OPT, DIM)
 -- Y =: movstd (..., "NANCOND")
 -- Y =: movstd (..., PROPERTY, VALUE)
     Calculate the moving standard deviation over a sliding window of length WLEN on data X.

     If WLEN is a scalar, the function movstd is applied to a moving window of length WLEN.  When WLEN is an odd number the window is symmetric and includes (WLEN - 1) / 2 elements on either side of the central element.  For example, when calculating the output at index 5 with a window length of 3, movstd uses data elements [4, 5, 6].  If WLEN is an even number, the window is asymmetric and has WLEN/2 elements to the left of the central element and WLEN/2 - 1 elements to the right of the central element.  For example, when calculating the output at index 5 with a window length of 4, movstd uses data elements [3, 4, 5, 6].

     If WLEN is an array with two elements [NB, NA], the function is applied to a moving window -NB:NA.  This window includes NB number of elements before the current element and NA number of elements after the current element. The current element is always included.  For example, given WLEN = [3, 0], the data used to calculate index 5 is [2, 3, 4, 5].

     The optional argument OPT determines the type of normalization to use. Valid values are

     0:
          normalize with N-1, provides the square root of the best unbiased estimator of the variance [default]

     1:
          normalize with N, this provides the square root of the second moment around the mean

     If the optional argument DIM is given, operate along this dimension. The normalization argument OPT must be given before the dimension. To use the default value for OPT you may pass an empty input argument [].

     The optional string argument "NANCOND" controls whether NaN and NA values should be included ("includenan"), or excluded ("omitnan"), from the data passed to movstd.  The default is "includenan".  Caution: the "omitnan" option is not yet implemented.

     The calculation can be controlled by specifying PROPERTY/VALUE pairs.  Valid properties are

     "Endpoints"

          This property controls how results are calculated at the boundaries (endpoints) of the window.  Possible values are:

          "shrink"  (default)
               The window is truncated at the beginning and end of the array to exclude elements for which there is no source data.  For example, with a window of length 3, Y(1) = movstd (X(1:2)), and Y(end) = movstd (X(end-1:end)).

          "discard"
               Any Y values that use a window extending beyond the original data array are deleted.  For example, with a 10-element data vector and a window of length 3, the output will contain only 8 elements.  The first element would require calculating the function over indices [0, 1, 2] and is therefore discarded.  The last element would require calculating the function over indices [9, 10, 11] and is therefore discarded.

          "fill"
               Any window elements outside the data array are replaced by NaN.  For example, with a window of length 3, Y(1) = movstd ([NaN, X(1:2)]), and Y(end) = movstd ([X(end-1:end), NaN]). This option usually results in Y having NaN values at the boundaries, although it is influenced by how movstd handles NaN, and also by the property "nancond".

          USER_VALUE
               Any window elements outside the data array are replaced by the specified value USER_VALUE which must be a numeric scalar.  For example, with a window of length 3, Y(1) = movstd ([USER_VALUE, X(1:2)]), and Y(end) = movstd ([X(end-1:end), USER_VALUE]). A common choice for USER_VALUE is 0.

          "same"
               Any window elements outside the data array are replaced by the value of X at the boundary.  For example, with a window of length 3, Y(1) = movstd ([X(1), X(1:2)]), and Y(end) = movstd ([X(end-1:end), X(end)]).

          "periodic"
               The window is wrapped so that any missing data elements are taken from the other side of the data.  For example, with a window of length 3, Y(1) = movstd ([X(end), X(1:2)]), and Y(end) = movstd ([X(end-1:end), X(1)]).

     "SamplePoints"
          Caution: This option is not yet implemented.

     Programming Note: This function is a wrapper which calls movfun. For additional options and documentation, see movfun.

     movfun, movslice, movmad, movmax, movmean, movmedian, movmin, movprod, movsum, movvar

 -- Y =: movsum (X, WLEN)
 -- Y =: movsum (X, [NB, NA])
 -- Y =: movsum (..., DIM)
 -- Y =: movsum (..., "NANCOND")
 -- Y =: movsum (..., PROPERTY, VALUE)
     Calculate the moving sum over a sliding window of length WLEN on data X.

     If WLEN is a scalar, the function movsum is applied to a moving window of length WLEN.  When WLEN is an odd number the window is symmetric and includes (WLEN - 1) / 2 elements on either side of the central element.  For example, when calculating the output at index 5 with a window length of 3, movsum uses data elements [4, 5, 6].  If WLEN is an even number, the window is asymmetric and has WLEN/2 elements to the left of the central element and WLEN/2 - 1 elements to the right of the central element.  For example, when calculating the output at index 5 with a window length of 4, movsum uses data elements [3, 4, 5, 6].

     If WLEN is an array with two elements [NB, NA], the function is applied to a moving window -NB:NA.  This window includes NB number of elements before the current element and NA number of elements after the current element. The current element is always included.  For example, given WLEN = [3, 0], the data used to calculate index 5 is [2, 3, 4, 5].

     If the optional argument DIM is given, operate along this dimension.

     The optional string argument "NANCOND" controls whether NaN and NA values should be included ("includenan"), or excluded ("omitnan"), from the data passed to movsum.  The default is "includenan".  Caution: the "omitnan" option is not yet implemented.

     The calculation can be controlled by specifying PROPERTY/VALUE pairs.  Valid properties are

     "Endpoints"

          This property controls how results are calculated at the boundaries (endpoints) of the window.  Possible values are:

          "shrink"  (default)
               The window is truncated at the beginning and end of the array to exclude elements for which there is no source data.  For example, with a window of length 3, Y(1) = movsum (X(1:2)), and Y(end) = movsum (X(end-1:end)).

          "discard"
               Any Y values that use a window extending beyond the original data array are deleted.  For example, with a 10-element data vector and a window of length 3, the output will contain only 8 elements.  The first element would require calculating the function over indices [0, 1, 2] and is therefore discarded.  The last element would require calculating the function over indices [9, 10, 11] and is therefore discarded.

          "fill"
               Any window elements outside the data array are replaced by NaN.  For example, with a window of length 3, Y(1) = movsum ([NaN, X(1:2)]), and Y(end) = movsum ([X(end-1:end), NaN]). This option usually results in Y having NaN values at the boundaries, although it is influenced by how movsum handles NaN, and also by the property "nancond".

          USER_VALUE
               Any window elements outside the data array are replaced by the specified value USER_VALUE which must be a numeric scalar.  For example, with a window of length 3, Y(1) = movsum ([USER_VALUE, X(1:2)]), and Y(end) = movsum ([X(end-1:end), USER_VALUE]). A common choice for USER_VALUE is 0.

          "same"
               Any window elements outside the data array are replaced by the value of X at the boundary.  For example, with a window of length 3, Y(1) = movsum ([X(1), X(1:2)]), and Y(end) = movsum ([X(end-1:end), X(end)]).

          "periodic"
               The window is wrapped so that any missing data elements are taken from the other side of the data.  For example, with a window of length 3, Y(1) = movsum ([X(end), X(1:2)]), and Y(end) = movsum ([X(end-1:end), X(1)]).

     "SamplePoints"
          Caution: This option is not yet implemented.

     Programming Note: This function is a wrapper which calls movfun. For additional options and documentation, see movfun.

     movfun, movslice, movmad, movmax, movmean, movmedian, movmin, movprod, movstd, movvar

 -- Y =: movvar (X, WLEN)
 -- Y =: movvar (X, [NB, NA])
 -- Y =: movvar (..., OPT)
 -- Y =: movvar (..., OPT, DIM)
 -- Y =: movvar (..., "NANCOND")
 -- Y =: movvar (..., PROPERTY, VALUE)
     Calculate the moving variance over a sliding window of length WLEN on data X.

     If WLEN is a scalar, the function var is applied to a moving window of length WLEN.  When WLEN is an odd number the window is symmetric and includes (WLEN - 1) / 2 elements on either side of the central element.  For example, when calculating the output at index 5 with a window length of 3, movvar uses data elements [4, 5, 6].  If WLEN is an even number, the window is asymmetric and has WLEN/2 elements to the left of the central element and WLEN/2 - 1 elements to the right of the central element.  For example, when calculating the output at index 5 with a window length of 4, movvar uses data elements [3, 4, 5, 6].

     If WLEN is an array with two elements [NB, NA], the function is applied to a moving window -NB:NA.  This window includes NB number of elements before the current element and NA number of elements after the current element. The current element is always included.  For example, given WLEN = [3, 0], the data used to calculate index 5 is [2, 3, 4, 5].

     The optional argument OPT determines the type of normalization to use. Valid values are:

     0:
          normalize with N-1, provides the best unbiased estimator of the variance [default]

     1:
          normalizes with N, this provides the second moment around the mean

     If the optional argument DIM is given, operate along this dimension. The normalization argument OPT must be given before the dimension. To use the default value for OPT you may pass an empty input argument [].

     The optional string argument "NANCOND" controls whether NaN and NA values should be included ("includenan"), or excluded ("omitnan"), from the data passed to var.  The default is "includenan".  Caution: the "omitnan" option is not yet implemented.

     The calculation can be controlled by specifying PROPERTY/VALUE pairs.  Valid properties are

     "Endpoints"

          This property controls how results are calculated at the boundaries (endpoints) of the window.  Possible values are:

          "shrink"  (default)
               The window is truncated at the beginning and end of the array to exclude elements for which there is no source data.  For example, with a window of length 3, Y(1) = var (X(1:2)), and Y(end) = var (X(end-1:end)).

          "discard"
               Any Y values that use a window extending beyond the original data array are deleted.  For example, with a 10-element data vector and a window of length 3, the output will contain only 8 elements.  The first element would require calculating the function over indices [0, 1, 2] and is therefore discarded.  The last element would require calculating the function over indices [9, 10, 11] and is therefore discarded.

          "fill"
               Any window elements outside the data array are replaced by NaN.  For example, with a window of length 3, Y(1) = var ([NaN, X(1:2)]), and Y(end) = var ([X(end-1:end), NaN]). This option usually results in Y having NaN values at the boundaries, although it is influenced by how var handles NaN, and also by the property "nancond".

          USER_VALUE
               Any window elements outside the data array are replaced by the specified value USER_VALUE which must be a numeric scalar.  For example, with a window of length 3, Y(1) = var ([USER_VALUE, X(1:2)]), and Y(end) = var ([X(end-1:end), USER_VALUE]). A common choice for USER_VALUE is 0.

          "same"
               Any window elements outside the data array are replaced by the value of X at the boundary.  For example, with a window of length 3, Y(1) = var ([X(1), X(1:2)]), and Y(end) = var ([X(end-1:end), X(end)]).

          "periodic"
               The window is wrapped so that any missing data elements are taken from the other side of the data.  For example, with a window of length 3, Y(1) = var ([X(end), X(1:2)]), and Y(end) = var ([X(end-1:end), X(1)]).

     "SamplePoints"
          Caution: This option is not yet implemented.

     Programming Note: This function is a wrapper which calls movfun. For additional options and documentation, see movfun.

     movfun, movslice, movmad, movmax, movmean, movmedian, movmin, movprod, movstd, movsum

Basic Statistical Functions
===========================

Octave supports various helpful statistical functions.  Many are useful as initial steps to prepare a data set for further analysis.  Others provide different measures from those of the basic descriptive statistics.

 -- Y =: center (X)
 -- Y =: center (X, DIM)
     Center data by subtracting its mean.

     If X is a vector, subtract its mean.

     If X is a matrix, do the above for each column.

     If the optional argument DIM is given, operate along this dimension.

     Programming Note: center has obvious application for normalizing statistical data.  It is also useful for improving the precision of general numerical calculations.  Whenever there is a large value that is common to a batch of data, the mean can be subtracted off, the calculation performed, and then the mean added back to obtain the final answer. zscore

 -- Z =: zscore (X)
 -- Z =: zscore (X, OPT)
 -- Z =: zscore (X, OPT, DIM)
 -- [Z, MU, SIGMA] =: zscore (...)
     Compute the Z score of X.

     If X is a vector, subtract its mean and divide by its standard deviation.  If the standard deviation is zero, divide by 1 instead.

     The optional parameter OPT determines the normalization to use when computing the standard deviation and has the same definition as the corresponding parameter for std.

     If X is a matrix, calculate along the first non-singleton dimension. If the third optional argument DIM is given, operate along this dimension.

     The optional outputs MU and SIGMA contain the mean and standard deviation.

     mean, std, center

 -- Z =: normalize (X)
 -- Z =: normalize (X, DIM)
 -- Z =: normalize (..., METHOD)
 -- Z =: normalize (..., METHOD, OPTION)
 -- Z =: normalize (..., SCALE, SCALEOPTION, CENTER, CENTEROPTION)
 -- [Z, C, S] =: normalize (...)
     Return a normalization of the data in X using one of several available scaling and centering methods.

     normalize by default will return the zscore of X, defined as the number of standard deviations each element is from the mean of X.  This is equivalent to centering at the mean of the data and scaling by the standard deviation.

     The returned value Z will have the same size as X.  The optional return variables C and S are the centering and scaling factors used in the normalization such that:

            Z = (X - C) ./ S

     If X is a vector, normalize will operate on the data in X.

     If X is a matrix, normalize will operate independently on each column in X.

     If X is an N-dimensional array, normalize will operate independently on the first non-singleton dimension in X.

     If the optional second argument DIM is given, operate along this dimension.

     normalize ignores NaN values is X similar to the behavior of the omitnan option in std, mean, and median.

     The optional inputs METHOD and OPTION can be used to specify the type of normalization performed on X.  Note that only the scale and center options may be specified together using any of the methods defined below.  Valid normalization methods are:

     zscore
          (Default) Normalizes the elements in X to the scaled distance from a central value.  Valid Options:

          std
               (Default) Data is centered at mean (X) and scaled by the standard deviation.

          robust
               Data is centered at median (X) and scaled by the median absolute deviation.

     norm
          Z is the general vector norm of X, with OPTION being the normalization factor P that determines the vector norm type according to:

                 Z = [sum (abs (X) .^ P)] ^ (1/P)

          P can be any positive scalar, specific values being:

          P = 1
               X is normalized by sum (abs (X)).

          P = 2
               (Default) X is normalized by the Euclidian norm, or vector magnitude, of the elements.

          P = Inf
               X is normalized by max (abs (X)).

     scale
          X is scaled by a factor determined by OPTION, which can be a numeric scalar or one of the following:

          std
               (Default) X is scaled by its standard deviation.

          mad
               X is scaled by its median absolute deviation.

          first
               X is scaled by its first element.

          iqr
               X is scaled by its interquartile range.

     range
          X is scaled to fit the range specified by OPTION as a two element scalar row vector.  The default range is [0, 1].

     center
          X is shifted by an amount determined by OPTION, which can be a numeric scalar or one of the following:

          mean
               (Default) X is shifted by mean (X).

          median
               X is shifted by median (X).

     medianiqr
          X is shifted by median (X) and scaled by the interquartile range.

     Known MATLAB incompatibilities:

      1. The option DataVariables is not yet implemented for Table class X inputs.

     zscore, iqr, norm, rescale, std, median, mean, mad

 -- N =: histc (X, EDGES)
 -- N =: histc (X, EDGES, DIM)
 -- [N, IDX] =: histc (...)
     Compute histogram counts.

     When X is a vector, the function counts the number of elements of X that fall in the histogram bins defined by EDGES.  This must be a vector of monotonically increasing values that define the edges of the histogram bins. N(k) contains the number of elements in X for which EDGES(k) <= X < EDGES(k+1). The final element of N contains the number of elements of X exactly equal to the last element of EDGES.

     When X is an N-dimensional array, the computation is carried out along dimension DIM.  If not specified DIM defaults to the first non-singleton dimension.

     When a second output argument is requested an index matrix is also returned. The IDX matrix has the same size as X.  Each element of IDX contains the index of the histogram bin in which the corresponding element of X was counted. hist

unique function documented at unique is often useful for statistics.

 -- C =: nchoosek (N, K)
 -- C =: nchoosek (SET, K)

     Compute the binomial coefficient of N or list all possible combinations of a SET of items.

     If N is a scalar then calculate the binomial coefficient of N and K which is defined as

           /   \
           | n |    n (n-1) (n-2) ... (n-k+1)       n!
           |   |  = ------------------------- =  ---------
           | k |               k!                k! (n-k)!
           \   /

     This is the number of combinations of N items taken in groups of size K.

     If the first argument is a vector, SET, then generate all combinations of the elements of SET, taken K at a time, with one row per combination.  The result C has K columns and nchoosek (length (SET), K) rows.

     For example:

     How many ways can three items be grouped into pairs?

          nchoosek (3, 2)
             => 3

     What are the possible pairs?

          nchoosek (1:3, 2)
             =>  1   2
                 1   3
                 2   3

     Programming Note: When calculating the binomial coefficient nchoosek works only for non-negative, integer arguments.  Use bincoeff for non-integer and negative scalar arguments, or for computing many binomial coefficients at once with vector inputs for N or K.

     bincoeff, perms

 -- P =: perms (V)
 -- P =: perms (V, "unique")
     Generate all permutations of vector V with one row per permutation.

     Results are returned in reverse lexicographic order if V is in ascending order.  If V is in a different permutation, then the result is permuted that way too.  Consequently, an input in descending order yields a result in normal lexicographic order.  The result has size factorial (N) * N, where N is the length of V. Any repeated elements are included in the output.

     If the optional argument "unique" is given then only unique permutations are returned, using less memory and taking less time than calling unique (perms (V), "rows").

     Example 1

          perms ([1, 2, 3])
          =>
          3   2   1
          3   1   2
          2   3   1
          2   1   3
          1   3   2
          1   2   3

     Example 2

          perms ([1, 1, 2, 2], "unique")
          =>
          2   2   1   1
          2   1   2   1
          2   1   1   2
          1   2   2   1
          1   2   1   2
          1   1   2   2

     Programming Note: If the "unique" option is not used, the length of V should be no more than 10-12 to limit memory consumption.  Even with "unique", there should be no more than 10-12 unique elements in V. permute, randperm, nchoosek

 -- Y =: ranks (X)
 -- Y =: ranks (X, DIM)
 -- Y =: ranks (X, DIM, RTYPE)
     Return the ranks (in the sense of order statistics) of X along the first non-singleton dimension adjusted for ties.

     If the optional DIM argument is given, operate along this dimension.

     The optional parameter RTYPE determines how ties are handled.  All examples below assume an input of [ 1, 2, 2, 4 ].

     0 or "fractional" (default) for fractional ranking (1, 2.5,
          2.5, 4);

     1 or "competition" for competition ranking (1, 2, 2, 4);

     2 or "modified" for modified competition ranking (1, 3, 3, 4);

     3 or "ordinal" for ordinal ranking (1, 2, 3, 4);

     4 or "dense" for dense ranking (1, 2, 2, 3).

     spearman, kendall

 -- CNT =: run_count (X, N)
 -- CNT =: run_count (X, N, DIM)
     Count the upward runs along the first non-singleton dimension of X of length 1, 2, ..., N-1 and greater than or equal to N.

     If the optional argument DIM is given then operate along this dimension. runlength

 -- count =: runlength (X)
 -- [count, value] =: runlength (X)
     Find the lengths of all sequences of common values.

     COUNT is a vector with the lengths of each repeated value.

     The optional output VALUE contains the value that was repeated in the sequence.

          runlength ([2, 2, 0, 4, 4, 4, 0, 1, 1, 1, 1])
          =>   2   1   3   1   4

     run_count

Correlation and Regression Analysis
===================================

 -- C =: cov (X)
 -- C =: cov (X, Y)
 -- C =: cov (..., OPT)
 -- C =: cov (..., NANFLAG)
     Compute the covariance matrix.

     The covariance between two variable vectors A and  B is calculated as:

          cov (A,B) = 1/(N-1) * SUM_i (A(i) - mean (A)) * (B(i) - mean (B))

     where N is the length of the vectors A and B.

     If called with one argument, compute cov (X, X).  If X is a vector, this is the scalar variance of X.  If X is a matrix, each row of X is treated as an observation, and each column as a variable, and the (I, J)-th entry of cov (X) is the covariance between the I-th and J-th columns in X.  If X has dimensions n x m, the output C will be a m x m square covariance matrix.

     If called with two arguments, compute cov (X, Y), the covariance between two random variables X and Y.  X and Y must have the same number of elements, and will be treated as vectors with the covariance computed as cov (X(:), Y(:)).  The output will be a 2 x 2 covariance matrix.

     The optional argument OPT determines the type of normalization to use.  Valid values are

     0 [default]:
          Normalize with N-1.  This provides the best unbiased estimator of the covariance.

     1:
          Normalize with N.  This provides the second moment around the mean.  OPT is set to 1 for N = 1.

     The optional argument NANFLAG must appear last in the argument list and controls how NaN values are handled by cov.  The three valid values are:

     includenan [default]:
          Leave NaN values in X and Y.  Output will follow the normal rules for handling NaN values in arithmetic operations.

     omitrows:
          Rows containing NaN values are trimmed from both X and Y prior to calculating the covariance.  A NaN in one variable will remove that row from both X and Y.

     partialrows:
          Rows containing NaN values are ignored from both X and Y independently for each I-th and J-th covariance calculation.  This may result in a different number of observations, N, being used to calculated each element of the covariance matrix.

     Compatibility Note: Previous versions of cov treated rows X and Y as multivariate random variables.  This version attempts to maintain full compatibility with MATLAB by treating X and Y as two univariate distributions regardless of shape, resulting in a 2x2 output matrix.  Code relying on Octave's previous definition will need to be modified when running this newer version of cov.  The previous behavior can be obtained by using the NaN package's covm function as covm (X, Y, "D"). corr

 -- R =: corr (X)
 -- R =: corr (X, Y)
     Compute matrix of correlation coefficients.

     If each row of X and Y is an observation and each column is a variable, then the (I, J)-th entry of corr (X, Y) is the correlation between the I-th variable in X and the J-th variable in Y. X and Y must have the same number of rows (observations).

          corr (X,Y) = cov (X,Y) / (std (X) * std (Y))

     If called with one argument, compute corr (X, X), the correlation between the columns of X. cov

 -- R =: corrcoef (X)
 -- R =: corrcoef (X, Y)
 -- R =: corrcoef (..., PARAM, VALUE, ...)
 -- [R, P] =: corrcoef (...)
 -- [R, P, LCI, HCI] =: corrcoef (...)
     Compute a matrix of correlation coefficients.

     X is an array where each column contains a variable and each row is an observation.

     If a second input Y (of the same size as X) is given then calculate the correlation coefficients between X and Y.

     PARAM, VALUE are optional pairs of parameters and values which modify the calculation.  Valid options are:

     "alpha"
          Confidence level used for the bounds of the confidence interval, LCI and HCI.  Default is 0.05, i.e., 95% confidence interval.

     "rows"
          Determine processing of NaN values.  Acceptable values are "all", "complete", and "pairwise".  Default is "all". With "complete", only the rows without NaN values are considered. With "pairwise", the selection of NaN-free rows is made for each pair of variables.

     Output R is a matrix of Pearson's product moment correlation coefficients for each pair of variables.

     Output P is a matrix of pair-wise p-values testing for the null hypothesis of a correlation coefficient of zero.

     Outputs LCI and HCI are matrices containing, respectively, the lower and higher bounds of the 95% confidence interval of each correlation coefficient. corr, cov, std

 -- RHO =: spearman (X)
 -- RHO =: spearman (X, Y)
     Compute Spearman's rank correlation coefficient RHO.

     For two data vectors X and Y, Spearman's RHO is the correlation coefficient of the ranks of X and Y.

     If X and Y are drawn from independent distributions, RHO has zero mean and variance 1 / (N - 1), where N is the length of the X and Y vectors, and is asymptotically normally distributed.

     spearman (X) is equivalent to spearman (X, X). ranks, kendall

 -- TAU =: kendall (X)
 -- TAU =: kendall (X, Y)
     Compute Kendall's TAU.

     For two data vectors X, Y of common length N, Kendall's TAU is the correlation of the signs of all rank differences of X and Y; i.e., if both X and Y have distinct entries, then

                   1
          TAU = -------   SUM sign (Q(i) - Q(j)) * sign (R(i) - R(j))
                N (N-1)   i,j

     in which the Q(i) and R(i) are the ranks of X and Y, respectively.

     If X and Y are drawn from independent distributions, Kendall's TAU is asymptotically normal with mean 0 and variance (2 * (2N+5)) / (9 * N * (N-1)).

     kendall (X) is equivalent to kendall (X, X). ranks, spearman

Distributions
=============

Octave has functions for computing the Probability Density Function (PDF), the Cumulative Distribution function (CDF), and the quantile (the inverse of the CDF) for arbitrary user-defined distributions (discrete) and for experimental data (empirical).

The following table summarizes the supported distributions (in alphabetical order).

Distribution
PDF
CDF
Quantile
Univariate Discrete Distribution
discrete_pdf
discrete_cdf
discrete_inv
Empirical Distribution
empirical_pdf
empirical_cdf
empirical_inv

 -- PDF =: discrete_pdf (X, V, P)
     For each element of X, compute the probability density function (PDF) at X of a univariate discrete distribution which assumes the values in V with probabilities P.

 -- CDF =: discrete_cdf (X, V, P)
     For each element of X, compute the cumulative distribution function (CDF) at X of a univariate discrete distribution which assumes the values in V with probabilities P.

 -- Q =: discrete_inv (X, V, P)
     For each element of X, compute the quantile (the inverse of the CDF) at X of the univariate distribution which assumes the values in V with probabilities P.

 -- PDF =: empirical_pdf (X, DATA)
     For each element of X, compute the probability density function (PDF) at X of the empirical distribution obtained from the univariate sample DATA.

 -- CDF =: empirical_cdf (X, DATA)
     For each element of X, compute the cumulative distribution function (CDF) at X of the empirical distribution obtained from the univariate sample DATA.

 -- Q =: empirical_inv (X, DATA)
     For each element of X, compute the quantile (the inverse of the CDF) at X of the empirical distribution obtained from the univariate sample DATA.

Random Number Generation
========================

Octave can generate random numbers from a large number of distributions.  The random number generators are based on the random number generators described in Special Utility Matrices.

The following table summarizes the available random number generators (in alphabetical order).

Distribution                    Function
Univariate Discrete Distribution    discrete_rnd
Empirical Distribution              empirical_rnd
Exponential Distribution            rande
Gamma Distribution                  randg
Poisson Distribution                randp
Standard Normal Distribution        randn
Uniform Distribution                rand
Uniform Distribution (integers)     randi

 -- RND =: discrete_rnd (V, P)
 -- RND =: discrete_rnd (V, P, R)
 -- RND =: discrete_rnd (V, P, R, C, ...)
 -- RND =: discrete_rnd (V, P, [SZ])
     Return a matrix of random samples from the univariate distribution which assumes the values in V with probabilities P.

     When called with a single size argument, return a square matrix with the dimension specified.  When called with more than one scalar argument the first two arguments are taken as the number of rows and columns and any further arguments specify additional matrix dimensions.  The size may also be specified with a vector of dimensions SZ.

     If no size arguments are given then the result matrix is the common size of V and P.

 -- RND =: empirical_rnd (DATA)
 -- RND =: empirical_rnd (DATA, R)
 -- RND =: empirical_rnd (DATA, R, C, ...)
 -- RND =: empirical_rnd (DATA, [SZ])
     Return a matrix of random samples from the empirical distribution obtained from the univariate sample DATA.

     When called with a single size argument, return a square matrix with the dimension specified.  When called with more than one scalar argument the first two arguments are taken as the number of rows and columns and any further arguments specify additional matrix dimensions.  The size may also be specified with a vector of dimensions SZ.

     If no size arguments are given then the result matrix is a random ordering of the sample DATA.