ParaMonte Fortran 2.0.0
Parallel Monte Carlo and Machine Learning Library
See the latest version documentation.
pm_distCov::getCovRand Interface Reference

Generate and return a random positive-definite (correlation or covariance) matrix using the Gram method.
More...

Detailed Description

Generate and return a random positive-definite (correlation or covariance) matrix using the Gram method.

See the documentation of pm_distCov for details.
See also setCovRand for generating random covariance matrices using the method of Dvine or Onion.

Parameters
[in]mold: The input scalar of,
  1. type complex of kind any supported by the processor (e.g., CK, CK32, CK64, or CK128),
  2. type real of kind any supported by the processor (e.g., RK, RK32, RK64, or RK128),
whose type and kind determines the type and kind of the output rand.
The value of mold is ignored entirely within the algorithm.
[in]ndim: The input positive scalar of type integer of default kind IK, representing the rank of the matrix (the number of dimensions) of shape (ndim, ndim).
(optional. It must be present if and only if the input scale argument is missing or is a scalar.)
[in]scale: The input scalar or contiguous vector of size ndim of type real of the same kind as the output argument rand, representing the scale of the matrix (e.g., the standard deviation of a covariance matrix) along each dimension.
(optional. default = 1.. It can be present if and only if it is a scalar or, it is a vector and the input argument ndim is missing.)
Returns
rand : The output matrix of shape (1:ndim, 1:ndim) of the same type and kind as the input argument mold, containing a random positive-definite matrix.


Possible calling interfaces

rand(1:ndim, 1:ndim) = getCovRand(mold, scale(1:ndim))
rand(1:ndim, 1:ndim) = getCovRand(mold, ndim, scale = scale)
Generate and return a random positive-definite (correlation or covariance) matrix using the Gram meth...
Definition: pm_distCov.F90:394
This module contains classes and procedures for generating random matrices distributed on the space o...
Definition: pm_distCov.F90:72
Warning
The condition all([0 < scale]) must hold for the corresponding input arguments.
These conditions are verified only if the library is built with the preprocessor macro CHECK_ENABLED=1.
Note
Unlike the case for setCovRand, when the input argument scale is missing, the diagonal elements of the output correlation matrix are strictly enforced to 1.


Example usage

1program example
2
3 use pm_kind, only: SK, IK, LK
4 use pm_io, only: display_type
5 use pm_distUnif, only: getUnifRand
6 use pm_matrixChol, only: setChoLow
9 use pm_distCov, only: getCovRand
10
11 implicit none
12
13 integer(IK) :: itry, ndim
14
15 type(display_type) :: disp
16 disp = display_type(file = "main.out.F90")
17
18 block
19 use pm_kind, only: TKG => RKS ! all real kinds are supported.
20 real(TKG), allocatable :: rand(:,:)
21 real(TKG) :: mold
22 real(TKG) :: scale
23 do itry = 1, 5
24 call disp%skip()
25 call disp%show("ndim = getUnifRand(2, 3)")
26 ndim = getUnifRand(2, 3)
27 call disp%show("scale = getUnifRand(1, 10)")
28 scale = getUnifRand(1, 10)
29 call disp%show("rand = getCovRand(mold, ndim)")
30 rand = getCovRand(mold, ndim)
31 call disp%show("rand")
32 call disp%show( rand )
33 call disp%show("isMatClass(rand, posdefmat)")
34 call disp%show( isMatClass(rand, posdefmat) )
35 call disp%show("rand = getCovRand(mold, ndim, scale)")
36 rand = getCovRand(mold, ndim, scale)
37 call disp%show("rand")
38 call disp%show( rand )
39 call disp%show("isMatClass(rand, posdefmat)")
40 call disp%show( isMatClass(rand, posdefmat) )
41 call disp%skip()
42 call disp%skip()
43 call disp%show("rand = getCovRand(mold = 0._TKG, scale = [(real(itry, TKG), itry = 1, 5)])")
44 rand = getCovRand(mold = 0._TKG, scale = [(real(itry, TKG), itry = 1, 5)])
45 call disp%show("rand")
46 call disp%show( rand )
47 call disp%show("isMatClass(rand, posdefmat)")
48 call disp%show( isMatClass(rand, posdefmat) )
49 call disp%skip()
50 end do
51 end block
52
53 block
54 use pm_kind, only: TKG => CKS ! all real kinds are supported.
55 complex(TKG), allocatable :: rand(:,:)
56 complex(TKG) :: mold
57 real(TKG) :: scale
58 do itry = 1, 5
59 call disp%skip()
60 call disp%show("ndim = getUnifRand(2, 3)")
61 ndim = getUnifRand(2, 3)
62 call disp%show("scale = getUnifRand(1, 10)")
63 scale = getUnifRand(1, 10)
64 call disp%show("rand = getCovRand(mold, ndim)")
65 rand = getCovRand(mold, ndim)
66 call disp%show("rand")
67 call disp%show( rand )
68 call disp%show("isMatClass(rand, posdefmat)")
69 call disp%show( isMatClass(rand, posdefmat) )
70 call disp%show("rand = getCovRand(mold, ndim, scale)")
71 rand = getCovRand(mold, ndim, scale)
72 call disp%show("rand")
73 call disp%show( rand )
74 call disp%show("isMatClass(rand, posdefmat)")
75 call disp%show( isMatClass(rand, posdefmat) )
76 call disp%skip()
77 call disp%skip()
78 call disp%show("rand = getCovRand(mold = 0._TKG, scale = [(real(itry, TKG), itry = 1, 5)])")
79 rand = getCovRand(mold = 0._TKG, scale = [(real(itry, TKG), itry = 1, 5)])
80 call disp%show("rand")
81 call disp%show( rand )
82 call disp%show("isMatClass(rand, posdefmat)")
83 call disp%show( isMatClass(rand, posdefmat) )
84 call disp%skip()
85 end do
86 end block
87
88end program example
Generate and return a (collection) of random vector(s) of size ndim from the ndim-dimensional MultiVa...
Generate and return a scalar or a contiguous array of rank 1 of length s1 of randomly uniformly distr...
This is a generic method of the derived type display_type with pass attribute.
Definition: pm_io.F90:11726
This is a generic method of the derived type display_type with pass attribute.
Definition: pm_io.F90:11508
[LEGACY code] Return the lower-triangle of the Cholesky factorization of the symmetric positive-def...
Generate and return .true. if and only if the input matrix is of the specified input class.
This module contains classes and procedures for computing various statistical quantities related to t...
This module contains classes and procedures for computing various statistical quantities related to t...
This module contains classes and procedures for input/output (IO) or generic display operations on st...
Definition: pm_io.F90:252
type(display_type) disp
This is a scalar module variable an object of type display_type for general display.
Definition: pm_io.F90:11393
This module defines the relevant Fortran kind type-parameters frequently used in the ParaMonte librar...
Definition: pm_kind.F90:268
integer, parameter LK
The default logical kind in the ParaMonte library: kind(.true.) in Fortran, kind(....
Definition: pm_kind.F90:541
integer, parameter CKS
The single-precision complex kind in Fortran mode. On most platforms, this is a 32-bit real kind.
Definition: pm_kind.F90:570
integer, parameter IK
The default integer kind in the ParaMonte library: int32 in Fortran, c_int32_t in C-Fortran Interoper...
Definition: pm_kind.F90:540
integer, parameter SK
The default character kind in the ParaMonte library: kind("a") in Fortran, c_char in C-Fortran Intero...
Definition: pm_kind.F90:539
integer, parameter RKS
The single-precision real kind in Fortran mode. On most platforms, this is an 32-bit real kind.
Definition: pm_kind.F90:567
This module contains procedures and generic interfaces for computing the Cholesky factorization of po...
This module contains abstract and concrete derived types that are required for compile-time resolutio...
type(posdefmat_type), parameter posdefmat
This is a scalar parameter object of type hermitian_type that is exclusively used to signify the Herm...
Generate and return an object of type display_type.
Definition: pm_io.F90:10282

Example Unix compile command via Intel ifort compiler
1#!/usr/bin/env sh
2rm main.exe
3ifort -fpp -standard-semantics -O3 -Wl,-rpath,../../../lib -I../../../inc main.F90 ../../../lib/libparamonte* -o main.exe
4./main.exe

Example Windows Batch compile command via Intel ifort compiler
1del main.exe
2set PATH=..\..\..\lib;%PATH%
3ifort /fpp /standard-semantics /O3 /I:..\..\..\include main.F90 ..\..\..\lib\libparamonte*.lib /exe:main.exe
4main.exe

Example Unix / MinGW compile command via GNU gfortran compiler
1#!/usr/bin/env sh
2rm main.exe
3gfortran -cpp -ffree-line-length-none -O3 -Wl,-rpath,../../../lib -I../../../inc main.F90 ../../../lib/libparamonte* -o main.exe
4./main.exe

Example output
1
2ndim = getUnifRand(2, 3)
3scale = getUnifRand(1, 10)
4rand = getCovRand(mold, ndim)
5rand
6+1.00000000, -0.524043262
7-0.524043262, +1.00000000
9T
10rand = getCovRand(mold, ndim, scale)
11rand
12+16.0000000, +15.4355307
13+15.4355307, +15.9999981
15T
16
17
18rand = getCovRand(mold = 0._TKG, scale = [(real(itry, TKG), itry = 1, 5)])
19rand
20+0.999999881, +1.76461983, -0.623302639, -1.77094305, -0.371090323
21+1.76461983, +4.00000000, +1.30485487, -6.34312534, +2.20254874
22-0.623302639, +1.30485487, +9.00000000, -7.44724178, +11.8881130
23-1.77094305, -6.34312534, -7.44724178, +15.9999990, -6.73858166
24-0.371090323, +2.20254874, +11.8881130, -6.73858166, +24.9999962
26T
27
28
29ndim = getUnifRand(2, 3)
30scale = getUnifRand(1, 10)
31rand = getCovRand(mold, ndim)
32rand
33+1.00000000, -0.797388732
34-0.797388732, +1.00000000
36T
37rand = getCovRand(mold, ndim, scale)
38rand
39+4.00000000, +3.58696795
40+3.58696795, +4.00000048
42T
43
44
45rand = getCovRand(mold = 0._TKG, scale = [(real(itry, TKG), itry = 1, 5)])
46rand
47+1.00000000, +1.66718113, -0.281988472, -2.71459675, +2.36408281
48+1.66718113, +4.00000000, -0.181999475, -6.50552034, +5.91791439
49-0.281988472, -0.181999475, +8.99999905, +3.72426367, -8.28747845
50-2.71459675, -6.50552034, +3.72426367, +16.0000000, -11.6185513
51+2.36408281, +5.91791439, -8.28747845, -11.6185513, +25.0000019
53T
54
55
56ndim = getUnifRand(2, 3)
57scale = getUnifRand(1, 10)
58rand = getCovRand(mold, ndim)
59rand
60+1.00000000, +0.310113430
61+0.310113430, +1.00000000
63T
64rand = getCovRand(mold, ndim, scale)
65rand
66+4.00000000, -1.22807634
67-1.22807634, +4.00000048
69T
70
71
72rand = getCovRand(mold = 0._TKG, scale = [(real(itry, TKG), itry = 1, 5)])
73rand
74+1.00000000, +1.64157057, +0.275721699, -1.20346296, -3.76897526
75+1.64157057, +4.00000048, -2.91801953, -1.37555671, -9.67187881
76+0.275721699, -2.91801953, +9.00000000, -0.858826399, +8.06080151
77-1.20346296, -1.37555671, -0.858826399, +15.9999990, -0.250711441
78-3.76897526, -9.67187881, +8.06080151, -0.250711441, +24.9999981
80T
81
82
83ndim = getUnifRand(2, 3)
84scale = getUnifRand(1, 10)
85rand = getCovRand(mold, ndim)
86rand
87+1.00000000, +0.776790857, -0.935969293
88+0.776790857, +1.00000000, -0.947810173
89-0.935969293, -0.947810173, +1.00000000
91T
92rand = getCovRand(mold, ndim, scale)
93rand
94+24.9999962, +6.79447794, -18.1401615
95+6.79447794, +24.9999981, -16.6497517
96-18.1401615, -16.6497517, +24.9999981
98T
99
100
101rand = getCovRand(mold = 0._TKG, scale = [(real(itry, TKG), itry = 1, 5)])
102rand
103+1.00000000, -0.700695872, -0.297288060, +2.49565697, +1.86308646
104-0.700695872, +4.00000048, -3.89484358, +1.35728788, +3.22167921
105-0.297288060, -3.89484358, +9.00000191, -8.32446861, -11.4707747
106+2.49565697, +1.35728788, -8.32446861, +16.0000019, +16.1395988
107+1.86308646, +3.22167921, -11.4707747, +16.1395988, +25.0000000
109T
110
111
112ndim = getUnifRand(2, 3)
113scale = getUnifRand(1, 10)
114rand = getCovRand(mold, ndim)
115rand
116+1.00000000, +0.342483401, +0.630225241
117+0.342483401, +1.00000000, -7.10024166
118+0.630225241, -7.10024166, +1.00000000
120F
121rand = getCovRand(mold, ndim, scale)
122rand
123+4.00000000, +2.57892585, +0.963036656
124+2.57892585, +1112.85693, +52.1254807
125+0.963036656, +52.1254807, +4.00000000
127T
128
129
130rand = getCovRand(mold = 0._TKG, scale = [(real(itry, TKG), itry = 1, 5)])
131rand
132+0.999999881, -1.99783361, -2.03061795, -2.07726645, -0.552124619
133-1.99783361, +4.00000048, +4.10076046, +3.95902753, +0.820038557
134-2.03061795, +4.10076046, +9.00000191, -0.800278187, +4.46266890
135-2.07726645, +3.95902753, -0.800278187, +16.0000000, +9.42004013
136-0.552124619, +0.820038557, +4.46266890, +9.42004013, +25.0000000
138T
139
140
141ndim = getUnifRand(2, 3)
142scale = getUnifRand(1, 10)
143rand = getCovRand(mold, ndim)
144rand
145(+1.00000000, +0.00000000), (+0.627244771, +0.640735149)
146(+0.627244771, -0.640735149), (+1.00000000, +0.00000000)
148T
149rand = getCovRand(mold, ndim, scale)
150rand
151(+8.99999905, +0.00000000), (+5.18754864, +6.44286060)
152(+5.18754864, -6.44286060), (+8.99999905, +0.00000000)
154T
155
156
157rand = getCovRand(mold = 0._TKG, scale = [(real(itry, TKG), itry = 1, 5)])
158rand
159(+1.00000000, +0.00000000), (+0.334222376, +0.00000000), (-1.18526649, +0.00000000), (-2.28484535, +0.00000000), (+2.81603003, +0.00000000)
160(+0.334222376, +0.00000000), (+4.00000048, +0.00000000), (+2.86596727, +0.00000000), (-2.06729579, +0.00000000), (-3.07763791, +0.00000000)
161(-1.18526649, +0.00000000), (+2.86596727, +0.00000000), (+8.99999905, +0.00000000), (+7.84128141, +0.00000000), (-4.79801416, +0.00000000)
162(-2.28484535, +0.00000000), (-2.06729579, +0.00000000), (+7.84128141, +0.00000000), (+15.9999990, +0.00000000), (-6.78592443, +0.00000000)
163(+2.81603003, +0.00000000), (-3.07763791, +0.00000000), (-4.79801416, +0.00000000), (-6.78592443, +0.00000000), (+25.0000019, +0.00000000)
165T
166
167
168ndim = getUnifRand(2, 3)
169scale = getUnifRand(1, 10)
170rand = getCovRand(mold, ndim)
171rand
172(+1.00000000, +0.00000000), (+101.273308, -164.486343), (-2.97160268, +3.67032337)
173(+101.273308, +164.486343), (+1.00000000, +0.00000000), (-9.40372849, -1.74395859)
174(-2.97160268, -3.67032337), (-9.40372849, +1.74395859), (+1.00000000, +0.00000000)
176F
177rand = getCovRand(mold, ndim, scale)
178rand
179(+383.480499, +0.00000000), (+400.990356, -658.393127), (+25.0472527, +11.5311289)
180(+400.990356, +658.393127), (+1567.96289, +0.00000000), (+4.56126595, +49.6938362)
181(+25.0472527, -11.5311289), (+4.56126595, -49.6938362), (+3.99999976, +0.00000000)
183T
184
185
186rand = getCovRand(mold = 0._TKG, scale = [(real(itry, TKG), itry = 1, 5)])
187rand
188(+1.00000000, +0.00000000), (+1.10453939, +0.00000000), (+1.85971785, +0.00000000), (-2.36060834, +0.00000000), (-0.200386032, +0.00000000)
189(+1.10453939, +0.00000000), (+4.00000048, +0.00000000), (+4.96809816, +0.00000000), (+1.04312468, +0.00000000), (-4.53579235, +0.00000000)
190(+1.85971785, +0.00000000), (+4.96809816, +0.00000000), (+9.00000191, +0.00000000), (-0.323873609E-1, +0.00000000), (-0.401068479, +0.00000000)
191(-2.36060834, +0.00000000), (+1.04312468, +0.00000000), (-0.323873609E-1, +0.00000000), (+16.0000000, +0.00000000), (-1.61214328, +0.00000000)
192(-0.200386032, +0.00000000), (-4.53579235, +0.00000000), (-0.401068479, +0.00000000), (-1.61214328, +0.00000000), (+25.0000000, +0.00000000)
194T
195
196
197ndim = getUnifRand(2, 3)
198scale = getUnifRand(1, 10)
199rand = getCovRand(mold, ndim)
200rand
201(+1.00000000, +0.00000000), (+0.945548594, +0.380004644E-1)
202(+0.945548594, -0.380004644E-1), (+1.00000000, +0.00000000)
204T
205rand = getCovRand(mold, ndim, scale)
206rand
207(+81.0000000, +0.00000000), (+51.5637627, -39.3834076)
208(+51.5637627, +39.3834076), (+80.9999924, +0.00000000)
210T
211
212
213rand = getCovRand(mold = 0._TKG, scale = [(real(itry, TKG), itry = 1, 5)])
214rand
215(+1.00000000, +0.00000000), (+1.00697184, +0.00000000), (-2.10672903, +0.00000000), (+2.73027420, +0.00000000), (-3.62604022, +0.00000000)
216(+1.00697184, +0.00000000), (+4.00000000, +0.00000000), (-5.18921804, +0.00000000), (-2.00285459, +0.00000000), (-8.36297607, +0.00000000)
217(-2.10672903, +0.00000000), (-5.18921804, +0.00000000), (+8.99999905, +0.00000000), (+0.261529684, +0.00000000), (+13.0373650, +0.00000000)
218(+2.73027420, +0.00000000), (-2.00285459, +0.00000000), (+0.261529684, +0.00000000), (+16.0000019, +0.00000000), (-2.51530981, +0.00000000)
219(-3.62604022, +0.00000000), (-8.36297607, +0.00000000), (+13.0373650, +0.00000000), (-2.51530981, +0.00000000), (+24.9999981, +0.00000000)
221T
222
223
224ndim = getUnifRand(2, 3)
225scale = getUnifRand(1, 10)
226rand = getCovRand(mold, ndim)
227rand
228(+1.00000000, +0.00000000), (-0.300385416, +0.254224598)
229(-0.300385416, -0.254224598), (+1.00000000, +0.00000000)
231T
232rand = getCovRand(mold, ndim, scale)
233rand
234(+4.00000000, +0.00000000), (-3.33755755, +1.11532593)
235(-3.33755755, -1.11532593), (+3.99999976, +0.00000000)
237T
238
239
240rand = getCovRand(mold = 0._TKG, scale = [(real(itry, TKG), itry = 1, 5)])
241rand
242(+0.999999881, +0.00000000), (+1.96369064, +0.00000000), (+1.90418851, +0.00000000), (-1.38460803, +0.00000000), (+4.33383417, +0.00000000)
243(+1.96369064, +0.00000000), (+4.00000000, +0.00000000), (+4.57652998, +0.00000000), (-3.22220087, +0.00000000), (+8.92168617, +0.00000000)
244(+1.90418851, +0.00000000), (+4.57652998, +0.00000000), (+9.00000000, +0.00000000), (-3.59543681, +0.00000000), (+11.8303909, +0.00000000)
245(-1.38460803, +0.00000000), (-3.22220087, +0.00000000), (-3.59543681, +0.00000000), (+16.0000000, +0.00000000), (-3.41339946, +0.00000000)
246(+4.33383417, +0.00000000), (+8.92168617, +0.00000000), (+11.8303909, +0.00000000), (-3.41339946, +0.00000000), (+25.0000019, +0.00000000)
248T
249
250
251ndim = getUnifRand(2, 3)
252scale = getUnifRand(1, 10)
253rand = getCovRand(mold, ndim)
254rand
255(+1.00000000, +0.00000000), (+687.797241, -1192.57373), (-2.96659255, +15.1212234)
256(+687.797241, +1192.57373), (+1.00000000, +0.00000000), (-25.7348957, +9.24851418)
257(-2.96659255, -15.1212234), (-25.7348957, -9.24851418), (+1.00000000, +0.00000000)
259F
260rand = getCovRand(mold, ndim, scale)
261rand
262(+12181.3086, +0.00000000), (+10992.7461, -19067.9453), (-174.310532, -211.353210)
263(+10992.7461, +19067.9453), (+39860.5195, +0.00000000), (+160.786758, -451.156982)
264(-174.310532, +211.353210), (+160.786758, +451.156982), (+15.9999981, +0.00000000)
266T
267
268
269rand = getCovRand(mold = 0._TKG, scale = [(real(itry, TKG), itry = 1, 5)])
270rand
271(+1.00000000, +0.00000000), (-1.72044694, +0.00000000), (-1.67484498, +0.00000000), (+2.07150126, +0.00000000), (+0.399662018, +0.00000000)
272(-1.72044694, +0.00000000), (+3.99999952, +0.00000000), (+0.880984783, +0.00000000), (-6.19967222, +0.00000000), (-3.93389726, +0.00000000)
273(-1.67484498, +0.00000000), (+0.880984783, +0.00000000), (+9.00000000, +0.00000000), (-1.83047819, +0.00000000), (+8.52118683, +0.00000000)
274(+2.07150126, +0.00000000), (-6.19967222, +0.00000000), (-1.83047819, +0.00000000), (+16.0000000, +0.00000000), (+4.82812881, +0.00000000)
275(+0.399662018, +0.00000000), (-3.93389726, +0.00000000), (+8.52118683, +0.00000000), (+4.82812881, +0.00000000), (+25.0000000, +0.00000000)
277T
278
279
Test:
test_pm_distCov


Final Remarks


If you believe this algorithm or its documentation can be improved, we appreciate your contribution and help to edit this page's documentation and source file on GitHub.
For details on the naming abbreviations, see this page.
For details on the naming conventions, see this page.
This software is distributed under the MIT license with additional terms outlined below.

  1. If you use any parts or concepts from this library to any extent, please acknowledge the usage by citing the relevant publications of the ParaMonte library.
  2. If you regenerate any parts/ideas from this library in a programming environment other than those currently supported by this ParaMonte library (i.e., other than C, C++, Fortran, MATLAB, Python, R), please also ask the end users to cite this original ParaMonte library.

This software is available to the public under a highly permissive license.
Help us justify its continued development and maintenance by acknowledging its benefit to society, distributing it, and contributing to it.

Author:
Amir Shahmoradi, Monday March 6, 2017, 3:22 pm, Institute for Computational Engineering and Sciences (ICES), The University of Texas at Austin.

Definition at line 394 of file pm_distCov.F90.


The documentation for this interface was generated from the following file: