Skip to content
24 changes: 16 additions & 8 deletions src/stdlib_experimental_stats.fypp
Original file line number Diff line number Diff line change
Expand Up @@ -213,9 +213,10 @@ module stdlib_experimental_stats
#:for k1, t1 in RC_KINDS_TYPES
#:for rank in RANKS
#:set RName = rname("moment_all",rank, t1, k1)
module function ${RName}$(x, order, mask) result(res)
module function ${RName}$(x, order, center, mask) result(res)
${t1}$, intent(in) :: x${ranksuffix(rank)}$
integer, intent(in) :: order
${t1}$, intent(in), optional :: center
logical, intent(in), optional :: mask
${t1}$ :: res
end function ${RName}$
Expand All @@ -225,9 +226,10 @@ module stdlib_experimental_stats
#:for k1, t1 in INT_KINDS_TYPES
#:for rank in RANKS
#:set RName = rname("moment_all",rank, t1, k1, 'dp')
module function ${RName}$(x, order, mask) result(res)
module function ${RName}$(x, order, center, mask) result(res)
${t1}$, intent(in) :: x${ranksuffix(rank)}$
integer, intent(in) :: order
real(dp), intent(in), optional :: center
logical, intent(in), optional :: mask
real(dp) :: res
end function ${RName}$
Expand All @@ -237,10 +239,11 @@ module stdlib_experimental_stats
#:for k1, t1 in RC_KINDS_TYPES
#:for rank in RANKS
#:set RName = rname("moment",rank, t1, k1)
module function ${RName}$(x, order, dim, mask) result(res)
module function ${RName}$(x, order, dim, center, mask) result(res)
${t1}$, intent(in) :: x${ranksuffix(rank)}$
integer, intent(in) :: order
integer, intent(in) :: dim
${t1}$, intent(in), optional :: center${reduced_shape('x', rank, 'dim')}$
logical, intent(in), optional :: mask
${t1}$ :: res${reduced_shape('x', rank, 'dim')}$
end function ${RName}$
Expand All @@ -250,10 +253,11 @@ module stdlib_experimental_stats
#:for k1, t1 in INT_KINDS_TYPES
#:for rank in RANKS
#:set RName = rname("moment",rank, t1, k1, 'dp')
module function ${RName}$(x, order, dim, mask) result(res)
module function ${RName}$(x, order, dim, center, mask) result(res)
${t1}$, intent(in) :: x${ranksuffix(rank)}$
integer, intent(in) :: order
integer, intent(in) :: dim
real(dp),intent(in), optional :: center${reduced_shape('x', rank, 'dim')}$
logical, intent(in), optional :: mask
real(dp) :: res${reduced_shape('x', rank, 'dim')}$
end function ${RName}$
Expand All @@ -263,9 +267,10 @@ module stdlib_experimental_stats
#:for k1, t1 in RC_KINDS_TYPES
#:for rank in RANKS
#:set RName = rname("moment_mask_all",rank, t1, k1)
module function ${RName}$(x, order, mask) result(res)
module function ${RName}$(x, order, center, mask) result(res)
${t1}$, intent(in) :: x${ranksuffix(rank)}$
integer, intent(in) :: order
${t1}$, intent(in), optional :: center
logical, intent(in) :: mask${ranksuffix(rank)}$
${t1}$ :: res
end function ${RName}$
Expand All @@ -275,9 +280,10 @@ module stdlib_experimental_stats
#:for k1, t1 in INT_KINDS_TYPES
#:for rank in RANKS
#:set RName = rname("moment_mask_all",rank, t1, k1, 'dp')
module function ${RName}$(x, order, mask) result(res)
module function ${RName}$(x, order, center, mask) result(res)
${t1}$, intent(in) :: x${ranksuffix(rank)}$
integer, intent(in) :: order
real(dp),intent(in), optional :: center
logical, intent(in) :: mask${ranksuffix(rank)}$
real(dp) :: res
end function ${RName}$
Expand All @@ -287,10 +293,11 @@ module stdlib_experimental_stats
#:for k1, t1 in RC_KINDS_TYPES
#:for rank in RANKS
#:set RName = rname("moment_mask",rank, t1, k1)
module function ${RName}$(x, order, dim, mask) result(res)
module function ${RName}$(x, order, dim, center, mask) result(res)
${t1}$, intent(in) :: x${ranksuffix(rank)}$
integer, intent(in) :: order
integer, intent(in) :: dim
${t1}$, intent(in), optional :: center${reduced_shape('x', rank, 'dim')}$
logical, intent(in) :: mask${ranksuffix(rank)}$
${t1}$ :: res${reduced_shape('x', rank, 'dim')}$
end function ${RName}$
Expand All @@ -300,10 +307,11 @@ module stdlib_experimental_stats
#:for k1, t1 in INT_KINDS_TYPES
#:for rank in RANKS
#:set RName = rname("moment_mask",rank, t1, k1, 'dp')
module function ${RName}$(x, order, dim, mask) result(res)
module function ${RName}$(x, order, dim, center, mask) result(res)
${t1}$, intent(in) :: x${ranksuffix(rank)}$
integer, intent(in) :: order
integer, intent(in) :: dim
real(dp), intent(in), optional :: center${reduced_shape('x', rank, 'dim')}$
logical, intent(in) :: mask${ranksuffix(rank)}$
real(dp) :: res${reduced_shape('x', rank, 'dim')}$
end function ${RName}$
Expand Down
45 changes: 30 additions & 15 deletions src/stdlib_experimental_stats.md
Original file line number Diff line number Diff line change
@@ -1,7 +1,7 @@
# Descriptive statistics

* [`mean` - mean of array elements](#mean---mean-of-array-elements)
* [`moment` - central moment of array elements](#moment---central-moment-of-array-elements)
* [`moment` - raw and central moments of array elements](#moment---raw-and-central-moments-of-array-elements)
* [`var` - variance of array elements](#var---variance-of-array-elements)

## `mean` - mean of array elements
Expand Down Expand Up @@ -48,25 +48,34 @@ program demo_mean
end program demo_mean
```

## `moment` - central moment of array elements
## `moment` - raw and central moments of array elements

### Description

Returns the _k_-th order central moment of all the elements of `array`, or of the elements of `array` along dimension `dim` if provided, and if the corresponding element in `mask` is `true`.
Returns the _k_-th order raw moment of all the elements of `array`, or of the elements of `array` along dimension `dim` if provided, and if the corresponding element in `mask` is `true`.

The _k_-th order central moment is defined as :
If an array `center` is provided, the function returns the _k_-th order central moment of all the elements of `array`, or of the elements of `array` along dimension `dim` if provided, and if the corresponding element in `mask` is `true`.


The _k_-th order raw moment is defined as :

```
moment(array) = 1/n sum_i (array(i) - mean(array))^k
moment(array) = 1/n sum_i (array(i))^k
```

where n is the number of elements.

The _k_-th order central moment is defined as :

```
moment(array) = 1/n sum_i (array(i) - center)^k
```

### Syntax

`result = moment(array, order [, mask])`
`result = moment(array, order [, center [, mask]])`

`result = moment(array, order, dim [, mask])`
`result = moment(array, order, dim [, center [, mask]])`

### Arguments

Expand All @@ -76,29 +85,35 @@ where n is the number of elements.

`dim`: Shall be a scalar of type `integer` with a value in the range from 1 to n, where n is the rank of `array`.

`center` (optional): Shall be a scalar of the same type of `array` if `array` is `real` or `complex`, or of type `real(dp)` if `array` is of type `integer`. If `dim` is provided, `center` shall be an array (with a shape similar to that of `array` with dimension `dim` dropped) of the same type of `array` if `array` is `real` or `complex`, or of type `real(dp)` if `array` is of type `integer`.

`mask` (optional): Shall be of type `logical` and either by a scalar or an array of the same shape as `array`.

### Return value

If `array` is of type `real` or `complex`, the result is of the same type as `array`.
If `array` is of type `integer`, the result is of type `real(dp)`.

If `dim` is absent, a scalar with the _k_-th central moment of all elements in `array` is returned. Otherwise, an array of rank n-1, where n equals the rank of `array`, and a shape similar to that of `array` with dimension `dim` dropped is returned.
If `dim` is absent, a scalar with the _k_-th raw (or central if `center` is provided) moment of all elements in `array` is returned. Otherwise, an array of rank n-1, where n equals the rank of `array`, and a shape similar to that of `array` with dimension `dim` dropped is returned.

If `mask` is specified, the result is the _k_-th central moment of all elements of `array` corresponding to `true` elements of `mask`. If every element of `mask` is `false`, the result is IEEE `NaN`.
If `mask` is specified, the result is the _k_-th raw (or central if `center` is provided) moment of all elements of `array` corresponding to `true` elements of `mask`. If every element of `mask` is `false`, the result is IEEE `NaN`.

### Example

```fortran
program demo_moment
use stdlib_experimental_stats, only: moment
use stdlib_experimental_stats, only: mean, moment
implicit none
real :: x(1:6) = [ 1., 2., 3., 4., 5., 6. ]
print *, moment(x, 2) !returns 2.9167
print *, moment( reshape(x, [ 2, 3 ] ), 2) !returns 2.9167
print *, moment( reshape(x, [ 2, 3 ] ), 2, 1) !returns [0.25, 0.25, 0.25]
print *, moment( reshape(x, [ 2, 3 ] ), 2, 1,&
reshape(x, [ 2, 3 ] ) > 3.) !returns [NaN, 0., 0.25]
real :: y(1:2, 1:3) = reshape([ 1., 2., 3., 4., 5., 6. ], [ 2, 3])
print *, moment(x, 2, center = mean(x)) !returns 2.9167
print *, moment( y, 2,&
center = mean(y)) !returns 2.9167
print *, moment( y, 2, 1,&
center = mean(y, 1)) !returns [0.25, 0.25, 0.25]
print *, moment( y, 2, 1,&
center = mean(y, 1, y > 3.),&
mask = (y > 3.)) !returns [NaN, 0., 0.25]
end program demo_moment
```

Expand Down
Loading