rvs_normal - normal distribution random variatesExperimental
A normal continuous random variate distribution, also known as Gaussian, or Gauss or Laplace-Gauss distribution.
The location loc specifies the mean or expectation (). The scale specifies the standard deviation ().
Without argument, the function returns a standard normal distributed random variate .
With two arguments, the function returns a normal distributed random variate . For complex arguments, the real and imaginary parts are independent of each other.
With three arguments, the function returns a rank-1 array of normal distributed random variates.
Note
The algorithm used for generating exponential random variates is fundamentally limited to double precision.1
result = rvs_normal ([loc, scale] [[, array_size]])
Elemental function (passing both loc and scale).
loc: optional argument has intent(in) and is a scalar of type real or complex.
scale: optional argument has intent(in) and is a positive scalar of type real or complex.
array_size: optional argument has intent(in) and is a scalar of type integer.
loc and scale arguments must be of the same type.
The result is a scalar or rank-1 array, with a size of array_size, and the same type as scale and loc. If scale is non-positive, the result is NaN.
program example_normal_rvs
use stdlib_random, only: random_seed
use stdlib_stats_distribution_normal, only: norm => rvs_normal
implicit none
real :: a(2, 3, 4), b(2, 3, 4)
complex :: loc, scale
integer :: seed_put, seed_get
seed_put = 1234567
call random_seed(seed_put, seed_get)
print *, norm() !single standard normal random variate
! 0.563655198
print *, norm(1.0, 2.0)
!normal random variate mu=1.0, sigma=2.0
! -0.633261681
print *, norm(0.0, 1.0, 10) !an array of 10 standard norml random variates
! -3.38123664E-02 -0.190365672 0.220678389 -0.424612164 -0.249541596
! 0.865260184 1.11086845 -0.328349441 1.10873628 1.27049923
a(:, :, :) = 1.0
b(:, :, :) = 1.0
print *, norm(a, b) ! a rank 3 random variates array
!0.152776539 -7.51764774E-02 1.47208166 0.180561781 1.32407105
! 1.20383692 0.123445868 -0.455737948 -0.469808221 1.60750175
! 1.05748117 0.720934749 0.407810807 1.48165631 2.31749439
! 0.414566994 3.06084275 1.86505437 1.36338580 7.26878643E-02
! 0.178585172 1.39557445 0.828021586 0.872084975
loc = (-1.0, 2.0)
scale = (2.0, 1.0)
print *, norm(loc, scale)
!single complex normal random variate with real part of mu=-1, sigma=2;
!imagainary part of mu=2.0 and sigma=1.0
! (1.22566295,2.12518454)
end program example_normal_rvs
pdf_normal - normal distribution probability density functionExperimental
The probability density function (pdf) of the single real variable normal distribution:
For a complex varible with independent real and imaginary parts, the joint probability density function is the product of the the corresponding real and imaginary marginal pdfs:2
result = pdf_normal (x, loc, scale)
Elemental function
x: has intent(in) and is a scalar of type real or complex.
loc: has intent(in) and is a scalar of type real or complex.
scale: has intent(in) and is a positive scalar of type real or complex.
All three arguments must have the same type.
The result is a scalar or an array, with a shape conformable to the arguments, and the same type as the input arguments. If scale is non-positive, the result is NaN.
program example_normal_pdf
use stdlib_random, only: random_seed
use stdlib_stats_distribution_normal, only: norm_pdf => pdf_normal, &
norm => rvs_normal
implicit none
real, dimension(3, 4, 5) :: x, mu, sigma
real :: xsum
complex :: loc, scale
integer :: seed_put, seed_get, i
seed_put = 1234567
call random_seed(seed_put, seed_get)
! probability density at x=1.0 in standard normal
print *, norm_pdf(1.0, 0., 1.)
! 0.241970733
! probability density at x=2.0 with mu=-1.0 and sigma=2.0
print *, norm_pdf(2.0, -1.0, 2.0)
! 6.47588000E-02
! probability density at x=1.0 with mu=1.0 and sigma=-1.0 (out of range)
print *, norm_pdf(1.0, 1.0, -1.0)
! NaN
! standard normal random variates array
x = reshape(norm(0.0, 1.0, 60), [3, 4, 5])
! standard normal probability density array
mu(:, :, :) = 0.0
sigma(:, :, :) = 1.0
print *, norm_pdf(x, mu, sigma)
! 0.340346158 0.285823315 0.398714304 0.391778737 0.389345556
! 0.364551932 0.386712372 0.274370432 0.215250477 0.378006011
! 0.215760440 0.177990928 0.278640658 0.223813817 0.356875211
! 0.285167664 0.378533930 0.390739858 0.271684974 0.138273031
! 0.135456234 0.331718773 0.398283750 0.383706540
! probability density array where sigma<=0.0 for certain elements
print *, norm_pdf([1.0, 1.0, 1.0], [1.0, 1.0, 1.0], [1.0, 0.0, -1.0])
! 0.398942292 NaN NaN
! `pdf_normal` is pure and, thus, can be called concurrently
xsum = 0.0
do concurrent (i=1:size(x,3))
xsum = xsum + sum(norm_pdf(x(:,:,i), mu(:,:,i), sigma(:,:,i)))
end do
print *, xsum
! 18.0754433
! complex normal probability density function at x=(1.5,1.0) with real part
! of mu=1.0, sigma=1.0 and imaginary part of mu=-0.5, sigma=2.0
loc = (1.0, -0.5)
scale = (1.0, 2.)
print *, norm_pdf((1.5, 1.0), loc, scale)
! 5.30100204E-02
end program example_normal_pdf
cdf_normal - normal distribution cumulative distribution functionExperimental
Cumulative distribution function of the single real variable normal distribution:
For the complex variable with independent real and imaginary parts, the joint cumulative distribution function is the product of the corresponding real and imaginary marginal cdfs:2
result = cdf_normal (x, loc, scale)
Elemental function
x: has intent(in) and is a scalar of type real or complex.
loc: has intent(in) and is a scalar of type real or complex.
scale: has intent(in) and is a positive scalar of type real or complex.
All three arguments must have the same type.
The result is a scalar or an array, with a shape conformable to the arguments, and the same type as the input arguments. If scale is non-positive, the result is NaN.
program example_normal_cdf
use stdlib_random, only: random_seed
use stdlib_stats_distribution_normal, only: norm_cdf => cdf_normal, &
norm => rvs_normal
implicit none
real, dimension(2, 3, 4) :: x, mu, sigma
real :: xsum
complex :: loc, scale
integer :: seed_put, seed_get, i
seed_put = 1234567
call random_seed(seed_put, seed_get)
! standard normal cumulative probability at x=0.0
print *, norm_cdf(0.0, 0.0, 1.0)
! 0.500000000
! cumulative probability at x=2.0 with mu=-1.0 sigma=2.0
print *, norm_cdf(2.0, -1.0, 2.0)
! 0.933192849
! cumulative probability at x=1.0 with mu=1.0 and sigma=-1.0 (out of range)
print *, norm_cdf(1.0, 1.0, -1.0)
! NaN
! standard normal random variates array
x = reshape(norm(0.0, 1.0, 24), [2, 3, 4])
! standard normal cumulative array
mu(:, :, :) = 0.0
sigma(:, :, :) = 1.0
print *, norm_cdf(x, mu, sigma)
! 0.713505626 0.207069695 0.486513376 0.424511284 0.587328553
! 0.335559726 0.401470929 0.806552052 0.866687536 0.371323735
! 0.866228044 0.898046613 0.198435277 0.141147852 0.681565762
! 0.206268221 0.627057910 0.580759525 0.190364420 7.27325380E-02
! 7.08068311E-02 0.728241026 0.522919059 0.390097380
! cumulative probability array where sigma<=0.0 for certain elements
print *, norm_cdf([1.0, 1.0, 1.0], [1.0, 1.0, 1.0], [1.0, 0.0, -1.0])
! 0.500000000 NaN NaN
! `cdf_normal` is pure and, thus, can be called concurrently
xsum = 0.0
do concurrent (i=1:size(x,3))
xsum = xsum + sum(norm_cdf(x(:,:,i), mu(:,:,i), sigma(:,:,i)))
end do
print *, xsum
! 11.3751936
! complex normal cumulative distribution at x=(0.5,-0.5) with real part of
! mu=1.0, sigma=0.5 and imaginary part of mu=0.0, sigma=1.0
loc = (1.0, 0.0)
scale = (0.5, 1.0)
print *, norm_cdf((0.5, -0.5), loc, scale)
! 4.89511043E-02
end program example_normal_cdf
Marsaglia, George, and Wai Wan Tsang. "The ziggurat method for generating random variables." Journal of statistical software 5 (2000): 1-7. ↩
Miller, Scott, and Donald Childers. Probability and random processes: With applications to signal processing and communications. Academic Press, 2012 (p. 197). ↩↩