Skip to content

Probability Distribution and Statistical Functions -- Uniform Distribution Module #593

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Merged
merged 3 commits into from
Dec 12, 2021
Merged

Probability Distribution and Statistical Functions -- Uniform Distribution Module #593

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
3 commits
Select commit Hold shift + click to select a range
35a3ecc
add pdf/cdf with corresponding precisions
Dec 11, 2021
fe14e4c
correction for kinds
Dec 11, 2021
e0df0e2
add description of variable loc for pdf/cdf
Dec 12, 2021
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
36 changes: 24 additions & 12 deletions doc/specs/stdlib_stats_distribution_uniform.md
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -36,7 +36,7 @@ Return a randomized rank one array of the input type.

```fortran
program demo_shuffle
use stdlib_stats_distribution_PRNG, only : random_seed
use stdlib_random, only : random_seed
use stdlib_stats_distribution_uniform, only : shuffle
implicit none
integer :: seed_put, seed_get, i
Expand Down Expand Up @@ -79,14 +79,16 @@ Without argument the function returns a scalar standard uniformly distributed va

With single argument `scale` of `integer` type the function returns a scalar uniformly distributed variate of `integer` type on [0,scale]. This is the standard Rectangular distribution.

With single argument `scale` of `real` or `complex` type the function returns a scalar uniformly distributed variate of `real` or `complex` type on [0, scale].
With single argument `scale` of `real` or `complex` type the function returns a scalar uniformly distributed variate of `real` type on [0, scale] or `complex` type on [(0, 0i), (scale, i(scale))].

With double arguments `loc` and `scale` the function returns a scalar uniformly distributed random variates of `integer`, `real` or `complex` type on [loc, loc + scale] dependent of input type.
With double arguments `loc` and `scale` the function returns a scalar uniformly distributed random variates of `integer` or `real` type on [loc, loc + scale], or `complex` type on [(loc, i(loc)), ((loc + scale), i(loc + scale))], dependent of input type.

With triple arguments `loc`, `scale` and `array_size` the function returns a rank one array of uniformly distributed variates of `integer`, `real` or `complex` type with an array size of `array_size`.

For `complex` type, the real part and imaginary part are independent of each other.

Note: the algorithm used for generating uniform random variates is fundamentally limited to double precision.

### Syntax

`result = [[stdlib_stats_distribution_uniform(module):rvs_uniform(interface)]]([[loc,] scale] [[[,array_size]]])`
Expand All @@ -101,19 +103,19 @@ Elemental function (without the third argument).

`scale`: optional argument has `intent(in)` and is a scalar of type `integer`, `real` or `complex`.

`array_size`: optional argument has `intent(in)` and is a scalar of type `integer`.
`array_size`: optional argument has `intent(in)` and is a scalar of type `integer` with default kind.

`loc` and `scale` must have the same type and kind when both are present.

### Return value

The result is a scalar or a rank one array, with size of `array_size`, of type `integer`, `real` or `complex` depending on the input type.
The result is a scalar or a rank one array with size of `array_size`, of type `integer`, `real` or `complex` depending on the input type.

### Example

```fortran
program demo_uniform_rvs
use stdlib_stats_distribution_PRNG, only:random_seed
use stdlib_random, only:random_seed
use stdlib_stats_distribution_uniform, only:uni=> rvs_uniform

implicit none
Expand Down Expand Up @@ -193,15 +195,19 @@ program demo_uniform_rvs
end program demo_uniform_rvs
```

## `pdf_uniform` - Uniform probability density function
## `pdf_uniform` - Uniform distribution probability density function

### Status

Experimental

### Description

The probability density function of the uniform distribution.
The probability density function of the uniform distribution:

f(x) = 0 x < loc or x > loc + scale for all types uniform distributions

For random variable x in [loc, loc + scale]:

f(x) = 1 / (scale + 1); for discrete uniform distribution.

Expand Down Expand Up @@ -235,7 +241,7 @@ The result is a scalar or an array, with a shape conformable to arguments, of ty

```fortran
program demo_uniform_pdf
use stdlib_stats_distribution_PRNG, only : random_seed
use stdlib_random, only : random_seed
use stdlib_stats_distribution_uniform, only : uni_pdf => pdf_uniform, &
uni => rvs_uniform

Expand Down Expand Up @@ -286,15 +292,21 @@ end program demo_uniform_pdf

```

## `cdf_uniform` - Uniform cumulative distribution function
## `cdf_uniform` - Uniform distribution cumulative distribution function

### Status

Experimental

### Description

Cumulative distribution function of the uniform distribution
Cumulative distribution function of the uniform distribution:

F(x) = 0 x < loc for all types uniform distributions

F(x) = 1 x > loc + scale for all types uniform distributions

For random variable x in [loc, loc + scale]:

F(x) = (x - loc + 1) / (scale + 1); for discrete uniform distribution.

Expand Down Expand Up @@ -328,7 +340,7 @@ The result is a scalar or an array, with a shape conformable to arguments, of ty

```fortran
program demo_uniform_cdf
use stdlib_stats_distribution_PRNG, only : random_seed
use stdlib_random, only : random_seed
use stdlib_stats_distribution_uniform, only : uni_cdf => cdf_uniform, &
uni => rvs_uniform

Expand Down
51 changes: 27 additions & 24 deletions src/stdlib_stats_distribution_uniform.fypp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -2,7 +2,7 @@
#:set RC_KINDS_TYPES = REAL_KINDS_TYPES + CMPLX_KINDS_TYPES
#:set ALL_KINDS_TYPES = INT_KINDS_TYPES + RC_KINDS_TYPES
module stdlib_stats_distribution_uniform
use stdlib_kinds
use stdlib_kinds, only : sp, dp, xdp, qp, int8, int16, int32, int64
use stdlib_error, only : error_stop
use stdlib_random, only : dist_rand

Expand Down Expand Up @@ -383,14 +383,15 @@ contains
elemental function pdf_unif_${t1[0]}$${k1}$(x, loc, scale) result(res)

${t1}$, intent(in) :: x, loc, scale
real :: res
${t1}$ :: res
${t1}$, parameter :: zero = 0.0_${k1}$, one = 1.0_${k1}$

if(scale == 0.0_${k1}$) then
res = 0.0
if(scale == zero) then
res = zero
else if(x < loc .or. x > (loc + scale)) then
res = 0.0
res = zero
else
res = 1.0 / scale
res = one / scale
end if
end function pdf_unif_${t1[0]}$${k1}$

Expand All @@ -402,17 +403,17 @@ contains
elemental function pdf_unif_${t1[0]}$${k1}$(x, loc, scale) result(res)

${t1}$, intent(in) :: x, loc, scale
real :: res
real(${k1}$) :: tr, ti
real(${k1}$) :: res, tr, ti
real(${k1}$), parameter :: zero = 0.0_${k1}$, one = 1.0_${k1}$

tr = loc % re + scale % re; ti = loc % im + scale % im
if(scale == (0.0_${k1}$,0.0_${k1}$)) then
res = 0.0
if(scale == (zero, zero)) then
res = zero
else if((x % re >= loc % re .and. x % re <= tr) .and. &
(x % im >= loc % im .and. x % im <= ti)) then
res = 1.0 / (scale % re * scale % im)
res = one / (scale % re * scale % im)
else
res = 0.0
res = zero
end if
end function pdf_unif_${t1[0]}$${k1}$

Expand Down Expand Up @@ -445,16 +446,17 @@ contains
elemental function cdf_unif_${t1[0]}$${k1}$(x, loc, scale) result(res)

${t1}$, intent(in) :: x, loc, scale
real :: res
${t1}$ :: res
${t1}$, parameter :: zero = 0.0_${k1}$, one = 1.0_${k1}$

if(scale == 0.0_${k1}$) then
res = 0.0
if(scale == zero) then
res = zero
else if(x < loc) then
res = 0.0
res = zero
else if(x >= loc .and. x <= (loc + scale)) then
res = (x - loc) / scale
else
res = 1.0
res = one
end if
end function cdf_unif_${t1[0]}$${k1}$

Expand All @@ -466,29 +468,30 @@ contains
elemental function cdf_unif_${t1[0]}$${k1}$(x, loc, scale) result(res)

${t1}$, intent(in) :: x, loc, scale
real :: res
real(${k1}$) :: res
real(${k1}$), parameter :: zero = 0.0_${k1}$, one = 1.0_${k1}$
logical :: r1, r2, i1, i2

if(scale == (0.0_${k1}$,0.0_${k1}$)) then
res = 0.0
if(scale == (zero, zero)) then
res = zero
return
end if
r1 = x % re < loc % re
r2 = x % re > (loc % re + scale % re)
i1 = x % im < loc % im
i2 = x % im > (loc % im + scale % im)
if(r1 .or. i1) then
res = 0.0
res = zero
else if((.not. r1) .and. (.not. r2) .and. i2) then
res = (x % re - loc % re) / scale % re
else if((.not. i1) .and. (.not. i2) .and. r2) then
res = (x % im - loc % im) / scale % im
else if((.not. r1) .and. (.not. r2) .and. (.not. i1) .and. (.not. i2)) &
then
res = (x % re - loc % re) * (x % im - loc % im) / &
(scale % re * scale % im)
res = ((x % re - loc % re) / scale % re) * ((x % im - loc % im) / &
scale % im)
else if(r2 .and. i2)then
res = 1.0
res = one
end if
end function cdf_unif_${t1[0]}$${k1}$

Expand Down
72 changes: 57 additions & 15 deletions src/tests/stats/test_distribution_uniform.fypp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -207,16 +207,18 @@ contains
subroutine test_uni_pdf_${t1[0]}$${k1}$
${t1}$ :: x1, x2(3,4), loc, scale
integer :: seed, get, i
real :: res(3,5)
#:if t1[0] == "i"
#! for integer type
real :: res(3, 5)
real, parameter :: ans(15) = [(1.96078438E-02, i=1,15)]
#:elif t1[0] == "r"
#! for real type
real, parameter :: ans(15) = [(0.5, i=1,15)]
${t1}$ :: res(3, 5)
${t1}$, parameter :: ans(15) = [(0.5_${k1}$, i=1,15)]
#:else
#! for complex type
real, parameter :: ans(15) = [(1.0, i=1,15)]
real(${k1}$) :: res(3, 5)
real(${k1}$), parameter :: ans(15) = [(1.0_${k1}$, i=1,15)]
#:endif

print *, "Test pdf_uniform_${t1[0]}$${k1}$"
Expand All @@ -236,8 +238,15 @@ contains
x2 = reshape(uni_rvs(loc, scale, 12), [3,4])
res(:,1) = uni_pdf(x1, loc, scale)
res(:, 2:5) = uni_pdf(x2, loc, scale)
#:if t1[0] == "i"
#! for integer type
call check(all(abs(res - reshape(ans,[3,5])) < sptol), &
msg = "pdf_uniform_${t1[0]}$${k1}$ failed", warn=warn)
#:else
#! for real and complex types
call check(all(abs(res - reshape(ans,[3,5])) < ${k1}$tol), &
msg = "pdf_uniform_${t1[0]}$${k1}$ failed", warn=warn)
#:endif
end subroutine test_uni_pdf_${t1[0]}$${k1}$

#:endfor
Expand All @@ -247,47 +256,73 @@ contains
#:for k1, t1 in ALL_KINDS_TYPES
subroutine test_uni_cdf_${t1[0]}$${k1}$
${t1}$ :: x1, x2(3,4), loc, scale
real :: res(3,5)
integer :: seed, get
#:if k1 == "int8"
real :: res(3,5)
real, parameter :: ans(15) = [0.435643554, 0.435643554, 0.435643554, &
0.702970326, 0.653465331, 0.485148519, &
0.386138618, 0.386138618, 0.336633652, &
0.277227730, 0.237623766, 0.524752498, &
0.732673287, 0.534653485, 0.415841579]
#:elif k1 == "int16"
real :: res(3,5)
real, parameter :: ans(15) = [0.178217828, 0.178217828, 0.178217828, &
0.465346545, 0.673267305, 0.247524753, &
0.158415839, 0.792079210, 0.742574275, &
0.574257433, 0.881188095, 0.663366318, &
0.524752498, 0.623762369, 0.178217828]
#:elif k1 == "int32"
real :: res(3,5)
real, parameter :: ans(15) = [0.732673287, 0.732673287, 0.732673287, &
0.722772300, 0.792079210, 5.94059415E-02,&
0.841584146, 0.405940592, 0.960396051, &
0.534653485, 0.782178223, 0.861386120, &
0.564356446, 0.613861382, 0.306930691]
#:elif k1 == "int64"
real :: res(3,5)
real, parameter :: ans(15) = [0.455445558, 0.455445558, 0.455445558, &
0.277227730, 0.455445558, 0.930693090, &
0.851485133, 0.623762369, 5.94059415E-02,&
0.693069279, 0.544554472, 0.207920790, &
0.306930691, 0.356435657, 0.128712878]
#:elif t1[0] == "r"
#! for real type
real, parameter :: ans(15) = [0.170192942, 0.170192942, 0.170192942, &
0.276106149, 0.754930079, 0.406620681, &
0.187742814, 0.651605546, 0.934733927, &
0.151271492, 0.987674534, 0.130533904, &
0.106271908, 9.27578658E-02, 0.203399420]
${t1}$ :: res(3,5)
${t1}$, parameter :: ans(15) = &
[0.170192944297557408050991512027394492_${k1}$, &
0.170192944297557408050991512027394492_${k1}$, &
0.170192944297557408050991512027394492_${k1}$, &
0.276106146274646191418611351764411665_${k1}$, &
0.754930097473875072466853453079238534_${k1}$, &
0.406620682573118008562573777453508228_${k1}$, &
0.187742819191801080247472555129206739_${k1}$, &
0.651605526090507591874256831943057477_${k1}$, &
0.934733949732104885121941606485052034_${k1}$, &
0.151271491851613287815681019310432021_${k1}$, &
0.987674522797719611766353864368284121_${k1}$, &
0.130533899463404684526679488953959662_${k1}$, &
0.106271905921876880229959283497009892_${k1}$, &
9.27578652240113182836367400341259781E-0002_${k1}$, &
0.203399426853420439709196898547816090_${k1}$]
#:else
#! for complex type
real, parameter :: ans(15) = [4.69913185E-02, 4.69913185E-02, &
4.69913185E-02, 0.306970179, 0.122334257,&
0.141398594, 0.128925011, 9.85755492E-03,&
8.16527531E-02, 0.163921610, &
7.81712309E-02, 0.446415812, &
5.31753292E-04, 0.101455867, 0.155276477]
real(${k1}$) :: res(3,5)
real(${k1}$), parameter :: ans(15) = &
[4.69913179731340971083526490627998168E-0002_${k1}$, &
4.69913179731340971083526490627998168E-0002_${k1}$, &
4.69913179731340971083526490627998168E-0002_${k1}$, &
0.306970191529817593217448363707416739_${k1}$, &
0.122334258469188588238756489506609443_${k1}$, &
0.141398599060326408705075175176932616_${k1}$, &
0.128925006861443729884744412460848140_${k1}$, &
9.85755512660026594506599410104817938E-0003_${k1}$, &
8.16527497645585445208592050401597260E-0002_${k1}$, &
0.163921605454974749736935624944263178_${k1}$, &
7.81712317416218284294000447064256003E-0002_${k1}$, &
0.446415807686727375005224206895756087_${k1}$, &
5.31753272901435018841591264266743165E-0004_${k1}$, &
0.101455865191097416942685556683943046_${k1}$, &
0.155276470981788516449112374966730510_${k1}$]
#:endif

print *, "Test cdf_uniform_${t1[0]}$${k1}$"
Expand All @@ -307,8 +342,15 @@ contains
x2 = reshape(uni_rvs(loc, scale, 12), [3,4])
res(:,1) = uni_cdf(x1, loc, scale)
res(:, 2:5) = uni_cdf(x2, loc, scale)
#:if t1[0] == "i"
#! for integer type
call check(all(abs(res - reshape(ans,[3,5])) < sptol), &
msg = "cdf_uniform_${t1[0]}$${k1}$ failed", warn=warn)
#:else
#! for real and complex types
call check(all(abs(res - reshape(ans,[3,5])) < ${k1}$tol), &
msg = "cdf_uniform_${t1[0]}$${k1}$ failed", warn=warn)
#:endif
end subroutine test_uni_cdf_${t1[0]}$${k1}$

#:endfor
Expand Down