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

Return the covariance matrix corresponding to the input (potentially weighted) correlation matrix or return the biased sample covariance matrix of the input array of shape (ndim, nsam) or (nsam, ndim) or a pair of (potentially weighted) time series x(1:nsam) and y(1:nsam) where ndim is the number of data dimensions (the number of data attributes) and nsam is the number of data points.
More...

Detailed Description

Return the covariance matrix corresponding to the input (potentially weighted) correlation matrix or return the biased sample covariance matrix of the input array of shape (ndim, nsam) or (nsam, ndim) or a pair of (potentially weighted) time series x(1:nsam) and y(1:nsam) where ndim is the number of data dimensions (the number of data attributes) and nsam is the number of data points.

This generic interface performs either of the following two computational tasks:

  1. compute the covariance matrix corresponding to an input correlation matrix and vector of standard deviations in arbitrary ndim dimensions.
  2. compute the sample covariance matrix of a random sample of nsam observations in arbitrary ndim dimensional space.

See the documentation of the parent module pm_sampleCov for algorithmic details and sample covariance matrix definition.

Parameters
[in,out]cov: The output positive semi-definite square matrix of shape (1:ndim, 1:ndim) 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),
representing the covariance matrix corresponding to the input sample or cor, whichever is present.
On output,
  1. If the input arguments x and y are present, then ndim == 2 holds by definition and the output covariance is of the form

    \begin{equation} \ms{cov} = \begin{bmatrix} \sigma_{xx} && \sigma_{yx} \\ \sigma_{xy} && \sigma_{yy} \end{bmatrix} \end{equation}

  2. Otherwise, is subset is present, then only the specified input subset will be overwritten with the covariance matrix.
    Any elements not in the specified input subset remain intact.
[in]subset: The input scalar constant argument that can be any of the following:
  1. The constant lowDia, implying that only the lower-diagonal subset of the output covariance matrix must be computed.
  2. The constant uppDia, implying that only the upper-diagonal subset of the output covariance matrix must be computed.
This input argument is merely serves to resolve the different procedures of this generic interface from each other at compile-time.
(optional. It must be present if and only if the input arguments x and y are both missing.)
[in]cor: The input positive semi-definite square matrix of shape (1:ndim, 1:ndim) of the same type and kind as the input cov, representing the correlation matrix based upon which the output covariance matrix cov is to be computed.
(optional. It must be present if and only if the input arguments std and subsetr are present and the rest are missing.)
[in]subsetr: The input scalar constant argument that can be any of the following:
  1. The constant lowDia, implying that only the lower-diagonal subset of the input correlation matrix must be used.
  2. The constant uppDia, implying that only the upper-diagonal subset of the input correlation matrix must be used.
This input argument is merely serves to resolve the different procedures of this generic interface from each other at compile-time.
Although the allowed subset constants imply the use of the diagonal elements, they are by definition assumed to be 1 and therefore are not referenced in the algorithm.
(optional. It must be present if and only if the input arguments subset, cor, and std are present and the rest are missing.)
[in]std: The input positive vector of shape (1:ndim) of type real of the same kind as the input cov, containing the standard deviations of the ndim data attributes based upon which the output covariance matrix cov is to be computed.
(optional. It must be present if and only if the input arguments cor and subsetr are present and the rest are missing.)
[in]mean: The input contiguous vector of shape (ndim) of the same type and kind as the output cov containing the sample mean.
  1. If the input sample is missing and x and y are present, then mean must be a vector of size ndim = 2.
  2. If the input sample is present and is a 2D array, then mean must be a vector of size ndim = size(sample, 3 - dim) (i.e., computed along the specified input dimension dim).
The mean vector can be readily computed via getMean or setMean.
(optional. default = getFilled(0., ndim) or [0., 0.] depending on the rank of the input sample and dim or the presence of x and y.)
[in]x: The input contiguous vector of shape (nsam) of the same type and kind as the output cov, containing the first attribute x of the observational sample, where nsam is the number of observations in the sample.
(optional. It must be present if and only if the input argument y is present and cor and sample are missing.)
[in]y: The input contiguous vector of shape (nsam) of the same type and kind as the output cov, containing the second attribute x of the observational sample, where nsam is the number of observations in the sample.
(optional. It must be present if and only if the input argument x is present and cor and sample are missing.)
[in]sample: The input contiguous array of shape (ndim, nsam) or (nsam, ndim) of the same type and kind as the output cov, containing the sample comprised of nsam observations each with ndim attributes.
If sample is a matrix, then the direction along which the covariance matrix is computed is dictated by the input argument dim.
(optional. It must be present if and only if the input argument dim is present and x and y are missing.)
[in]dim: The input scalar integer of default kind IK indicating the dimension of sample along which the covariance matrix must be computed.
  1. If dim = 1, the input sample is assumed to have the shape (nsam, ndim).
  2. If dim = 2, the input sample is assumed to have the shape (ndim, nsam).
(optional. It must be present if and only if the input argument sample is present.)
[in]weight: The contiguous vector of length nsam of,
  1. type integer of default kind IK, or
  2. type real of the same kind as the kind of the output cov,
containing the corresponding weights of individual nsam observations in sample.
(optional. default = getFilled(1, nsam). It can be present if and only if the input arguments sample or x and y are present.)
[in]weisum: The input scalar of the same type and kind as the input weight containing sum(weight).
This quantity is a byproduct of computing the mean of a sample and is automatically returned by setMean.
(optional. It must be present if and only if the weight argument is present and the input arguments cor, subsetr, and std are missing.)


Possible calling interfaces

use pm_sampleCov, only: setCov, lowDia, uppDia
! correlation matrix to covariance matrix.
call setCov(cov(1:ndim, 1:ndim), subset, cor(1:ndim, 1:ndim), subsetr, std(1:ndim))
! XY time series covariance matrix.
call setCov(cov(1:2, 1:2) , x(1:nsam), y(1:nsam))
call setCov(cov(1:2, 1:2), mean(1:2), x(1:nsam), y(1:nsam))
call setCov(cov(1:2, 1:2) , x(1:nsam), y(1:nsam), weight(1:nsam))
call setCov(cov(1:2, 1:2), mean(1:2), x(1:nsam), y(1:nsam), weight(1:nsam), weisum)
! sample covariance matrix.
call setCov(cov(1:ndim, 1:ndim), subset , sample(:,:), dim)
call setCov(cov(1:ndim, 1:ndim), subset, mean(1:ndim), sample(:,:), dim)
call setCov(cov(1:ndim, 1:ndim), subset, , sample(:,:), dim, weight(1:nsam))
call setCov(cov(1:ndim, 1:ndim), subset, mean(1:ndim), sample(:,:), dim, weight(1:nsam), weisum)
Return the covariance matrix corresponding to the input (potentially weighted) correlation matrix or ...
This module contains classes and procedures for computing the properties related to the covariance ma...
Warning
The condition any(x(1) /= x) must hold for the corresponding input arguments.
The condition any(y(1) /= x) must hold for the corresponding input arguments.
The input sample must contain at least two unique values per sample attribute.
The condition 0 < sum(weight) must hold for the corresponding input arguments.
The condition all(0. <= weight) must hold for the corresponding input arguments.
The condition 1 < size(sample, dim) must hold for the corresponding input arguments.
The condition 0 < size(sample, 3 - dim) must hold for the corresponding input arguments.
The condition all(shape(cov) == shape(cor)) must hold for the corresponding input arguments.
The condition 1 <= dim .and. dim <= rank(sample) must hold for the corresponding input arguments.
The condition size(sample, 3 - dim) == size(variance) must hold for the corresponding input arguments.
The condition size(sample, 3 - dim) == size(mean) must hold for the corresponding input arguments.
The condition size(sample, dim) == size(weight) must hold for the corresponding input arguments.
The condition size(x) == size(weight) must hold for the corresponding input arguments.
The condition size(x) == size(y) must hold for the corresponding input arguments.
These conditions are verified only if the library is built with the preprocessor macro CHECK_ENABLED=1.
The pure procedure(s) documented herein become impure when the ParaMonte library is compiled with preprocessor macro CHECK_ENABLED=1.
By default, these procedures are pure in release build and impure in debug and testing builds.
Note
Note the effects of bias-correction in computing the covariance matrix become noticeable only for very sample sample sizes (i.e., when nsam is small).
See also
getCor
setCor
getCov
setCov
getVar
setVar
getMean
setMean
getCovMerged
setCovMerged
setCovMeanMerged
getVarCorrection
getMeanMerged
setMeanMerged
getVarMerged
setVarMerged


Example usage

1program example
2
3 use pm_kind, only: SK, IK, LK
4 use pm_sampleMean, only: getMean
5 use pm_sampleMean, only: setMean
9 use pm_distUnif, only: getUnifRand
10 use pm_arrayFill, only: getFilled
11 use pm_sampleCov, only: uppDia
12 use pm_sampleCov, only: lowDia
13 use pm_sampleCov, only: setCov
14 use pm_io, only: display_type
15 use pm_io, only: getFormat
16
17 implicit none
18
19 integer(IK) :: itry, ntry = 10
20 type(display_type) :: disp
21 character(:), allocatable :: format
22 disp = display_type(file = "main.out.F90")
23
24 call disp%skip()
25 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
26 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
27 call disp%show("!Convert correlation matrix and standard deviation to covariance matrix.")
28 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
29 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
30 call disp%skip()
31
32 block
33 use pm_kind, only: TKG => RKS ! All other real types are also supported.
34 use pm_matrixCopy, only: setMatCopy, rdpack
35 use pm_distCov, only: getCovRand
36 integer(IK) :: ndim
37 real(TKG), allocatable :: cov(:,:), cor(:,:), std(:)
38 format = getFormat(mold = [0._TKG], ed = SK_"es", signed = .true._LK)
39 do itry = 1, ntry
40 call disp%skip()
41 call disp%show("ndim = getUnifRand(1, 7)")
42 ndim = getUnifRand(1, 7)
43 call disp%show("ndim")
44 call disp%show( ndim )
45 call disp%show("std = getUnifRand(1, 10, ndim)")
46 std = getUnifRand(1, 10, ndim)
47 call disp%show("std")
48 call disp%show( std , format = format )
49 call disp%show("cor = getCovRand(1., ndim)")
50 cor = getCovRand(1., ndim)
51 call disp%show("cor")
52 call disp%show( cor , format = format )
53 call disp%show("cov = getFilled(0._TKG, ndim, ndim + 1)")
54 cov = getFilled(0._TKG, ndim, ndim + 1)
55 call disp%show("call setMatCopy(cov(1:ndim, 1:ndim), rdpack, cor, rdpack, lowDia)")
56 call setMatCopy(cov(1:ndim, 1:ndim), rdpack, cor, rdpack, lowDia)
57 call disp%show("cov")
58 call disp%show( cov , format = format )
59 call disp%show("call setCov(cov(1:ndim, 2:ndim+1), uppDia, cov(1:ndim, 1:ndim), lowDia, std) ! convert lower-triangle correlation matrix to upper-triangle covariance matrix.")
60 call setCov(cov(1:ndim, 2:ndim+1), uppDia, cov(1:ndim, 1:ndim), lowDia, std)
61 call disp%show("cov")
62 call disp%show( cov , format = format )
63 call disp%skip()
64 call disp%show("cov = getFilled(0._TKG, ndim, ndim + 1)")
65 cov = getFilled(0._TKG, ndim, ndim + 1)
66 call disp%show("call setMatCopy(cov(1:ndim, 2:ndim+1), rdpack, cor, rdpack, uppDia)")
67 call setMatCopy(cov(1:ndim, 2:ndim+1), rdpack, cor, rdpack, uppDia)
68 call disp%show("cov")
69 call disp%show( cov , format = format )
70 call disp%show("call setCov(cov(1:ndim, 1:ndim), lowDia, cov(1:ndim, 2:ndim+1), uppDia, std) ! convert upper-triangle correlation matrix to lower-triangle covariance matrix.")
71 call setCov(cov(1:ndim, 1:ndim), lowDia, cov(1:ndim, 2:ndim+1), uppDia, std)
72 call disp%show("cov")
73 call disp%show( cov , format = format )
74 call disp%skip()
75 end do
76 end block
77
78 call disp%skip()
79 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
80 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
81 call disp%show("!Compute the covariance matrix of a 2-D sample.")
82 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
83 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
84 call disp%skip()
85
86 block
87 use pm_kind, only: TKG => RKS ! All other real types are also supported.
88 real(TKG), allocatable :: sample(:,:), samShifted(:,:), cov(:,:), mean(:)
89 integer(IK) :: ndim, nsam, dim
90 call disp%show("ndim = 2; nsam = 10; dim = 2")
91 ndim = 2; nsam = 10; dim = 2
92 call disp%show("sample = reshape(getUnifRand(1, 20, ndim * nsam), shape = [ndim, nsam], order = [2, 1])")
93 sample = reshape(getUnifRand(1, 20, ndim * nsam), shape = [ndim, nsam], order = [2, 1])
94 call disp%show("sample")
95 call disp%show( sample )
96 call disp%show("mean = getMean(sample, dim)")
97 mean = getMean(sample, dim)
98 call disp%show("mean")
99 call disp%show( mean )
100 call disp%show("samShifted = getShifted(sample, dim, -mean)")
101 samShifted = getShifted(sample, dim, -mean)
102 call disp%show("samShifted")
103 call disp%show( samShifted )
104 call disp%show("cov = getFilled(0., ndim, ndim)")
105 cov = getFilled(0., ndim, ndim)
106 call disp%show("call setCov(cov, uppDia, mean, sample, dim)")
107 call setCov(cov, uppDia, mean, sample, dim)
108 call disp%show("cov")
109 call disp%show( cov )
110 call disp%show("cov = getFilled(0., ndim, ndim)")
111 cov = getFilled(0., ndim, ndim)
112 call disp%show("call setCov(cov, uppDia, samShifted, dim) ! same result as above.")
113 call setCov(cov, uppDia, samShifted, dim)
114 call disp%show("cov")
115 call disp%show( cov )
116 call disp%skip()
117 call disp%show("cov = getFilled(0., ndim, ndim)")
118 cov = getFilled(0., ndim, ndim)
119 call disp%show("call setCov(cov, lowDia, mean, sample, dim)")
120 call setCov(cov, lowDia, mean, sample, dim)
121 call disp%show("cov")
122 call disp%show( cov )
123 call disp%show("cov = getFilled(0., ndim, ndim)")
124 cov = getFilled(0., ndim, ndim)
125 call disp%show("call setCov(cov, lowDia, samShifted, dim) ! same result as above.")
126 call setCov(cov, lowDia, samShifted, dim)
127 call disp%show("cov")
128 call disp%show( cov )
129 call disp%skip()
130 call disp%show("Compute the sample covariance along the first dimension.", deliml = SK_'''')
131 call disp%skip()
132 call disp%show("dim = 1")
133 dim = 1
134 call disp%show("cov = getFilled(0., ndim, ndim)")
135 cov = getFilled(0., ndim, ndim)
136 call disp%show("call setCov(cov, uppDia, mean, transpose(sample), dim)")
137 call setCov(cov, uppDia, mean, transpose(sample), dim)
138 call disp%show("cov")
139 call disp%show( cov )
140 call disp%show("cov = getFilled(0., ndim, ndim)")
141 cov = getFilled(0., ndim, ndim)
142 call disp%show("call setCov(cov, uppDia, transpose(samShifted), dim) ! same result as above.")
143 call setCov(cov, uppDia, transpose(samShifted), dim)
144 call disp%show("cov")
145 call disp%show( cov )
146 call disp%skip()
147 call disp%show("cov = getFilled(0., ndim, ndim)")
148 cov = getFilled(0., ndim, ndim)
149 call disp%show("call setCov(cov, lowDia, mean, transpose(sample), dim)")
150 call setCov(cov, lowDia, mean, transpose(sample), dim)
151 call disp%show("cov")
152 call disp%show( cov )
153 call disp%show("cov = getFilled(0., ndim, ndim)")
154 cov = getFilled(0., ndim, ndim)
155 call disp%show("call setCov(cov, lowDia, transpose(samShifted), dim) ! same result as above.")
156 call setCov(cov, lowDia, transpose(samShifted), dim)
157 call disp%show("cov")
158 call disp%show( cov )
159 call disp%skip()
160 call disp%show("Compute the full sample covariance for a pair of time series.", deliml = SK_'''')
161 call disp%skip()
162 call disp%show("cov = getFilled(0., ndim, ndim)")
163 cov = getFilled(0., ndim, ndim)
164 call disp%show("call setCov(cov, mean, sample(1,:), sample(2,:))")
165 call setCov(cov, mean, sample(1,:), sample(2,:))
166 call disp%show("cov")
167 call disp%show( cov )
168 call disp%show("cov = getFilled(0., ndim, ndim)")
169 cov = getFilled(0., ndim, ndim)
170 call disp%show("call setCov(cov, samShifted(1,:), samShifted(2,:)) ! same result as above.")
171 call setCov(cov, samShifted(1,:), samShifted(2,:))
172 call disp%show("cov")
173 call disp%show( cov )
174 call disp%skip()
175 end block
176
177 call disp%skip()
178 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
179 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
180 call disp%show("!Compute the biased covariance matrix of a weighted 2-D sample.")
181 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
182 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
183 call disp%skip()
184
185 block
186 use pm_kind, only: TKG => RKS ! All other real types are also supported.
187 real(TKG) :: rweisum
188 integer(IK) :: iweisum
189 real(TKG), allocatable :: rweight(:)
190 integer(IK), allocatable :: iweight(:)
191 real(TKG), allocatable :: sample(:,:), samShifted(:,:), cov(:,:), mean(:)
192 integer(IK) :: ndim, nsam, dim
193 call disp%show("ndim = 2; nsam = 10; dim = 2")
194 ndim = 2; nsam = 10; dim = 2
195 call disp%show("sample = reshape(getUnifRand(1, 20, ndim * nsam), shape = [ndim, nsam], order = [2, 1])")
196 sample = reshape(getUnifRand(1, 20, ndim * nsam), shape = [ndim, nsam], order = [2, 1])
197 call disp%show("sample")
198 call disp%show( sample )
199 call disp%show("call setResized(mean, ndim)")
200 call setResized(mean, ndim)
201 call disp%show("iweight = getUnifRand(1, 10, nsam) ! integer-valued weights.")
202 iweight = getUnifRand(1, 10, nsam) ! integer-valued weights.
203 call disp%show("iweight")
204 call disp%show( iweight )
205 call disp%show("call setMean(mean, sample, dim, iweight, iweisum)")
206 call setMean(mean, sample, dim, iweight, iweisum)
207 call disp%show("mean")
208 call disp%show( mean )
209 call disp%show("iweisum")
210 call disp%show( iweisum )
211 call disp%show("rweight = iweight ! or real-valued weights.")
212 rweight = iweight ! or real-valued weights.
213 call disp%show("iweight")
214 call disp%show( iweight )
215 call disp%show("call setMean(mean, sample, dim, rweight, rweisum)")
216 call setMean(mean, sample, dim, rweight, rweisum)
217 call disp%show("mean")
218 call disp%show( mean )
219 call disp%show("rweisum")
220 call disp%show( rweisum )
221
222 call disp%skip()
223 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
224 call disp%show("!Compute the covariance matrix integer weights.")
225 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
226 call disp%skip()
227
228 call disp%show("samShifted = getShifted(sample, dim, -mean)")
229 samShifted = getShifted(sample, dim, -mean)
230 call disp%show("samShifted")
231 call disp%show( samShifted )
232 call disp%show("cov = getFilled(0., ndim, ndim)")
233 cov = getFilled(0., ndim, ndim)
234 call disp%show("call setCov(cov, uppDia, mean, sample, dim, iweight, iweisum)")
235 call setCov(cov, uppDia, mean, sample, dim, iweight, iweisum)
236 call disp%show("cov")
237 call disp%show( cov )
238 call disp%show("cov = getFilled(0., ndim, ndim)")
239 cov = getFilled(0., ndim, ndim)
240 call disp%show("call setCov(cov, uppDia, samShifted, dim, iweight, sum(iweight)) ! same result as above.")
241 call setCov(cov, uppDia, samShifted, dim, iweight, sum(iweight))
242 call disp%show("cov")
243 call disp%show( cov )
244 call disp%skip()
245 call disp%show("cov = getFilled(0., ndim, ndim)")
246 cov = getFilled(0., ndim, ndim)
247 call disp%show("call setCov(cov, lowDia, mean, sample, dim, iweight, iweisum)")
248 call setCov(cov, lowDia, mean, sample, dim, iweight, iweisum)
249 call disp%show("cov")
250 call disp%show( cov )
251 call disp%show("cov = getFilled(0., ndim, ndim)")
252 cov = getFilled(0., ndim, ndim)
253 call disp%show("call setCov(cov, lowDia, samShifted, dim, iweight, sum(iweight)) ! same result as above.")
254 call setCov(cov, lowDia, samShifted, dim, iweight, sum(iweight))
255 call disp%show("cov")
256 call disp%show( cov )
257 call disp%skip()
258 call disp%show("Compute the sample covariance along the first dimension.", deliml = SK_'''')
259 call disp%skip()
260 call disp%show("dim = 1")
261 dim = 1
262 call disp%show("cov = getFilled(0., ndim, ndim)")
263 cov = getFilled(0., ndim, ndim)
264 call disp%show("call setCov(cov, uppDia, mean, transpose(sample), dim, iweight, iweisum)")
265 call setCov(cov, uppDia, mean, transpose(sample), dim, iweight, iweisum)
266 call disp%show("cov")
267 call disp%show( cov )
268 call disp%show("cov = getFilled(0., ndim, ndim)")
269 cov = getFilled(0., ndim, ndim)
270 call disp%show("call setCov(cov, uppDia, transpose(samShifted), dim, iweight, sum(iweight)) ! same result as above.")
271 call setCov(cov, uppDia, transpose(samShifted), dim, iweight, sum(iweight))
272 call disp%show("cov")
273 call disp%show( cov )
274 call disp%skip()
275 call disp%show("cov = getFilled(0., ndim, ndim)")
276 cov = getFilled(0., ndim, ndim)
277 call disp%show("call setCov(cov, lowDia, mean, transpose(sample), dim, iweight, iweisum)")
278 call setCov(cov, lowDia, mean, transpose(sample), dim, iweight, iweisum)
279 call disp%show("cov")
280 call disp%show( cov )
281 call disp%show("cov = getFilled(0., ndim, ndim)")
282 cov = getFilled(0., ndim, ndim)
283 call disp%show("call setCov(cov, lowDia, transpose(samShifted), dim, iweight, sum(iweight)) ! same result as above.")
284 call setCov(cov, lowDia, transpose(samShifted), dim, iweight, sum(iweight))
285 call disp%show("cov")
286 call disp%show( cov )
287 call disp%skip()
288 call disp%show("Compute the full sample covariance for a pair of time series.", deliml = SK_'''')
289 call disp%skip()
290 call disp%show("cov = getFilled(0., ndim, ndim)")
291 cov = getFilled(0., ndim, ndim)
292 call disp%show("call setCov(cov, mean, sample(1,:), sample(2,:), iweight, iweisum)")
293 call setCov(cov, mean, sample(1,:), sample(2,:), iweight, iweisum)
294 call disp%show("cov")
295 call disp%show( cov )
296 call disp%show("cov = getFilled(0., ndim, ndim)")
297 cov = getFilled(0., ndim, ndim)
298 call disp%show("call setCov(cov, samShifted(1,:), samShifted(2,:), iweight, sum(iweight)) ! same result as above.")
299 call setCov(cov, samShifted(1,:), samShifted(2,:), iweight, sum(iweight))
300 call disp%show("cov")
301 call disp%show( cov )
302 call disp%skip()
303
304 call disp%skip()
305 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
306 call disp%show("!Compute the covariance matrix real weights.")
307 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
308 call disp%skip()
309
310 call disp%show("dim = 2_IK")
311 dim = 2_IK
312 call disp%show("samShifted = getShifted(sample, dim, -mean)")
313 samShifted = getShifted(sample, dim, -mean)
314 call disp%show("samShifted")
315 call disp%show( samShifted )
316 call disp%show("cov = getFilled(0., ndim, ndim)")
317 cov = getFilled(0., ndim, ndim)
318 call disp%show("call setCov(cov, uppDia, mean, sample, dim, rweight, rweisum)")
319 call setCov(cov, uppDia, mean, sample, dim, rweight, rweisum)
320 call disp%show("cov")
321 call disp%show( cov )
322 call disp%show("cov = getFilled(0., ndim, ndim)")
323 cov = getFilled(0., ndim, ndim)
324 call disp%show("call setCov(cov, uppDia, samShifted, dim, rweight, sum(rweight)) ! same result as above.")
325 call setCov(cov, uppDia, samShifted, dim, rweight, sum(rweight))
326 call disp%show("cov")
327 call disp%show( cov )
328 call disp%skip()
329 call disp%show("cov = getFilled(0., ndim, ndim)")
330 cov = getFilled(0., ndim, ndim)
331 call disp%show("call setCov(cov, lowDia, mean, sample, dim, rweight, rweisum)")
332 call setCov(cov, lowDia, mean, sample, dim, rweight, rweisum)
333 call disp%show("cov")
334 call disp%show( cov )
335 call disp%show("cov = getFilled(0., ndim, ndim)")
336 cov = getFilled(0., ndim, ndim)
337 call disp%show("call setCov(cov, lowDia, samShifted, dim, rweight, sum(rweight)) ! same result as above.")
338 call setCov(cov, lowDia, samShifted, dim, rweight, sum(rweight))
339 call disp%show("cov")
340 call disp%show( cov )
341 call disp%skip()
342 call disp%show("Compute the sample covariance along the first dimension.", deliml = SK_'''')
343 call disp%skip()
344 call disp%show("dim = 1")
345 dim = 1
346 call disp%show("cov = getFilled(0., ndim, ndim)")
347 cov = getFilled(0., ndim, ndim)
348 call disp%show("call setCov(cov, uppDia, mean, transpose(sample), dim, rweight, rweisum)")
349 call setCov(cov, uppDia, mean, transpose(sample), dim, rweight, rweisum)
350 call disp%show("cov")
351 call disp%show( cov )
352 call disp%show("cov = getFilled(0., ndim, ndim)")
353 cov = getFilled(0., ndim, ndim)
354 call disp%show("call setCov(cov, uppDia, transpose(samShifted), dim, rweight, sum(rweight)) ! same result as above.")
355 call setCov(cov, uppDia, transpose(samShifted), dim, rweight, sum(rweight))
356 call disp%show("cov")
357 call disp%show( cov )
358 call disp%skip()
359 call disp%show("cov = getFilled(0., ndim, ndim)")
360 cov = getFilled(0., ndim, ndim)
361 call disp%show("call setCov(cov, lowDia, mean, transpose(sample), dim, rweight, rweisum)")
362 call setCov(cov, lowDia, mean, transpose(sample), dim, rweight, rweisum)
363 call disp%show("cov")
364 call disp%show( cov )
365 call disp%show("cov = getFilled(0., ndim, ndim)")
366 cov = getFilled(0., ndim, ndim)
367 call disp%show("call setCov(cov, lowDia, transpose(samShifted), dim, rweight, sum(rweight)) ! same result as above.")
368 call setCov(cov, lowDia, transpose(samShifted), dim, rweight, sum(rweight))
369 call disp%show("cov")
370 call disp%show( cov )
371 call disp%skip()
372 call disp%show("Compute the full sample covariance for a pair of time series.", deliml = SK_'''')
373 call disp%skip()
374 call disp%show("cov = getFilled(0., ndim, ndim)")
375 cov = getFilled(0., ndim, ndim)
376 call disp%show("call setCov(cov, mean, sample(1,:), sample(2,:), rweight, rweisum)")
377 call setCov(cov, mean, sample(1,:), sample(2,:), rweight, rweisum)
378 call disp%show("cov")
379 call disp%show( cov )
380 call disp%show("cov = getFilled(0., ndim, ndim)")
381 cov = getFilled(0., ndim, ndim)
382 call disp%show("call setCov(cov, samShifted(1,:), samShifted(2,:), rweight, sum(rweight)) ! same result as above.")
383 call setCov(cov, samShifted(1,:), samShifted(2,:), rweight, sum(rweight))
384 call disp%show("cov")
385 call disp%show( cov )
386 call disp%skip()
387 end block
388
389end program example
Generate and return an array of the specified rank and shape of arbitrary intrinsic type and kind wit...
Allocate or resize (shrink or expand) an input allocatable scalar string or array of rank 1....
Generate count evenly spaced points over the interval [x1, x2] if x1 < x2, or [x2,...
Generate and return a random positive-definite (correlation or covariance) matrix using the Gram meth...
Definition: pm_distCov.F90:394
Generate and return a scalar or a contiguous array of rank 1 of length s1 of randomly uniformly distr...
Generate and return a generic or type/kind-specific IO format with the requested specifications that ...
Definition: pm_io.F90:18485
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
Copy a desired subset of the input source matrix of arbitrary shape (:) or (:,:) to the target subset...
Generate and return the (weighted) mean of an input sample of nsam observations with ndim = 1 or 2 at...
Return the (weighted) mean of a pair of time series or of an input sample of nsam observations with n...
Generate a sample of shape (nsam), or (ndim, nsam) or (nsam, ndim) that is shifted by the specified i...
This module contains procedures and generic interfaces for convenient allocation and filling of array...
This module contains procedures and generic interfaces for resizing allocatable arrays of various typ...
This module contains procedures and generic interfaces for generating arrays with linear or logarithm...
This module contains classes and procedures for generating random matrices distributed on the space o...
Definition: pm_distCov.F90:72
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 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 relevant to copying (diagonal or upper/lower t...
This module contains classes and procedures for computing the first moment (i.e., the statistical mea...
This module contains classes and procedures for shifting univariate or multivariate samples by arbitr...
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
2!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
3!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
4!Convert correlation matrix and standard deviation to covariance matrix.
5!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
6!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
7
8
9ndim = getUnifRand(1, 7)
10ndim
11+2
12std = getUnifRand(1, 10, ndim)
13std
14+4.000000E+00, +2.000000E+00
15cor = getCovRand(1., ndim)
16cor
17+1.000000E+00, -7.630376E-01
18-7.630376E-01, +1.000000E+00
19cov = getFilled(0._TKG, ndim, ndim + 1)
20call setMatCopy(cov(1:ndim, 1:ndim), rdpack, cor, rdpack, lowDia)
21cov
22+1.000000E+00, +0.000000E+00, +0.000000E+00
23-7.630376E-01, +1.000000E+00, +0.000000E+00
24call setCov(cov(1:ndim, 2:ndim+1), uppDia, cov(1:ndim, 1:ndim), lowDia, std) ! convert lower-triangle correlation matrix to upper-triangle covariance matrix.
25cov
26+1.000000E+00, +1.600000E+01, -6.104300E+00
27-7.630376E-01, +1.000000E+00, +4.000000E+00
28
29cov = getFilled(0._TKG, ndim, ndim + 1)
30call setMatCopy(cov(1:ndim, 2:ndim+1), rdpack, cor, rdpack, uppDia)
31cov
32+0.000000E+00, +1.000000E+00, -7.630376E-01
33+0.000000E+00, +0.000000E+00, +1.000000E+00
34call setCov(cov(1:ndim, 1:ndim), lowDia, cov(1:ndim, 2:ndim+1), uppDia, std) ! convert upper-triangle correlation matrix to lower-triangle covariance matrix.
35cov
36+1.600000E+01, +1.000000E+00, -7.630376E-01
37-6.104300E+00, +4.000000E+00, +1.000000E+00
38
39
40ndim = getUnifRand(1, 7)
41ndim
42+4
43std = getUnifRand(1, 10, ndim)
44std
45+5.000000E+00, +5.000000E+00, +1.000000E+00, +1.000000E+01
46cor = getCovRand(1., ndim)
47cor
48+1.000000E+00, +2.005533E-01, -1.883450E-01, +1.922954E-01
49+2.005533E-01, +1.000000E+00, -6.863858E-01, +2.038796E-02
50-1.883450E-01, -6.863858E-01, +1.000000E+00, -7.081194E-01
51+1.922954E-01, +2.038796E-02, -7.081194E-01, +1.000000E+00
52cov = getFilled(0._TKG, ndim, ndim + 1)
53call setMatCopy(cov(1:ndim, 1:ndim), rdpack, cor, rdpack, lowDia)
54cov
55+1.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00
56+2.005533E-01, +1.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00
57-1.883450E-01, -6.863858E-01, +1.000000E+00, +0.000000E+00, +0.000000E+00
58+1.922954E-01, +2.038796E-02, -7.081194E-01, +1.000000E+00, +0.000000E+00
59call setCov(cov(1:ndim, 2:ndim+1), uppDia, cov(1:ndim, 1:ndim), lowDia, std) ! convert lower-triangle correlation matrix to upper-triangle covariance matrix.
60cov
61+1.000000E+00, +2.500000E+01, +5.013832E+00, -9.417253E-01, +9.614773E+00
62+2.005533E-01, +1.000000E+00, +2.500000E+01, -3.431929E+00, +1.019398E+00
63-1.883450E-01, -6.863858E-01, +1.000000E+00, +1.000000E+00, -7.081194E+00
64+1.922954E-01, +2.038796E-02, -7.081194E-01, +1.000000E+00, +1.000000E+02
65
66cov = getFilled(0._TKG, ndim, ndim + 1)
67call setMatCopy(cov(1:ndim, 2:ndim+1), rdpack, cor, rdpack, uppDia)
68cov
69+0.000000E+00, +1.000000E+00, +2.005533E-01, -1.883450E-01, +1.922954E-01
70+0.000000E+00, +0.000000E+00, +1.000000E+00, -6.863858E-01, +2.038796E-02
71+0.000000E+00, +0.000000E+00, +0.000000E+00, +1.000000E+00, -7.081194E-01
72+0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +1.000000E+00
73call setCov(cov(1:ndim, 1:ndim), lowDia, cov(1:ndim, 2:ndim+1), uppDia, std) ! convert upper-triangle correlation matrix to lower-triangle covariance matrix.
74cov
75+2.500000E+01, +1.000000E+00, +2.005533E-01, -1.883450E-01, +1.922954E-01
76+5.013832E+00, +2.500000E+01, +1.000000E+00, -6.863858E-01, +2.038796E-02
77-9.417253E-01, -3.431929E+00, +1.000000E+00, +1.000000E+00, -7.081194E-01
78+9.614773E+00, +1.019398E+00, -7.081194E+00, +1.000000E+02, +1.000000E+00
79
80
81ndim = getUnifRand(1, 7)
82ndim
83+3
84std = getUnifRand(1, 10, ndim)
85std
86+5.000000E+00, +2.000000E+00, +4.000000E+00
87cor = getCovRand(1., ndim)
88cor
89+1.000000E+00, +5.682755E-02, -2.933485E-01
90+5.682755E-02, +1.000000E+00, +6.507880E-01
91-2.933485E-01, +6.507880E-01, +1.000000E+00
92cov = getFilled(0._TKG, ndim, ndim + 1)
93call setMatCopy(cov(1:ndim, 1:ndim), rdpack, cor, rdpack, lowDia)
94cov
95+1.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00
96+5.682755E-02, +1.000000E+00, +0.000000E+00, +0.000000E+00
97-2.933485E-01, +6.507880E-01, +1.000000E+00, +0.000000E+00
98call setCov(cov(1:ndim, 2:ndim+1), uppDia, cov(1:ndim, 1:ndim), lowDia, std) ! convert lower-triangle correlation matrix to upper-triangle covariance matrix.
99cov
100+1.000000E+00, +2.500000E+01, +5.682755E-01, -5.866969E+00
101+5.682755E-02, +1.000000E+00, +4.000000E+00, +5.206304E+00
102-2.933485E-01, +6.507880E-01, +1.000000E+00, +1.600000E+01
103
104cov = getFilled(0._TKG, ndim, ndim + 1)
105call setMatCopy(cov(1:ndim, 2:ndim+1), rdpack, cor, rdpack, uppDia)
106cov
107+0.000000E+00, +1.000000E+00, +5.682755E-02, -2.933485E-01
108+0.000000E+00, +0.000000E+00, +1.000000E+00, +6.507880E-01
109+0.000000E+00, +0.000000E+00, +0.000000E+00, +1.000000E+00
110call setCov(cov(1:ndim, 1:ndim), lowDia, cov(1:ndim, 2:ndim+1), uppDia, std) ! convert upper-triangle correlation matrix to lower-triangle covariance matrix.
111cov
112+2.500000E+01, +1.000000E+00, +5.682755E-02, -2.933485E-01
113+5.682755E-01, +4.000000E+00, +1.000000E+00, +6.507880E-01
114-5.866969E+00, +5.206304E+00, +1.600000E+01, +1.000000E+00
115
116
117ndim = getUnifRand(1, 7)
118ndim
119+2
120std = getUnifRand(1, 10, ndim)
121std
122+7.000000E+00, +8.000000E+00
123cor = getCovRand(1., ndim)
124cor
125+1.000000E+00, +6.787896E-01
126+6.787896E-01, +1.000000E+00
127cov = getFilled(0._TKG, ndim, ndim + 1)
128call setMatCopy(cov(1:ndim, 1:ndim), rdpack, cor, rdpack, lowDia)
129cov
130+1.000000E+00, +0.000000E+00, +0.000000E+00
131+6.787896E-01, +1.000000E+00, +0.000000E+00
132call setCov(cov(1:ndim, 2:ndim+1), uppDia, cov(1:ndim, 1:ndim), lowDia, std) ! convert lower-triangle correlation matrix to upper-triangle covariance matrix.
133cov
134+1.000000E+00, +4.900000E+01, +3.801222E+01
135+6.787896E-01, +1.000000E+00, +6.400000E+01
136
137cov = getFilled(0._TKG, ndim, ndim + 1)
138call setMatCopy(cov(1:ndim, 2:ndim+1), rdpack, cor, rdpack, uppDia)
139cov
140+0.000000E+00, +1.000000E+00, +6.787896E-01
141+0.000000E+00, +0.000000E+00, +1.000000E+00
142call setCov(cov(1:ndim, 1:ndim), lowDia, cov(1:ndim, 2:ndim+1), uppDia, std) ! convert upper-triangle correlation matrix to lower-triangle covariance matrix.
143cov
144+4.900000E+01, +1.000000E+00, +6.787896E-01
145+3.801222E+01, +6.400000E+01, +1.000000E+00
146
147
148ndim = getUnifRand(1, 7)
149ndim
150+5
151std = getUnifRand(1, 10, ndim)
152std
153+1.000000E+01, +7.000000E+00, +1.000000E+01, +9.000000E+00, +1.000000E+00
154cor = getCovRand(1., ndim)
155cor
156+1.000000E+00, +1.674978E-01, +3.232417E-01, -2.522079E-02, +4.693460E-01
157+1.674978E-01, +1.000000E+00, -7.061831E-01, -5.579556E-01, +1.891168E-01
158+3.232417E-01, -7.061831E-01, +1.000000E+00, +6.764400E-01, -2.016856E-01
159-2.522079E-02, -5.579556E-01, +6.764400E-01, +1.000000E+00, -7.950259E-01
160+4.693460E-01, +1.891168E-01, -2.016856E-01, -7.950259E-01, +1.000000E+00
161cov = getFilled(0._TKG, ndim, ndim + 1)
162call setMatCopy(cov(1:ndim, 1:ndim), rdpack, cor, rdpack, lowDia)
163cov
164+1.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00
165+1.674978E-01, +1.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00
166+3.232417E-01, -7.061831E-01, +1.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00
167-2.522079E-02, -5.579556E-01, +6.764400E-01, +1.000000E+00, +0.000000E+00, +0.000000E+00
168+4.693460E-01, +1.891168E-01, -2.016856E-01, -7.950259E-01, +1.000000E+00, +0.000000E+00
169call setCov(cov(1:ndim, 2:ndim+1), uppDia, cov(1:ndim, 1:ndim), lowDia, std) ! convert lower-triangle correlation matrix to upper-triangle covariance matrix.
170cov
171+1.000000E+00, +1.000000E+02, +1.172485E+01, +3.232417E+01, -2.269871E+00, +4.693460E+00
172+1.674978E-01, +1.000000E+00, +4.900000E+01, -4.943282E+01, -3.515120E+01, +1.323818E+00
173+3.232417E-01, -7.061831E-01, +1.000000E+00, +1.000000E+02, +6.087960E+01, -2.016856E+00
174-2.522079E-02, -5.579556E-01, +6.764400E-01, +1.000000E+00, +8.100000E+01, -7.155233E+00
175+4.693460E-01, +1.891168E-01, -2.016856E-01, -7.950259E-01, +1.000000E+00, +1.000000E+00
176
177cov = getFilled(0._TKG, ndim, ndim + 1)
178call setMatCopy(cov(1:ndim, 2:ndim+1), rdpack, cor, rdpack, uppDia)
179cov
180+0.000000E+00, +1.000000E+00, +1.674978E-01, +3.232417E-01, -2.522079E-02, +4.693460E-01
181+0.000000E+00, +0.000000E+00, +1.000000E+00, -7.061831E-01, -5.579556E-01, +1.891168E-01
182+0.000000E+00, +0.000000E+00, +0.000000E+00, +1.000000E+00, +6.764400E-01, -2.016856E-01
183+0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +1.000000E+00, -7.950259E-01
184+0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +1.000000E+00
185call setCov(cov(1:ndim, 1:ndim), lowDia, cov(1:ndim, 2:ndim+1), uppDia, std) ! convert upper-triangle correlation matrix to lower-triangle covariance matrix.
186cov
187+1.000000E+02, +1.000000E+00, +1.674978E-01, +3.232417E-01, -2.522079E-02, +4.693460E-01
188+1.172485E+01, +4.900000E+01, +1.000000E+00, -7.061831E-01, -5.579556E-01, +1.891168E-01
189+3.232417E+01, -4.943282E+01, +1.000000E+02, +1.000000E+00, +6.764400E-01, -2.016856E-01
190-2.269871E+00, -3.515120E+01, +6.087960E+01, +8.100000E+01, +1.000000E+00, -7.950259E-01
191+4.693460E+00, +1.323818E+00, -2.016856E+00, -7.155233E+00, +1.000000E+00, +1.000000E+00
192
193
194ndim = getUnifRand(1, 7)
195ndim
196+7
197std = getUnifRand(1, 10, ndim)
198std
199+1.000000E+00, +2.000000E+00, +4.000000E+00, +7.000000E+00, +1.000000E+01, +5.000000E+00, +6.000000E+00
200cor = getCovRand(1., ndim)
201cor
202+1.000000E+00, -8.015262E-01, +3.629454E-01, -2.950583E-01, +5.747038E-01, -5.639211E-01, -2.825486E-02
203-8.015262E-01, +1.000000E+00, -3.324285E-01, +4.862410E-01, -7.673784E-01, +4.132438E-01, -9.036963E-02
204+3.629454E-01, -3.324285E-01, +1.000000E+00, -2.069316E-01, +7.475872E-01, -7.230764E-01, +5.369414E-01
205-2.950583E-01, +4.862410E-01, -2.069316E-01, +1.000000E+00, -1.723250E-01, +5.715978E-01, -2.411085E-01
206+5.747038E-01, -7.673784E-01, +7.475872E-01, -1.723250E-01, +1.000000E+00, -4.669105E-01, +2.909592E-01
207-5.639211E-01, +4.132438E-01, -7.230764E-01, +5.715978E-01, -4.669105E-01, +1.000000E+00, -1.626980E-01
208-2.825486E-02, -9.036963E-02, +5.369414E-01, -2.411085E-01, +2.909592E-01, -1.626980E-01, +1.000000E+00
209cov = getFilled(0._TKG, ndim, ndim + 1)
210call setMatCopy(cov(1:ndim, 1:ndim), rdpack, cor, rdpack, lowDia)
211cov
212+1.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00
213-8.015262E-01, +1.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00
214+3.629454E-01, -3.324285E-01, +1.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00
215-2.950583E-01, +4.862410E-01, -2.069316E-01, +1.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00
216+5.747038E-01, -7.673784E-01, +7.475872E-01, -1.723250E-01, +1.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00
217-5.639211E-01, +4.132438E-01, -7.230764E-01, +5.715978E-01, -4.669105E-01, +1.000000E+00, +0.000000E+00, +0.000000E+00
218-2.825486E-02, -9.036963E-02, +5.369414E-01, -2.411085E-01, +2.909592E-01, -1.626980E-01, +1.000000E+00, +0.000000E+00
219call setCov(cov(1:ndim, 2:ndim+1), uppDia, cov(1:ndim, 1:ndim), lowDia, std) ! convert lower-triangle correlation matrix to upper-triangle covariance matrix.
220cov
221+1.000000E+00, +1.000000E+00, -1.603052E+00, +1.451782E+00, -2.065408E+00, +5.747038E+00, -2.819605E+00, -1.695291E-01
222-8.015262E-01, +1.000000E+00, +4.000000E+00, -2.659428E+00, +6.807374E+00, -1.534757E+01, +4.132438E+00, -1.084435E+00
223+3.629454E-01, -3.324285E-01, +1.000000E+00, +1.600000E+01, -5.794084E+00, +2.990349E+01, -1.446153E+01, +1.288659E+01
224-2.950583E-01, +4.862410E-01, -2.069316E-01, +1.000000E+00, +4.900000E+01, -1.206275E+01, +2.000592E+01, -1.012656E+01
225+5.747038E-01, -7.673784E-01, +7.475872E-01, -1.723250E-01, +1.000000E+00, +1.000000E+02, -2.334552E+01, +1.745755E+01
226-5.639211E-01, +4.132438E-01, -7.230764E-01, +5.715978E-01, -4.669105E-01, +1.000000E+00, +2.500000E+01, -4.880940E+00
227-2.825486E-02, -9.036963E-02, +5.369414E-01, -2.411085E-01, +2.909592E-01, -1.626980E-01, +1.000000E+00, +3.600000E+01
228
229cov = getFilled(0._TKG, ndim, ndim + 1)
230call setMatCopy(cov(1:ndim, 2:ndim+1), rdpack, cor, rdpack, uppDia)
231cov
232+0.000000E+00, +1.000000E+00, -8.015262E-01, +3.629454E-01, -2.950583E-01, +5.747038E-01, -5.639211E-01, -2.825486E-02
233+0.000000E+00, +0.000000E+00, +1.000000E+00, -3.324285E-01, +4.862410E-01, -7.673784E-01, +4.132438E-01, -9.036963E-02
234+0.000000E+00, +0.000000E+00, +0.000000E+00, +1.000000E+00, -2.069316E-01, +7.475872E-01, -7.230764E-01, +5.369414E-01
235+0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +1.000000E+00, -1.723250E-01, +5.715978E-01, -2.411085E-01
236+0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +1.000000E+00, -4.669105E-01, +2.909592E-01
237+0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +1.000000E+00, -1.626980E-01
238+0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +1.000000E+00
239call setCov(cov(1:ndim, 1:ndim), lowDia, cov(1:ndim, 2:ndim+1), uppDia, std) ! convert upper-triangle correlation matrix to lower-triangle covariance matrix.
240cov
241+1.000000E+00, +1.000000E+00, -8.015262E-01, +3.629454E-01, -2.950583E-01, +5.747038E-01, -5.639211E-01, -2.825486E-02
242-1.603052E+00, +4.000000E+00, +1.000000E+00, -3.324285E-01, +4.862410E-01, -7.673784E-01, +4.132438E-01, -9.036963E-02
243+1.451782E+00, -2.659428E+00, +1.600000E+01, +1.000000E+00, -2.069316E-01, +7.475872E-01, -7.230764E-01, +5.369414E-01
244-2.065408E+00, +6.807374E+00, -5.794084E+00, +4.900000E+01, +1.000000E+00, -1.723250E-01, +5.715978E-01, -2.411085E-01
245+5.747038E+00, -1.534757E+01, +2.990349E+01, -1.206275E+01, +1.000000E+02, +1.000000E+00, -4.669105E-01, +2.909592E-01
246-2.819605E+00, +4.132438E+00, -1.446153E+01, +2.000592E+01, -2.334552E+01, +2.500000E+01, +1.000000E+00, -1.626980E-01
247-1.695291E-01, -1.084435E+00, +1.288659E+01, -1.012656E+01, +1.745755E+01, -4.880940E+00, +3.600000E+01, +1.000000E+00
248
249
250ndim = getUnifRand(1, 7)
251ndim
252+6
253std = getUnifRand(1, 10, ndim)
254std
255+9.000000E+00, +5.000000E+00, +6.000000E+00, +3.000000E+00, +5.000000E+00, +6.000000E+00
256cor = getCovRand(1., ndim)
257cor
258+1.000000E+00, -8.654518E-01, -7.282782E-01, +9.898201E-02, +1.921818E-01, +5.303325E-01
259-8.654518E-01, +1.000000E+00, +6.867865E-01, -5.515494E-01, -1.489851E-02, -2.569512E-01
260-7.282782E-01, +6.867865E-01, +1.000000E+00, -3.788385E-01, -1.680465E-01, -6.524359E-01
261+9.898201E-02, -5.515494E-01, -3.788385E-01, +1.000000E+00, -3.473597E-01, -2.055989E-01
262+1.921818E-01, -1.489851E-02, -1.680465E-01, -3.473597E-01, +1.000000E+00, +7.154284E-01
263+5.303325E-01, -2.569512E-01, -6.524359E-01, -2.055989E-01, +7.154284E-01, +1.000000E+00
264cov = getFilled(0._TKG, ndim, ndim + 1)
265call setMatCopy(cov(1:ndim, 1:ndim), rdpack, cor, rdpack, lowDia)
266cov
267+1.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00
268-8.654518E-01, +1.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00
269-7.282782E-01, +6.867865E-01, +1.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00
270+9.898201E-02, -5.515494E-01, -3.788385E-01, +1.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00
271+1.921818E-01, -1.489851E-02, -1.680465E-01, -3.473597E-01, +1.000000E+00, +0.000000E+00, +0.000000E+00
272+5.303325E-01, -2.569512E-01, -6.524359E-01, -2.055989E-01, +7.154284E-01, +1.000000E+00, +0.000000E+00
273call setCov(cov(1:ndim, 2:ndim+1), uppDia, cov(1:ndim, 1:ndim), lowDia, std) ! convert lower-triangle correlation matrix to upper-triangle covariance matrix.
274cov
275+1.000000E+00, +8.100000E+01, -3.894533E+01, -3.932702E+01, +2.672514E+00, +8.648179E+00, +2.863795E+01
276-8.654518E-01, +1.000000E+00, +2.500000E+01, +2.060359E+01, -8.273241E+00, -3.724627E-01, -7.708535E+00
277-7.282782E-01, +6.867865E-01, +1.000000E+00, +3.600000E+01, -6.819093E+00, -5.041395E+00, -2.348769E+01
278+9.898201E-02, -5.515494E-01, -3.788385E-01, +1.000000E+00, +9.000000E+00, -5.210396E+00, -3.700780E+00
279+1.921818E-01, -1.489851E-02, -1.680465E-01, -3.473597E-01, +1.000000E+00, +2.500000E+01, +2.146285E+01
280+5.303325E-01, -2.569512E-01, -6.524359E-01, -2.055989E-01, +7.154284E-01, +1.000000E+00, +3.600000E+01
281
282cov = getFilled(0._TKG, ndim, ndim + 1)
283call setMatCopy(cov(1:ndim, 2:ndim+1), rdpack, cor, rdpack, uppDia)
284cov
285+0.000000E+00, +1.000000E+00, -8.654518E-01, -7.282782E-01, +9.898201E-02, +1.921818E-01, +5.303325E-01
286+0.000000E+00, +0.000000E+00, +1.000000E+00, +6.867865E-01, -5.515494E-01, -1.489851E-02, -2.569512E-01
287+0.000000E+00, +0.000000E+00, +0.000000E+00, +1.000000E+00, -3.788385E-01, -1.680465E-01, -6.524359E-01
288+0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +1.000000E+00, -3.473597E-01, -2.055989E-01
289+0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +1.000000E+00, +7.154284E-01
290+0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +1.000000E+00
291call setCov(cov(1:ndim, 1:ndim), lowDia, cov(1:ndim, 2:ndim+1), uppDia, std) ! convert upper-triangle correlation matrix to lower-triangle covariance matrix.
292cov
293+8.100000E+01, +1.000000E+00, -8.654518E-01, -7.282782E-01, +9.898201E-02, +1.921818E-01, +5.303325E-01
294-3.894533E+01, +2.500000E+01, +1.000000E+00, +6.867865E-01, -5.515494E-01, -1.489851E-02, -2.569512E-01
295-3.932702E+01, +2.060359E+01, +3.600000E+01, +1.000000E+00, -3.788385E-01, -1.680465E-01, -6.524359E-01
296+2.672514E+00, -8.273241E+00, -6.819093E+00, +9.000000E+00, +1.000000E+00, -3.473597E-01, -2.055989E-01
297+8.648179E+00, -3.724627E-01, -5.041395E+00, -5.210396E+00, +2.500000E+01, +1.000000E+00, +7.154284E-01
298+2.863795E+01, -7.708535E+00, -2.348769E+01, -3.700780E+00, +2.146285E+01, +3.600000E+01, +1.000000E+00
299
300
301ndim = getUnifRand(1, 7)
302ndim
303+4
304std = getUnifRand(1, 10, ndim)
305std
306+9.000000E+00, +2.000000E+00, +1.000000E+00, +2.000000E+00
307cor = getCovRand(1., ndim)
308cor
309+1.000000E+00, -9.067490E-01, +2.879211E-01, -2.249667E-01
310-9.067490E-01, +1.000000E+00, -6.610098E-01, +8.290101E-02
311+2.879211E-01, -6.610098E-01, +1.000000E+00, +3.001972E-01
312-2.249667E-01, +8.290101E-02, +3.001972E-01, +1.000000E+00
313cov = getFilled(0._TKG, ndim, ndim + 1)
314call setMatCopy(cov(1:ndim, 1:ndim), rdpack, cor, rdpack, lowDia)
315cov
316+1.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00
317-9.067490E-01, +1.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00
318+2.879211E-01, -6.610098E-01, +1.000000E+00, +0.000000E+00, +0.000000E+00
319-2.249667E-01, +8.290101E-02, +3.001972E-01, +1.000000E+00, +0.000000E+00
320call setCov(cov(1:ndim, 2:ndim+1), uppDia, cov(1:ndim, 1:ndim), lowDia, std) ! convert lower-triangle correlation matrix to upper-triangle covariance matrix.
321cov
322+1.000000E+00, +8.100000E+01, -1.632148E+01, +2.591290E+00, -4.049401E+00
323-9.067490E-01, +1.000000E+00, +4.000000E+00, -1.322020E+00, +3.316040E-01
324+2.879211E-01, -6.610098E-01, +1.000000E+00, +1.000000E+00, +6.003944E-01
325-2.249667E-01, +8.290101E-02, +3.001972E-01, +1.000000E+00, +4.000000E+00
326
327cov = getFilled(0._TKG, ndim, ndim + 1)
328call setMatCopy(cov(1:ndim, 2:ndim+1), rdpack, cor, rdpack, uppDia)
329cov
330+0.000000E+00, +1.000000E+00, -9.067490E-01, +2.879211E-01, -2.249667E-01
331+0.000000E+00, +0.000000E+00, +1.000000E+00, -6.610098E-01, +8.290101E-02
332+0.000000E+00, +0.000000E+00, +0.000000E+00, +1.000000E+00, +3.001972E-01
333+0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +1.000000E+00
334call setCov(cov(1:ndim, 1:ndim), lowDia, cov(1:ndim, 2:ndim+1), uppDia, std) ! convert upper-triangle correlation matrix to lower-triangle covariance matrix.
335cov
336+8.100000E+01, +1.000000E+00, -9.067490E-01, +2.879211E-01, -2.249667E-01
337-1.632148E+01, +4.000000E+00, +1.000000E+00, -6.610098E-01, +8.290101E-02
338+2.591290E+00, -1.322020E+00, +1.000000E+00, +1.000000E+00, +3.001972E-01
339-4.049401E+00, +3.316040E-01, +6.003944E-01, +4.000000E+00, +1.000000E+00
340
341
342ndim = getUnifRand(1, 7)
343ndim
344+7
345std = getUnifRand(1, 10, ndim)
346std
347+1.000000E+00, +7.000000E+00, +4.000000E+00, +6.000000E+00, +9.000000E+00, +4.000000E+00, +8.000000E+00
348cor = getCovRand(1., ndim)
349cor
350+1.000000E+00, +8.171124E-02, +1.833007E-01, -5.452252E-01, +4.400409E-01, +4.064409E-01, +3.405425E-02
351+8.171124E-02, +1.000000E+00, +3.496451E-01, +3.681374E-01, -1.257106E-01, +6.347072E-01, -2.514751E-01
352+1.833007E-01, +3.496451E-01, +1.000000E+00, +4.623571E-01, +4.696944E-01, +6.954970E-01, -1.163978E-01
353-5.452252E-01, +3.681374E-01, +4.623571E-01, +1.000000E+00, -3.147775E-01, +3.116178E-01, -4.606161E-01
354+4.400409E-01, -1.257106E-01, +4.696944E-01, -3.147775E-01, +1.000000E+00, +3.610232E-01, -2.773589E-02
355+4.064409E-01, +6.347072E-01, +6.954970E-01, +3.116178E-01, +3.610232E-01, +1.000000E+00, -1.573637E-01
356+3.405425E-02, -2.514751E-01, -1.163978E-01, -4.606161E-01, -2.773589E-02, -1.573637E-01, +1.000000E+00
357cov = getFilled(0._TKG, ndim, ndim + 1)
358call setMatCopy(cov(1:ndim, 1:ndim), rdpack, cor, rdpack, lowDia)
359cov
360+1.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00
361+8.171124E-02, +1.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00
362+1.833007E-01, +3.496451E-01, +1.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00
363-5.452252E-01, +3.681374E-01, +4.623571E-01, +1.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00
364+4.400409E-01, -1.257106E-01, +4.696944E-01, -3.147775E-01, +1.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00
365+4.064409E-01, +6.347072E-01, +6.954970E-01, +3.116178E-01, +3.610232E-01, +1.000000E+00, +0.000000E+00, +0.000000E+00
366+3.405425E-02, -2.514751E-01, -1.163978E-01, -4.606161E-01, -2.773589E-02, -1.573637E-01, +1.000000E+00, +0.000000E+00
367call setCov(cov(1:ndim, 2:ndim+1), uppDia, cov(1:ndim, 1:ndim), lowDia, std) ! convert lower-triangle correlation matrix to upper-triangle covariance matrix.
368cov
369+1.000000E+00, +1.000000E+00, +5.719787E-01, +7.332030E-01, -3.271351E+00, +3.960368E+00, +1.625764E+00, +2.724340E-01
370+8.171124E-02, +1.000000E+00, +4.900000E+01, +9.790064E+00, +1.546177E+01, -7.919768E+00, +1.777180E+01, -1.408260E+01
371+1.833007E-01, +3.496451E-01, +1.000000E+00, +1.600000E+01, +1.109657E+01, +1.690900E+01, +1.112795E+01, -3.724731E+00
372-5.452252E-01, +3.681374E-01, +4.623571E-01, +1.000000E+00, +3.600000E+01, -1.699798E+01, +7.478827E+00, -2.210957E+01
373+4.400409E-01, -1.257106E-01, +4.696944E-01, -3.147775E-01, +1.000000E+00, +8.100000E+01, +1.299683E+01, -1.996984E+00
374+4.064409E-01, +6.347072E-01, +6.954970E-01, +3.116178E-01, +3.610232E-01, +1.000000E+00, +1.600000E+01, -5.035639E+00
375+3.405425E-02, -2.514751E-01, -1.163978E-01, -4.606161E-01, -2.773589E-02, -1.573637E-01, +1.000000E+00, +6.400000E+01
376
377cov = getFilled(0._TKG, ndim, ndim + 1)
378call setMatCopy(cov(1:ndim, 2:ndim+1), rdpack, cor, rdpack, uppDia)
379cov
380+0.000000E+00, +1.000000E+00, +8.171124E-02, +1.833007E-01, -5.452252E-01, +4.400409E-01, +4.064409E-01, +3.405425E-02
381+0.000000E+00, +0.000000E+00, +1.000000E+00, +3.496451E-01, +3.681374E-01, -1.257106E-01, +6.347072E-01, -2.514751E-01
382+0.000000E+00, +0.000000E+00, +0.000000E+00, +1.000000E+00, +4.623571E-01, +4.696944E-01, +6.954970E-01, -1.163978E-01
383+0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +1.000000E+00, -3.147775E-01, +3.116178E-01, -4.606161E-01
384+0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +1.000000E+00, +3.610232E-01, -2.773589E-02
385+0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +1.000000E+00, -1.573637E-01
386+0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +1.000000E+00
387call setCov(cov(1:ndim, 1:ndim), lowDia, cov(1:ndim, 2:ndim+1), uppDia, std) ! convert upper-triangle correlation matrix to lower-triangle covariance matrix.
388cov
389+1.000000E+00, +1.000000E+00, +8.171124E-02, +1.833007E-01, -5.452252E-01, +4.400409E-01, +4.064409E-01, +3.405425E-02
390+5.719787E-01, +4.900000E+01, +1.000000E+00, +3.496451E-01, +3.681374E-01, -1.257106E-01, +6.347072E-01, -2.514751E-01
391+7.332030E-01, +9.790064E+00, +1.600000E+01, +1.000000E+00, +4.623571E-01, +4.696944E-01, +6.954970E-01, -1.163978E-01
392-3.271351E+00, +1.546177E+01, +1.109657E+01, +3.600000E+01, +1.000000E+00, -3.147775E-01, +3.116178E-01, -4.606161E-01
393+3.960368E+00, -7.919768E+00, +1.690900E+01, -1.699798E+01, +8.100000E+01, +1.000000E+00, +3.610232E-01, -2.773589E-02
394+1.625764E+00, +1.777180E+01, +1.112795E+01, +7.478827E+00, +1.299683E+01, +1.600000E+01, +1.000000E+00, -1.573637E-01
395+2.724340E-01, -1.408260E+01, -3.724731E+00, -2.210957E+01, -1.996984E+00, -5.035639E+00, +6.400000E+01, +1.000000E+00
396
397
398ndim = getUnifRand(1, 7)
399ndim
400+7
401std = getUnifRand(1, 10, ndim)
402std
403+3.000000E+00, +1.000000E+01, +9.000000E+00, +8.000000E+00, +4.000000E+00, +4.000000E+00, +4.000000E+00
404cor = getCovRand(1., ndim)
405cor
406+1.000000E+00, +2.009096E-01, -6.263593E-01, +3.890767E-01, +3.159687E-01, +3.334182E-02, -5.066900E-01
407+2.009096E-01, +1.000000E+00, -8.316587E-01, -6.206609E-01, -4.078837E-01, +4.053683E-01, -7.023707E-01
408-6.263593E-01, -8.316587E-01, +1.000000E+00, +3.916197E-01, +2.127269E-01, -1.632651E-01, +7.485247E-01
409+3.890767E-01, -6.206609E-01, +3.916197E-01, +1.000000E+00, +7.647105E-01, -1.790522E-01, +3.488088E-01
410+3.159687E-01, -4.078837E-01, +2.127269E-01, +7.647105E-01, +1.000000E+00, -3.859851E-01, +3.717218E-01
411+3.334182E-02, +4.053683E-01, -1.632651E-01, -1.790522E-01, -3.859851E-01, +1.000000E+00, -6.094661E-01
412-5.066900E-01, -7.023707E-01, +7.485247E-01, +3.488088E-01, +3.717218E-01, -6.094661E-01, +1.000000E+00
413cov = getFilled(0._TKG, ndim, ndim + 1)
414call setMatCopy(cov(1:ndim, 1:ndim), rdpack, cor, rdpack, lowDia)
415cov
416+1.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00
417+2.009096E-01, +1.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00
418-6.263593E-01, -8.316587E-01, +1.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00
419+3.890767E-01, -6.206609E-01, +3.916197E-01, +1.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00
420+3.159687E-01, -4.078837E-01, +2.127269E-01, +7.647105E-01, +1.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00
421+3.334182E-02, +4.053683E-01, -1.632651E-01, -1.790522E-01, -3.859851E-01, +1.000000E+00, +0.000000E+00, +0.000000E+00
422-5.066900E-01, -7.023707E-01, +7.485247E-01, +3.488088E-01, +3.717218E-01, -6.094661E-01, +1.000000E+00, +0.000000E+00
423call setCov(cov(1:ndim, 2:ndim+1), uppDia, cov(1:ndim, 1:ndim), lowDia, std) ! convert lower-triangle correlation matrix to upper-triangle covariance matrix.
424cov
425+1.000000E+00, +9.000000E+00, +6.027289E+00, -1.691170E+01, +9.337840E+00, +3.791625E+00, +4.001018E-01, -6.080279E+00
426+2.009096E-01, +1.000000E+00, +1.000000E+02, -7.484928E+01, -4.965287E+01, -1.631535E+01, +1.621473E+01, -2.809483E+01
427-6.263593E-01, -8.316587E-01, +1.000000E+00, +8.100000E+01, +2.819662E+01, +7.658167E+00, -5.877543E+00, +2.694689E+01
428+3.890767E-01, -6.206609E-01, +3.916197E-01, +1.000000E+00, +6.400000E+01, +2.447074E+01, -5.729670E+00, +1.116188E+01
429+3.159687E-01, -4.078837E-01, +2.127269E-01, +7.647105E-01, +1.000000E+00, +1.600000E+01, -6.175761E+00, +5.947549E+00
430+3.334182E-02, +4.053683E-01, -1.632651E-01, -1.790522E-01, -3.859851E-01, +1.000000E+00, +1.600000E+01, -9.751458E+00
431-5.066900E-01, -7.023707E-01, +7.485247E-01, +3.488088E-01, +3.717218E-01, -6.094661E-01, +1.000000E+00, +1.600000E+01
432
433cov = getFilled(0._TKG, ndim, ndim + 1)
434call setMatCopy(cov(1:ndim, 2:ndim+1), rdpack, cor, rdpack, uppDia)
435cov
436+0.000000E+00, +1.000000E+00, +2.009096E-01, -6.263593E-01, +3.890767E-01, +3.159687E-01, +3.334182E-02, -5.066900E-01
437+0.000000E+00, +0.000000E+00, +1.000000E+00, -8.316587E-01, -6.206609E-01, -4.078837E-01, +4.053683E-01, -7.023707E-01
438+0.000000E+00, +0.000000E+00, +0.000000E+00, +1.000000E+00, +3.916197E-01, +2.127269E-01, -1.632651E-01, +7.485247E-01
439+0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +1.000000E+00, +7.647105E-01, -1.790522E-01, +3.488088E-01
440+0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +1.000000E+00, -3.859851E-01, +3.717218E-01
441+0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +1.000000E+00, -6.094661E-01
442+0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +1.000000E+00
443call setCov(cov(1:ndim, 1:ndim), lowDia, cov(1:ndim, 2:ndim+1), uppDia, std) ! convert upper-triangle correlation matrix to lower-triangle covariance matrix.
444cov
445+9.000000E+00, +1.000000E+00, +2.009096E-01, -6.263593E-01, +3.890767E-01, +3.159687E-01, +3.334182E-02, -5.066900E-01
446+6.027289E+00, +1.000000E+02, +1.000000E+00, -8.316587E-01, -6.206609E-01, -4.078837E-01, +4.053683E-01, -7.023707E-01
447-1.691170E+01, -7.484928E+01, +8.100000E+01, +1.000000E+00, +3.916197E-01, +2.127269E-01, -1.632651E-01, +7.485247E-01
448+9.337840E+00, -4.965287E+01, +2.819662E+01, +6.400000E+01, +1.000000E+00, +7.647105E-01, -1.790522E-01, +3.488088E-01
449+3.791625E+00, -1.631535E+01, +7.658167E+00, +2.447074E+01, +1.600000E+01, +1.000000E+00, -3.859851E-01, +3.717218E-01
450+4.001018E-01, +1.621473E+01, -5.877543E+00, -5.729670E+00, -6.175761E+00, +1.600000E+01, +1.000000E+00, -6.094661E-01
451-6.080279E+00, -2.809483E+01, +2.694689E+01, +1.116188E+01, +5.947549E+00, -9.751458E+00, +1.600000E+01, +1.000000E+00
452
453
454!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
455!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
456!Compute the covariance matrix of a 2-D sample.
457!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
458!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
459
460ndim = 2; nsam = 10; dim = 2
461sample = reshape(getUnifRand(1, 20, ndim * nsam), shape = [ndim, nsam], order = [2, 1])
462sample
463+20.0000000, +7.00000000, +15.0000000, +9.00000000, +10.0000000, +7.00000000, +5.00000000, +7.00000000, +15.0000000, +20.0000000
464+16.0000000, +5.00000000, +10.0000000, +16.0000000, +6.00000000, +13.0000000, +18.0000000, +11.0000000, +11.0000000, +3.00000000
465mean = getMean(sample, dim)
466mean
467+11.5000000, +10.9000006
468samShifted = getShifted(sample, dim, -mean)
469samShifted
470+8.50000000, -4.50000000, +3.50000000, -2.50000000, -1.50000000, -4.50000000, -6.50000000, -4.50000000, +3.50000000, +8.50000000
471+5.09999943, -5.90000057, -0.900000572, +5.09999943, -4.90000057, +2.09999943, +7.09999943, +0.999994278E-1, +0.999994278E-1, -7.90000057
472cov = getFilled(0., ndim, ndim)
473call setCov(cov, uppDia, mean, sample, dim)
474cov
475+28.0500011, -6.15000010
476+0.00000000, +22.8899994
477cov = getFilled(0., ndim, ndim)
478call setCov(cov, uppDia, samShifted, dim) ! same result as above.
479cov
480+28.0500011, -6.15000010
481+0.00000000, +22.8899994
482
483cov = getFilled(0., ndim, ndim)
484call setCov(cov, lowDia, mean, sample, dim)
485cov
486+28.0500011, +0.00000000
487-6.15000010, +22.8899994
488cov = getFilled(0., ndim, ndim)
489call setCov(cov, lowDia, samShifted, dim) ! same result as above.
490cov
491+28.0500011, +0.00000000
492-6.15000010, +22.8899994
493
494'Compute the sample covariance along the first dimension.'
495
496dim = 1
497cov = getFilled(0., ndim, ndim)
498call setCov(cov, uppDia, mean, transpose(sample), dim)
499cov
500+28.0500011, -6.15000010
501+0.00000000, +22.8899994
502cov = getFilled(0., ndim, ndim)
503call setCov(cov, uppDia, transpose(samShifted), dim) ! same result as above.
504cov
505+28.0500011, -6.15000010
506+0.00000000, +22.8899994
507
508cov = getFilled(0., ndim, ndim)
509call setCov(cov, lowDia, mean, transpose(sample), dim)
510cov
511+28.0500011, +0.00000000
512-6.15000010, +22.8899994
513cov = getFilled(0., ndim, ndim)
514call setCov(cov, lowDia, transpose(samShifted), dim) ! same result as above.
515cov
516+28.0500011, +0.00000000
517-6.15000010, +22.8899994
518
519'Compute the full sample covariance for a pair of time series.'
520
521cov = getFilled(0., ndim, ndim)
522call setCov(cov, mean, sample(1,:), sample(2,:))
523cov
524+28.0500011, -6.15000010
525-6.15000010, +22.8899994
526cov = getFilled(0., ndim, ndim)
527call setCov(cov, samShifted(1,:), samShifted(2,:)) ! same result as above.
528cov
529+28.0500011, -6.15000010
530-6.15000010, +22.8899994
531
532
533!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
534!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
535!Compute the biased covariance matrix of a weighted 2-D sample.
536!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
537!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
538
539ndim = 2; nsam = 10; dim = 2
540sample = reshape(getUnifRand(1, 20, ndim * nsam), shape = [ndim, nsam], order = [2, 1])
541sample
542+7.00000000, +20.0000000, +18.0000000, +14.0000000, +20.0000000, +7.00000000, +10.0000000, +19.0000000, +7.00000000, +6.00000000
543+19.0000000, +1.00000000, +4.00000000, +15.0000000, +10.0000000, +8.00000000, +14.0000000, +7.00000000, +5.00000000, +3.00000000
544call setResized(mean, ndim)
545iweight = getUnifRand(1, 10, nsam) ! integer-valued weights.
546iweight
547+3, +8, +2, +4, +1, +7, +8, +8, +8, +3
548call setMean(mean, sample, dim, iweight, iweisum)
549mean
550+12.4615393, +8.00000000
551iweisum
552+52
553rweight = iweight ! or real-valued weights.
554iweight
555+3, +8, +2, +4, +1, +7, +8, +8, +8, +3
556call setMean(mean, sample, dim, rweight, rweisum)
557mean
558+12.4615393, +8.00000000
559rweisum
560+52.0000000
561
562!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
563!Compute the covariance matrix integer weights.
564!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
565
566samShifted = getShifted(sample, dim, -mean)
567samShifted
568-5.46153927, +7.53846073, +5.53846073, +1.53846073, +7.53846073, -5.46153927, -2.46153927, +6.53846073, -5.46153927, -6.46153927
569+11.0000000, -7.00000000, -4.00000000, +7.00000000, +2.00000000, +0.00000000, +6.00000000, -1.00000000, -3.00000000, -5.00000000
570cov = getFilled(0., ndim, ndim)
571call setCov(cov, uppDia, mean, sample, dim, iweight, iweisum)
572cov
573+31.4408321, -10.2115393
574+0.00000000, +27.5000019
575cov = getFilled(0., ndim, ndim)
576call setCov(cov, uppDia, samShifted, dim, iweight, sum(iweight)) ! same result as above.
577cov
578+31.4408321, -10.2115393
579+0.00000000, +27.5000019
580
581cov = getFilled(0., ndim, ndim)
582call setCov(cov, lowDia, mean, sample, dim, iweight, iweisum)
583cov
584+31.4408321, +0.00000000
585-10.2115393, +27.5000019
586cov = getFilled(0., ndim, ndim)
587call setCov(cov, lowDia, samShifted, dim, iweight, sum(iweight)) ! same result as above.
588cov
589+31.4408321, +0.00000000
590-10.2115393, +27.5000019
591
592'Compute the sample covariance along the first dimension.'
593
594dim = 1
595cov = getFilled(0., ndim, ndim)
596call setCov(cov, uppDia, mean, transpose(sample), dim, iweight, iweisum)
597cov
598+31.4408321, -10.2115393
599+0.00000000, +27.5000019
600cov = getFilled(0., ndim, ndim)
601call setCov(cov, uppDia, transpose(samShifted), dim, iweight, sum(iweight)) ! same result as above.
602cov
603+31.4408321, -10.2115393
604+0.00000000, +27.5000019
605
606cov = getFilled(0., ndim, ndim)
607call setCov(cov, lowDia, mean, transpose(sample), dim, iweight, iweisum)
608cov
609+31.4408321, +0.00000000
610-10.2115393, +27.5000019
611cov = getFilled(0., ndim, ndim)
612call setCov(cov, lowDia, transpose(samShifted), dim, iweight, sum(iweight)) ! same result as above.
613cov
614+31.4408321, +0.00000000
615-10.2115393, +27.5000019
616
617'Compute the full sample covariance for a pair of time series.'
618
619cov = getFilled(0., ndim, ndim)
620call setCov(cov, mean, sample(1,:), sample(2,:), iweight, iweisum)
621cov
622+31.4408321, -10.2115393
623-10.2115393, +27.5000019
624cov = getFilled(0., ndim, ndim)
625call setCov(cov, samShifted(1,:), samShifted(2,:), iweight, sum(iweight)) ! same result as above.
626cov
627+31.4408321, -10.2115393
628-10.2115393, +27.5000019
629
630
631!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
632!Compute the covariance matrix real weights.
633!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
634
635dim = 2_IK
636samShifted = getShifted(sample, dim, -mean)
637samShifted
638-5.46153927, +7.53846073, +5.53846073, +1.53846073, +7.53846073, -5.46153927, -2.46153927, +6.53846073, -5.46153927, -6.46153927
639+11.0000000, -7.00000000, -4.00000000, +7.00000000, +2.00000000, +0.00000000, +6.00000000, -1.00000000, -3.00000000, -5.00000000
640cov = getFilled(0., ndim, ndim)
641call setCov(cov, uppDia, mean, sample, dim, rweight, rweisum)
642cov
643+31.4408321, -10.2115393
644+0.00000000, +27.5000019
645cov = getFilled(0., ndim, ndim)
646call setCov(cov, uppDia, samShifted, dim, rweight, sum(rweight)) ! same result as above.
647cov
648+31.4408321, -10.2115393
649+0.00000000, +27.5000019
650
651cov = getFilled(0., ndim, ndim)
652call setCov(cov, lowDia, mean, sample, dim, rweight, rweisum)
653cov
654+31.4408321, +0.00000000
655-10.2115393, +27.5000019
656cov = getFilled(0., ndim, ndim)
657call setCov(cov, lowDia, samShifted, dim, rweight, sum(rweight)) ! same result as above.
658cov
659+31.4408321, +0.00000000
660-10.2115393, +27.5000019
661
662'Compute the sample covariance along the first dimension.'
663
664dim = 1
665cov = getFilled(0., ndim, ndim)
666call setCov(cov, uppDia, mean, transpose(sample), dim, rweight, rweisum)
667cov
668+31.4408321, -10.2115393
669+0.00000000, +27.5000019
670cov = getFilled(0., ndim, ndim)
671call setCov(cov, uppDia, transpose(samShifted), dim, rweight, sum(rweight)) ! same result as above.
672cov
673+31.4408321, -10.2115393
674+0.00000000, +27.5000019
675
676cov = getFilled(0., ndim, ndim)
677call setCov(cov, lowDia, mean, transpose(sample), dim, rweight, rweisum)
678cov
679+31.4408321, +0.00000000
680-10.2115393, +27.5000019
681cov = getFilled(0., ndim, ndim)
682call setCov(cov, lowDia, transpose(samShifted), dim, rweight, sum(rweight)) ! same result as above.
683cov
684+31.4408321, +0.00000000
685-10.2115393, +27.5000019
686
687'Compute the full sample covariance for a pair of time series.'
688
689cov = getFilled(0., ndim, ndim)
690call setCov(cov, mean, sample(1,:), sample(2,:), rweight, rweisum)
691cov
692+31.4408321, -10.2115393
693-10.2115393, +27.5000019
694cov = getFilled(0., ndim, ndim)
695call setCov(cov, samShifted(1,:), samShifted(2,:), rweight, sum(rweight)) ! same result as above.
696cov
697+31.4408321, -10.2115393
698-10.2115393, +27.5000019
699
700
Benchmarks:


Benchmark :: The runtime performance of setCov along different sample dimensions.

1! Test the performance of Cholesky factorization computation using an assumed-shape interface vs. explicit-shape interface.
2program benchmark
3
4 use pm_kind, only: IK, LK, RKG => RKD, SK
5 use pm_sampleCov, only: uppDia
6 use pm_bench, only: bench_type
7
8 implicit none
9
10 integer(IK) :: itry, ntry
11 integer(IK) :: i
12 integer(IK) :: iarr
13 integer(IK) :: fileUnit
14 integer(IK) , parameter :: NARR = 18_IK
15 integer(IK) , allocatable :: rperm(:)
16 real(RKG) , allocatable :: samdim1(:,:)
17 real(RKG) , allocatable :: samdim2(:,:)
18 type(bench_type), allocatable :: bench(:)
19 integer(IK) , parameter :: nsammax = 2**NARR
20 integer(IK) , parameter :: ndim = 5_IK
21 real(RKG) :: mean(ndim), cov(ndim, ndim)
22 integer(IK) :: idim, jdim, isam, nsam
23 real(RKG) :: dumm
24
25 bench = [ bench_type(name = SK_"intrinsicDIM1", exec = intrinsicDIM1, overhead = setOverhead) &
26 , bench_type(name = SK_"intrinsicDIM2", exec = intrinsicDIM2, overhead = setOverhead) &
27 , bench_type(name = SK_"setCovDIM1", exec = setCovDIM1, overhead = setOverhead) &
28 , bench_type(name = SK_"setCovDIM2", exec = setCovDIM2, overhead = setOverhead) &
29 ]
30
31 write(*,"(*(g0,:,' '))")
32 write(*,"(*(g0,:,' '))") "sample covariance benchmarking..."
33 write(*,"(*(g0,:,' '))")
34
35 open(newunit = fileUnit, file = "main.out", status = "replace")
36
37 write(fileUnit, "(*(g0,:,','))") "nsam", (bench(i)%name, i = 1, size(bench))
38
39 dumm = 0._RKG
40 loopOverMatrixSize: do iarr = 1, NARR - 1
41
42 nsam = 2**iarr
43 ntry = nsammax / nsam
44 allocate(samdim1(nsam, ndim), samdim2(ndim, nsam))
45 write(*,"(*(g0,:,' '))") "Benchmarking setCovDIM1() vs. setCovDIM2()", nsam, ntry
46
47 do i = 1, size(bench)
48 bench(i)%timing = bench(i)%getTiming()
49 end do
50
51 write(fileUnit,"(*(g0,:,','))") nsam, (bench(i)%timing%mean / ntry, i = 1, size(bench))
52 deallocate(samdim1, samdim2)
53
54 end do loopOverMatrixSize
55 write(*,"(*(g0,:,' '))") dumm
56
57 close(fileUnit)
58
59contains
60
61 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
62 ! procedure wrappers.
63 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
64
65 subroutine setOverhead()
66 do itry = 1, ntry
67 call random_number(mean)
68 dumm = dumm - mean(1)
69 end do
70 end subroutine
71
72 subroutine setSamDIM1()
73 do isam = 1, nsam
74 call random_number(samdim1(isam, 1 : ndim))
75 end do
76 end subroutine
77
78 subroutine setSamDIM2()
79 do isam = 1, nsam
80 call random_number(samdim2(1 : ndim, isam))
81 end do
82 end subroutine
83
84 subroutine setCovDIM1()
85 block
86 use pm_sampleCov, only: setCov
87 do itry = 1, ntry
88 call setSamDIM1()
89 call setCov(cov, uppDia, samdim1, dim = 1_IK)
90 dumm = dumm + cov(1,1)
91 end do
92 end block
93 end subroutine
94
95 subroutine setCovDIM2()
96 block
97 use pm_sampleCov, only: setCov
98 do itry = 1, ntry
99 call setSamDIM2()
100 call setCov(cov, uppDia, samdim2, dim = 2_IK)
101 dumm = dumm + cov(1,1)
102 end do
103 end block
104 end subroutine
105
106 subroutine intrinsicDIM1()
107 do itry = 1, ntry
108 call setSamDIM1()
109 do jdim = 1, ndim
110 do idim = 1, jdim
111 cov(idim, jdim) = dot_product(samdim1(1 : nsam, idim), samdim1(1 : nsam, jdim)) / size(samdim1, dim = 1)
112 end do
113 end do
114 dumm = dumm + cov(1,1)
115 end do
116 end subroutine
117
118 subroutine intrinsicDIM2()
119 do itry = 1, ntry
120 call setSamDIM2()
121 do jdim = 1, ndim
122 do idim = 1, jdim
123 cov(idim, jdim) = dot_product(samdim2(idim, 1 : nsam), samdim2(jdim, 1 : nsam)) / size(samdim2, dim = 2)
124 end do
125 end do
126 dumm = dumm + cov(1,1)
127 end do
128 end subroutine
129
130end program benchmark
Generate and return an object of type timing_type containing the benchmark timing information and sta...
Definition: pm_bench.F90:574
This module contains abstract interfaces and types that facilitate benchmarking of different procedur...
Definition: pm_bench.F90:41
integer, parameter RKD
The double precision real kind in Fortran mode. On most platforms, this is an 64-bit real kind.
Definition: pm_kind.F90:568
This is the class for creating benchmark and performance-profiling objects.
Definition: pm_bench.F90:386

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

Postprocessing of the benchmark output
1#!/usr/bin/env python
2
3import matplotlib.pyplot as plt
4import pandas as pd
5import numpy as np
6
7import os
8dirname = os.path.basename(os.getcwd())
9
10fontsize = 14
11
12df = pd.read_csv("main.out", delimiter = ",")
13colnames = list(df.columns.values)
14
15
18
19ax = plt.figure(figsize = 1.25 * np.array([6.4,4.6]), dpi = 200)
20ax = plt.subplot()
21
22for colname in colnames[1:]:
23 plt.plot( df[colnames[0]].values
24 , df[colname].values
25 , linewidth = 2
26 )
27
28plt.xticks(fontsize = fontsize)
29plt.yticks(fontsize = fontsize)
30ax.set_xlabel(colnames[0], fontsize = fontsize)
31ax.set_ylabel("Runtime [ seconds ]", fontsize = fontsize)
32ax.set_title(" vs. ".join(colnames[1:])+"\nLower is better.", fontsize = fontsize)
33ax.set_xscale("log")
34ax.set_yscale("log")
35plt.minorticks_on()
36plt.grid(visible = True, which = "both", axis = "both", color = "0.85", linestyle = "-")
37ax.tick_params(axis = "y", which = "minor")
38ax.tick_params(axis = "x", which = "minor")
39ax.legend ( colnames[1:]
40 #, loc='center left'
41 #, bbox_to_anchor=(1, 0.5)
42 , fontsize = fontsize
43 )
44
45plt.tight_layout()
46plt.savefig("benchmark." + dirname + ".runtime.png")
47
48
51
52ax = plt.figure(figsize = 1.25 * np.array([6.4,4.6]), dpi = 200)
53ax = plt.subplot()
54
55plt.plot( df[colnames[0]].values
56 , np.ones(len(df[colnames[0]].values))
57 , linestyle = "--"
58 #, color = "black"
59 , linewidth = 2
60 )
61for colname in colnames[2:]:
62 plt.plot( df[colnames[0]].values
63 , df[colname].values / df[colnames[1]].values
64 , linewidth = 2
65 )
66
67plt.xticks(fontsize = fontsize)
68plt.yticks(fontsize = fontsize)
69ax.set_xlabel(colnames[0], fontsize = fontsize)
70ax.set_ylabel("Runtime compared to {}".format(colnames[1]), fontsize = fontsize)
71ax.set_title("Runtime Ratio Comparison. Lower means faster.\nLower than 1 means faster than {}().".format(colnames[1]), fontsize = fontsize)
72ax.set_xscale("log")
73ax.set_yscale("log")
74plt.minorticks_on()
75plt.grid(visible = True, which = "both", axis = "both", color = "0.85", linestyle = "-")
76ax.tick_params(axis = "y", which = "minor")
77ax.tick_params(axis = "x", which = "minor")
78ax.legend ( colnames[1:]
79 #, bbox_to_anchor = (1, 0.5)
80 #, loc = "center left"
81 , fontsize = fontsize
82 )
83
84plt.tight_layout()
85plt.savefig("benchmark." + dirname + ".runtime.ratio.png")

Visualization of the benchmark output

Benchmark moral
  1. The procedures under the generic interface setCov can compute the covariance under two different sample axes.
  2. Recall that C is a row-major language while Fortran is a column-major language.
  3. As such, one would expect the computations for a sample whose observations are stored along the second axis would be faster in the Fortran programming language.
  4. However, such an expectation does not appear to hold at all times and appears to depend significantly on the computing architecture and the number of data attributes involved.
  5. The higher the number of data attributes, the more likely the computations along the second axis of sample will be faster.
  6. Note that for small number of data attributes, the computations along the second data axis involve a small loop that has significant computational cost due to the implicit branching involved in the loop.
Test:
test_pm_sampleCov
Todo:
Normal Priority: The examples of this generic interface should be extended to corrected weighted covariance matrices.


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.
Fatemeh Bagheri, Monday 02:15 AM, September 27, 2021, Dallas, TX

Definition at line 1694 of file pm_sampleCov.F90.


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