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:

compute the covariance matrix corresponding to an input correlation matrix and vector of standard deviations in arbitrary ndim dimensions.
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, type `complex` of kind any supported by the processor (e.g., CK, CK32, CK64, or CK128), 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, 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} 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: The constant lowDia, implying that only the lower-diagonal subset of the output covariance matrix must be computed. 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: The constant lowDia, implying that only the lower-diagonal subset of the input correlation matrix must be used. 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. If the input `sample` is missing and `x` and `y` are present, then `mean` must be a vector of size `ndim = 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. If `dim = 1`, the input `sample` is assumed to have the shape `(nsam, ndim)`. 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, type `integer` of default kind IK, or 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)

pm_sampleCov::setCov
Return the covariance matrix corresponding to the input (potentially weighted) correlation matrix or ...
Definition: pm_sampleCov.F90:1694

pm_sampleCov
This module contains classes and procedures for computing the properties related to the covariance ma...
Definition: pm_sampleCov.F90:170

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

6 use pm_sampleShift, only: getShifted

7 use pm_arraySpace, only: getLinSpace

8 use pm_arrayResize, only: setResized

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

pm_arrayFill::getFilled
Generate and return an array of the specified rank and shape of arbitrary intrinsic type and kind wit...
Definition: pm_arrayFill.F90:132

pm_arrayResize::setResized
Allocate or resize (shrink or expand) an input allocatable scalar string or array of rank 1....
Definition: pm_arrayResize.F90:249

pm_arraySpace::getLinSpace
Generate count evenly spaced points over the interval [x1, x2] if x1 < x2, or [x2,...
Definition: pm_arraySpace.F90:108

pm_distCov::getCovRand
Generate and return a random positive-definite (correlation or covariance) matrix using the Gram meth...
Definition: pm_distCov.F90:394

pm_distUnif::getUnifRand
Generate and return a scalar or a contiguous array of rank 1 of length s1 of randomly uniformly distr...
Definition: pm_distUnif.F90:4150

pm_io::getFormat
Generate and return a generic or type/kind-specific IO format with the requested specifications that ...
Definition: pm_io.F90:18485

pm_io::show
This is a generic method of the derived type display_type with pass attribute.
Definition: pm_io.F90:11726

pm_io::skip
This is a generic method of the derived type display_type with pass attribute.
Definition: pm_io.F90:11508

pm_matrixCopy::setMatCopy
Copy a desired subset of the input source matrix of arbitrary shape (:) or (:,:) to the target subset...
Definition: pm_matrixCopy.F90:17748

pm_sampleMean::getMean
Generate and return the (weighted) mean of an input sample of nsam observations with ndim = 1 or 2 at...
Definition: pm_sampleMean.F90:229

pm_sampleMean::setMean
Return the (weighted) mean of a pair of time series or of an input sample of nsam observations with n...
Definition: pm_sampleMean.F90:2008

pm_sampleShift::getShifted
Generate a sample of shape (nsam), or (ndim, nsam) or (nsam, ndim) that is shifted by the specified i...
Definition: pm_sampleShift.F90:197

pm_arrayFill
This module contains procedures and generic interfaces for convenient allocation and filling of array...
Definition: pm_arrayFill.F90:42

pm_arrayResize
This module contains procedures and generic interfaces for resizing allocatable arrays of various typ...
Definition: pm_arrayResize.F90:81

pm_arraySpace
This module contains procedures and generic interfaces for generating arrays with linear or logarithm...
Definition: pm_arraySpace.F90:33

pm_distCov
This module contains classes and procedures for generating random matrices distributed on the space o...
Definition: pm_distCov.F90:72

pm_distUnif
This module contains classes and procedures for computing various statistical quantities related to t...
Definition: pm_distUnif.F90:274

pm_io
This module contains classes and procedures for input/output (IO) or generic display operations on st...
Definition: pm_io.F90:252

pm_io::disp
type(display_type) disp
This is a scalar module variable an object of type display_type for general display.
Definition: pm_io.F90:11393

pm_kind
This module defines the relevant Fortran kind type-parameters frequently used in the ParaMonte librar...
Definition: pm_kind.F90:268

pm_kind::LK
integer, parameter LK
The default logical kind in the ParaMonte library: kind(.true.) in Fortran, kind(....
Definition: pm_kind.F90:541

pm_kind::IK
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

pm_kind::SK
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

pm_kind::RKS
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

pm_matrixCopy
This module contains procedures and generic interfaces relevant to copying (diagonal or upper/lower t...
Definition: pm_matrixCopy.F90:113

pm_sampleMean
This module contains classes and procedures for computing the first moment (i.e., the statistical mea...
Definition: pm_sampleMean.F90:104

pm_sampleShift
This module contains classes and procedures for shifting univariate or multivariate samples by arbitr...
Definition: pm_sampleShift.F90:59

pm_io::display_type
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+4

12std = getUnifRand(1, 10, ndim)

13std

14+7.000000E+00, +9.000000E+00, +6.000000E+00, +4.000000E+00

15cor = getCovRand(1., ndim)

16cor

17+1.000000E+00, +6.544206E-01, +5.815119E-01, +1.371571E-01

18+6.544206E-01, +1.000000E+00, +1.711279E-02, -4.053138E-01

19+5.815119E-01, +1.711279E-02, +1.000000E+00, +6.762713E-02

20+1.371571E-01, -4.053138E-01, +6.762713E-02, +1.000000E+00

21cov = getFilled(0._TKG, ndim, ndim + 1)

22call setMatCopy(cov(1:ndim, 1:ndim), rdpack, cor, rdpack, lowDia)

23cov

24+1.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00

25+6.544206E-01, +1.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00

26+5.815119E-01, +1.711279E-02, +1.000000E+00, +0.000000E+00, +0.000000E+00

27+1.371571E-01, -4.053138E-01, +6.762713E-02, +1.000000E+00, +0.000000E+00

28call 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.

29cov

30+1.000000E+00, +4.900000E+01, +4.122850E+01, +2.442350E+01, +3.840398E+00

31+6.544206E-01, +1.000000E+00, +8.100000E+01, +9.240907E-01, -1.459130E+01

32+5.815119E-01, +1.711279E-02, +1.000000E+00, +3.600000E+01, +1.623051E+00

33+1.371571E-01, -4.053138E-01, +6.762713E-02, +1.000000E+00, +1.600000E+01

34

35cov = getFilled(0._TKG, ndim, ndim + 1)

36call setMatCopy(cov(1:ndim, 2:ndim+1), rdpack, cor, rdpack, uppDia)

37cov

38+0.000000E+00, +1.000000E+00, +6.544206E-01, +5.815119E-01, +1.371571E-01

39+0.000000E+00, +0.000000E+00, +1.000000E+00, +1.711279E-02, -4.053138E-01

40+0.000000E+00, +0.000000E+00, +0.000000E+00, +1.000000E+00, +6.762713E-02

41+0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +1.000000E+00

42call 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.

43cov

44+4.900000E+01, +1.000000E+00, +6.544206E-01, +5.815119E-01, +1.371571E-01

45+4.122850E+01, +8.100000E+01, +1.000000E+00, +1.711279E-02, -4.053138E-01

46+2.442350E+01, +9.240907E-01, +3.600000E+01, +1.000000E+00, +6.762713E-02

47+3.840398E+00, -1.459130E+01, +1.623051E+00, +1.600000E+01, +1.000000E+00

48

49

50ndim = getUnifRand(1, 7)

51ndim

52+6

53std = getUnifRand(1, 10, ndim)

54std

55+2.000000E+00, +3.000000E+00, +9.000000E+00, +4.000000E+00, +2.000000E+00, +2.000000E+00

56cor = getCovRand(1., ndim)

57cor

58+1.000000E+00, +4.023281E-02, +7.728843E-01, +4.914861E-01, +9.849682E-02, -5.807413E-01

59+4.023281E-02, +1.000000E+00, +6.530263E-01, +6.725434E-01, -7.959484E-01, +5.932847E-01

60+7.728843E-01, +6.530263E-01, +1.000000E+00, +8.497998E-01, -4.249785E-01, -7.603369E-02

61+4.914861E-01, +6.725434E-01, +8.497998E-01, +1.000000E+00, -6.386417E-01, +1.259757E-01

62+9.849682E-02, -7.959484E-01, -4.249785E-01, -6.386417E-01, +1.000000E+00, -6.656111E-01

63-5.807413E-01, +5.932847E-01, -7.603369E-02, +1.259757E-01, -6.656111E-01, +1.000000E+00

64cov = getFilled(0._TKG, ndim, ndim + 1)

65call setMatCopy(cov(1:ndim, 1:ndim), rdpack, cor, rdpack, lowDia)

66cov

67+1.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00

68+4.023281E-02, +1.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00

69+7.728843E-01, +6.530263E-01, +1.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00

70+4.914861E-01, +6.725434E-01, +8.497998E-01, +1.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00

71+9.849682E-02, -7.959484E-01, -4.249785E-01, -6.386417E-01, +1.000000E+00, +0.000000E+00, +0.000000E+00

72-5.807413E-01, +5.932847E-01, -7.603369E-02, +1.259757E-01, -6.656111E-01, +1.000000E+00, +0.000000E+00

73call 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.

74cov

75+1.000000E+00, +4.000000E+00, +2.413969E-01, +1.391192E+01, +3.931889E+00, +3.939873E-01, -2.322965E+00

76+4.023281E-02, +1.000000E+00, +9.000000E+00, +1.763171E+01, +8.070520E+00, -4.775691E+00, +3.559708E+00

77+7.728843E-01, +6.530263E-01, +1.000000E+00, +8.100000E+01, +3.059279E+01, -7.649613E+00, -1.368606E+00

78+4.914861E-01, +6.725434E-01, +8.497998E-01, +1.000000E+00, +1.600000E+01, -5.109134E+00, +1.007806E+00

79+9.849682E-02, -7.959484E-01, -4.249785E-01, -6.386417E-01, +1.000000E+00, +4.000000E+00, -2.662444E+00

80-5.807413E-01, +5.932847E-01, -7.603369E-02, +1.259757E-01, -6.656111E-01, +1.000000E+00, +4.000000E+00

81

82cov = getFilled(0._TKG, ndim, ndim + 1)

83call setMatCopy(cov(1:ndim, 2:ndim+1), rdpack, cor, rdpack, uppDia)

84cov

85+0.000000E+00, +1.000000E+00, +4.023281E-02, +7.728843E-01, +4.914861E-01, +9.849682E-02, -5.807413E-01

86+0.000000E+00, +0.000000E+00, +1.000000E+00, +6.530263E-01, +6.725434E-01, -7.959484E-01, +5.932847E-01

87+0.000000E+00, +0.000000E+00, +0.000000E+00, +1.000000E+00, +8.497998E-01, -4.249785E-01, -7.603369E-02

88+0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +1.000000E+00, -6.386417E-01, +1.259757E-01

89+0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +1.000000E+00, -6.656111E-01

90+0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +1.000000E+00

91call 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.

92cov

93+4.000000E+00, +1.000000E+00, +4.023281E-02, +7.728843E-01, +4.914861E-01, +9.849682E-02, -5.807413E-01

94+2.413969E-01, +9.000000E+00, +1.000000E+00, +6.530263E-01, +6.725434E-01, -7.959484E-01, +5.932847E-01

95+1.391192E+01, +1.763171E+01, +8.100000E+01, +1.000000E+00, +8.497998E-01, -4.249785E-01, -7.603369E-02

96+3.931889E+00, +8.070520E+00, +3.059279E+01, +1.600000E+01, +1.000000E+00, -6.386417E-01, +1.259757E-01

97+3.939873E-01, -4.775691E+00, -7.649613E+00, -5.109134E+00, +4.000000E+00, +1.000000E+00, -6.656111E-01

98-2.322965E+00, +3.559708E+00, -1.368606E+00, +1.007806E+00, -2.662444E+00, +4.000000E+00, +1.000000E+00

99

100

101ndim = getUnifRand(1, 7)

102ndim

103+6

104std = getUnifRand(1, 10, ndim)

105std

106+5.000000E+00, +3.000000E+00, +1.000000E+01, +6.000000E+00, +9.000000E+00, +6.000000E+00

107cor = getCovRand(1., ndim)

108cor

109+1.000000E+00, +5.897937E-01, -7.679903E-01, -4.463248E-03, +2.606820E-01, +1.045563E-01

110+5.897937E-01, +1.000000E+00, +5.799925E-02, +5.155175E-01, +5.788521E-01, -2.220319E-01

111-7.679903E-01, +5.799925E-02, +1.000000E+00, +4.797719E-01, +1.336182E-01, -2.744980E-01

112-4.463248E-03, +5.155175E-01, +4.797719E-01, +1.000000E+00, +1.287769E-01, +3.666028E-02

113+2.606820E-01, +5.788521E-01, +1.336182E-01, +1.287769E-01, +1.000000E+00, -1.377943E-01

114+1.045563E-01, -2.220319E-01, -2.744980E-01, +3.666028E-02, -1.377943E-01, +1.000000E+00

115cov = getFilled(0._TKG, ndim, ndim + 1)

116call setMatCopy(cov(1:ndim, 1:ndim), rdpack, cor, rdpack, lowDia)

117cov

118+1.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00

119+5.897937E-01, +1.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00

120-7.679903E-01, +5.799925E-02, +1.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00

121-4.463248E-03, +5.155175E-01, +4.797719E-01, +1.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00

122+2.606820E-01, +5.788521E-01, +1.336182E-01, +1.287769E-01, +1.000000E+00, +0.000000E+00, +0.000000E+00

123+1.045563E-01, -2.220319E-01, -2.744980E-01, +3.666028E-02, -1.377943E-01, +1.000000E+00, +0.000000E+00

124call 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.

125cov

126+1.000000E+00, +2.500000E+01, +8.846906E+00, -3.839951E+01, -1.338974E-01, +1.173069E+01, +3.136689E+00

127+5.897937E-01, +1.000000E+00, +9.000000E+00, +1.739978E+00, +9.279314E+00, +1.562901E+01, -3.996574E+00

128-7.679903E-01, +5.799925E-02, +1.000000E+00, +1.000000E+02, +2.878631E+01, +1.202564E+01, -1.646988E+01

129-4.463248E-03, +5.155175E-01, +4.797719E-01, +1.000000E+00, +3.600000E+01, +6.953952E+00, +1.319770E+00

130+2.606820E-01, +5.788521E-01, +1.336182E-01, +1.287769E-01, +1.000000E+00, +8.100000E+01, -7.440895E+00

131+1.045563E-01, -2.220319E-01, -2.744980E-01, +3.666028E-02, -1.377943E-01, +1.000000E+00, +3.600000E+01

132

133cov = getFilled(0._TKG, ndim, ndim + 1)

134call setMatCopy(cov(1:ndim, 2:ndim+1), rdpack, cor, rdpack, uppDia)

135cov

136+0.000000E+00, +1.000000E+00, +5.897937E-01, -7.679903E-01, -4.463248E-03, +2.606820E-01, +1.045563E-01

137+0.000000E+00, +0.000000E+00, +1.000000E+00, +5.799925E-02, +5.155175E-01, +5.788521E-01, -2.220319E-01

138+0.000000E+00, +0.000000E+00, +0.000000E+00, +1.000000E+00, +4.797719E-01, +1.336182E-01, -2.744980E-01

139+0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +1.000000E+00, +1.287769E-01, +3.666028E-02

140+0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +1.000000E+00, -1.377943E-01

141+0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +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+2.500000E+01, +1.000000E+00, +5.897937E-01, -7.679903E-01, -4.463248E-03, +2.606820E-01, +1.045563E-01

145+8.846906E+00, +9.000000E+00, +1.000000E+00, +5.799925E-02, +5.155175E-01, +5.788521E-01, -2.220319E-01

146-3.839951E+01, +1.739978E+00, +1.000000E+02, +1.000000E+00, +4.797719E-01, +1.336182E-01, -2.744980E-01

147-1.338974E-01, +9.279314E+00, +2.878631E+01, +3.600000E+01, +1.000000E+00, +1.287769E-01, +3.666028E-02

148+1.173069E+01, +1.562901E+01, +1.202564E+01, +6.953952E+00, +8.100000E+01, +1.000000E+00, -1.377943E-01

149+3.136689E+00, -3.996574E+00, -1.646988E+01, +1.319770E+00, -7.440895E+00, +3.600000E+01, +1.000000E+00

150

151

152ndim = getUnifRand(1, 7)

153ndim

154+4

155std = getUnifRand(1, 10, ndim)

156std

157+6.000000E+00, +9.000000E+00, +2.000000E+00, +3.000000E+00

158cor = getCovRand(1., ndim)

159cor

160+1.000000E+00, -8.879269E-01, +2.942848E-01, -6.429215E-01

161-8.879269E-01, +1.000000E+00, +4.709339E-02, +7.716588E-01

162+2.942848E-01, +4.709339E-02, +1.000000E+00, -1.943900E-02

163-6.429215E-01, +7.716588E-01, -1.943900E-02, +1.000000E+00

164cov = getFilled(0._TKG, ndim, ndim + 1)

165call setMatCopy(cov(1:ndim, 1:ndim), rdpack, cor, rdpack, lowDia)

166cov

167+1.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00

168-8.879269E-01, +1.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00

169+2.942848E-01, +4.709339E-02, +1.000000E+00, +0.000000E+00, +0.000000E+00

170-6.429215E-01, +7.716588E-01, -1.943900E-02, +1.000000E+00, +0.000000E+00

171call 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.

172cov

173+1.000000E+00, +3.600000E+01, -4.794806E+01, +3.531417E+00, -1.157259E+01

174-8.879269E-01, +1.000000E+00, +8.100000E+01, +8.476810E-01, +2.083479E+01

175+2.942848E-01, +4.709339E-02, +1.000000E+00, +4.000000E+00, -1.166340E-01

176-6.429215E-01, +7.716588E-01, -1.943900E-02, +1.000000E+00, +9.000000E+00

177

178cov = getFilled(0._TKG, ndim, ndim + 1)

179call setMatCopy(cov(1:ndim, 2:ndim+1), rdpack, cor, rdpack, uppDia)

180cov

181+0.000000E+00, +1.000000E+00, -8.879269E-01, +2.942848E-01, -6.429215E-01

182+0.000000E+00, +0.000000E+00, +1.000000E+00, +4.709339E-02, +7.716588E-01

183+0.000000E+00, +0.000000E+00, +0.000000E+00, +1.000000E+00, -1.943900E-02

184+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+3.600000E+01, +1.000000E+00, -8.879269E-01, +2.942848E-01, -6.429215E-01

188-4.794806E+01, +8.100000E+01, +1.000000E+00, +4.709339E-02, +7.716588E-01

189+3.531417E+00, +8.476810E-01, +4.000000E+00, +1.000000E+00, -1.943900E-02

190-1.157259E+01, +2.083479E+01, -1.166340E-01, +9.000000E+00, +1.000000E+00

191

192

193ndim = getUnifRand(1, 7)

194ndim

195+1

196std = getUnifRand(1, 10, ndim)

197std

198+8.000000E+00

199cor = getCovRand(1., ndim)

200cor

201+1.000000E+00

202cov = getFilled(0._TKG, ndim, ndim + 1)

203call setMatCopy(cov(1:ndim, 1:ndim), rdpack, cor, rdpack, lowDia)

204cov

205+1.000000E+00, +0.000000E+00

206call 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.

207cov

208+1.000000E+00, +6.400000E+01

209

210cov = getFilled(0._TKG, ndim, ndim + 1)

211call setMatCopy(cov(1:ndim, 2:ndim+1), rdpack, cor, rdpack, uppDia)

212cov

213+0.000000E+00, +1.000000E+00

214call 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.

215cov

216+6.400000E+01, +1.000000E+00

217

218

219ndim = getUnifRand(1, 7)

220ndim

221+5

222std = getUnifRand(1, 10, ndim)

223std

224+6.000000E+00, +6.000000E+00, +7.000000E+00, +1.000000E+01, +5.000000E+00

225cor = getCovRand(1., ndim)

226cor

227+1.000000E+00, -1.710908E-02, -5.197237E-01, +1.015936E-01, +4.044273E-01

228-1.710908E-02, +1.000000E+00, +6.244778E-01, +7.290409E-01, -4.117261E-01

229-5.197237E-01, +6.244778E-01, +1.000000E+00, +7.947495E-01, -7.498054E-01

230+1.015936E-01, +7.290409E-01, +7.947495E-01, +1.000000E+00, -5.567722E-01

231+4.044273E-01, -4.117261E-01, -7.498054E-01, -5.567722E-01, +1.000000E+00

232cov = getFilled(0._TKG, ndim, ndim + 1)

233call setMatCopy(cov(1:ndim, 1:ndim), rdpack, cor, rdpack, lowDia)

234cov

235+1.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00

236-1.710908E-02, +1.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00

237-5.197237E-01, +6.244778E-01, +1.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00

238+1.015936E-01, +7.290409E-01, +7.947495E-01, +1.000000E+00, +0.000000E+00, +0.000000E+00

239+4.044273E-01, -4.117261E-01, -7.498054E-01, -5.567722E-01, +1.000000E+00, +0.000000E+00

240call 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.

241cov

242+1.000000E+00, +3.600000E+01, -6.159270E-01, -2.182839E+01, +6.095614E+00, +1.213282E+01

243-1.710908E-02, +1.000000E+00, +3.600000E+01, +2.622807E+01, +4.374245E+01, -1.235178E+01

244-5.197237E-01, +6.244778E-01, +1.000000E+00, +4.900000E+01, +5.563247E+01, -2.624319E+01

245+1.015936E-01, +7.290409E-01, +7.947495E-01, +1.000000E+00, +1.000000E+02, -2.783861E+01

246+4.044273E-01, -4.117261E-01, -7.498054E-01, -5.567722E-01, +1.000000E+00, +2.500000E+01

247

248cov = getFilled(0._TKG, ndim, ndim + 1)

249call setMatCopy(cov(1:ndim, 2:ndim+1), rdpack, cor, rdpack, uppDia)

250cov

251+0.000000E+00, +1.000000E+00, -1.710908E-02, -5.197237E-01, +1.015936E-01, +4.044273E-01

252+0.000000E+00, +0.000000E+00, +1.000000E+00, +6.244778E-01, +7.290409E-01, -4.117261E-01

253+0.000000E+00, +0.000000E+00, +0.000000E+00, +1.000000E+00, +7.947495E-01, -7.498054E-01

254+0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +1.000000E+00, -5.567722E-01

255+0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +1.000000E+00

256call 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.

257cov

258+3.600000E+01, +1.000000E+00, -1.710908E-02, -5.197237E-01, +1.015936E-01, +4.044273E-01

259-6.159270E-01, +3.600000E+01, +1.000000E+00, +6.244778E-01, +7.290409E-01, -4.117261E-01

260-2.182839E+01, +2.622807E+01, +4.900000E+01, +1.000000E+00, +7.947495E-01, -7.498054E-01

261+6.095614E+00, +4.374245E+01, +5.563247E+01, +1.000000E+02, +1.000000E+00, -5.567722E-01

262+1.213282E+01, -1.235178E+01, -2.624319E+01, -2.783861E+01, +2.500000E+01, +1.000000E+00

263

264

265ndim = getUnifRand(1, 7)

266ndim

267+5

268std = getUnifRand(1, 10, ndim)

269std

270+7.000000E+00, +2.000000E+00, +8.000000E+00, +1.000000E+00, +6.000000E+00

271cor = getCovRand(1., ndim)

272cor

273+1.000000E+00, -6.294793E-01, +7.991309E-01, -7.867544E-01, -3.664601E-02

274-6.294793E-01, +1.000000E+00, -3.614026E-02, +6.684650E-01, -5.464165E-01

275+7.991309E-01, -3.614026E-02, +1.000000E+00, -4.956060E-01, -4.787447E-01

276-7.867544E-01, +6.684650E-01, -4.956060E-01, +1.000000E+00, -2.422076E-01

277-3.664601E-02, -5.464165E-01, -4.787447E-01, -2.422076E-01, +1.000000E+00

278cov = getFilled(0._TKG, ndim, ndim + 1)

279call setMatCopy(cov(1:ndim, 1:ndim), rdpack, cor, rdpack, lowDia)

280cov

281+1.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00

282-6.294793E-01, +1.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00

283+7.991309E-01, -3.614026E-02, +1.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00

284-7.867544E-01, +6.684650E-01, -4.956060E-01, +1.000000E+00, +0.000000E+00, +0.000000E+00

285-3.664601E-02, -5.464165E-01, -4.787447E-01, -2.422076E-01, +1.000000E+00, +0.000000E+00

286call 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.

287cov

288+1.000000E+00, +4.900000E+01, -8.812710E+00, +4.475133E+01, -5.507280E+00, -1.539132E+00

289-6.294793E-01, +1.000000E+00, +4.000000E+00, -5.782442E-01, +1.336930E+00, -6.556997E+00

290+7.991309E-01, -3.614026E-02, +1.000000E+00, +6.400000E+01, -3.964848E+00, -2.297975E+01

291-7.867544E-01, +6.684650E-01, -4.956060E-01, +1.000000E+00, +1.000000E+00, -1.453246E+00

292-3.664601E-02, -5.464165E-01, -4.787447E-01, -2.422076E-01, +1.000000E+00, +3.600000E+01

293

294cov = getFilled(0._TKG, ndim, ndim + 1)

295call setMatCopy(cov(1:ndim, 2:ndim+1), rdpack, cor, rdpack, uppDia)

296cov

297+0.000000E+00, +1.000000E+00, -6.294793E-01, +7.991309E-01, -7.867544E-01, -3.664601E-02

298+0.000000E+00, +0.000000E+00, +1.000000E+00, -3.614026E-02, +6.684650E-01, -5.464165E-01

299+0.000000E+00, +0.000000E+00, +0.000000E+00, +1.000000E+00, -4.956060E-01, -4.787447E-01

300+0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +1.000000E+00, -2.422076E-01

301+0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +1.000000E+00

302call 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.

303cov

304+4.900000E+01, +1.000000E+00, -6.294793E-01, +7.991309E-01, -7.867544E-01, -3.664601E-02

305-8.812710E+00, +4.000000E+00, +1.000000E+00, -3.614026E-02, +6.684650E-01, -5.464165E-01

306+4.475133E+01, -5.782442E-01, +6.400000E+01, +1.000000E+00, -4.956060E-01, -4.787447E-01

307-5.507280E+00, +1.336930E+00, -3.964848E+00, +1.000000E+00, +1.000000E+00, -2.422076E-01

308-1.539132E+00, -6.556997E+00, -2.297975E+01, -1.453246E+00, +3.600000E+01, +1.000000E+00

309

310

311ndim = getUnifRand(1, 7)

312ndim

313+2

314std = getUnifRand(1, 10, ndim)

315std

316+9.000000E+00, +8.000000E+00

317cor = getCovRand(1., ndim)

318cor

319+1.000000E+00, +9.729016E-01

320+9.729016E-01, +1.000000E+00

321cov = getFilled(0._TKG, ndim, ndim + 1)

322call setMatCopy(cov(1:ndim, 1:ndim), rdpack, cor, rdpack, lowDia)

323cov

324+1.000000E+00, +0.000000E+00, +0.000000E+00

325+9.729016E-01, +1.000000E+00, +0.000000E+00

326call 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.

327cov

328+1.000000E+00, +8.100000E+01, +7.004891E+01

329+9.729016E-01, +1.000000E+00, +6.400000E+01

330

331cov = getFilled(0._TKG, ndim, ndim + 1)

332call setMatCopy(cov(1:ndim, 2:ndim+1), rdpack, cor, rdpack, uppDia)

333cov

334+0.000000E+00, +1.000000E+00, +9.729016E-01

335+0.000000E+00, +0.000000E+00, +1.000000E+00

336call 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.

337cov

338+8.100000E+01, +1.000000E+00, +9.729016E-01

339+7.004891E+01, +6.400000E+01, +1.000000E+00

340

341

342ndim = getUnifRand(1, 7)

343ndim

344+4

345std = getUnifRand(1, 10, ndim)

346std

347+7.000000E+00, +1.000000E+00, +4.000000E+00, +7.000000E+00

348cor = getCovRand(1., ndim)

349cor

350+1.000000E+00, -3.329586E-01, +1.571624E-01, -6.655376E-01

351-3.329586E-01, +1.000000E+00, +5.734950E-01, +8.516300E-01

352+1.571624E-01, +5.734950E-01, +1.000000E+00, +4.781445E-01

353-6.655376E-01, +8.516300E-01, +4.781445E-01, +1.000000E+00

354cov = getFilled(0._TKG, ndim, ndim + 1)

355call setMatCopy(cov(1:ndim, 1:ndim), rdpack, cor, rdpack, lowDia)

356cov

357+1.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00

358-3.329586E-01, +1.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00

359+1.571624E-01, +5.734950E-01, +1.000000E+00, +0.000000E+00, +0.000000E+00

360-6.655376E-01, +8.516300E-01, +4.781445E-01, +1.000000E+00, +0.000000E+00

361call 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.

362cov

363+1.000000E+00, +4.900000E+01, -2.330710E+00, +4.400548E+00, -3.261134E+01

364-3.329586E-01, +1.000000E+00, +1.000000E+00, +2.293980E+00, +5.961410E+00

365+1.571624E-01, +5.734950E-01, +1.000000E+00, +1.600000E+01, +1.338805E+01

366-6.655376E-01, +8.516300E-01, +4.781445E-01, +1.000000E+00, +4.900000E+01

367

368cov = getFilled(0._TKG, ndim, ndim + 1)

369call setMatCopy(cov(1:ndim, 2:ndim+1), rdpack, cor, rdpack, uppDia)

370cov

371+0.000000E+00, +1.000000E+00, -3.329586E-01, +1.571624E-01, -6.655376E-01

372+0.000000E+00, +0.000000E+00, +1.000000E+00, +5.734950E-01, +8.516300E-01

373+0.000000E+00, +0.000000E+00, +0.000000E+00, +1.000000E+00, +4.781445E-01

374+0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +1.000000E+00

375call 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.

376cov

377+4.900000E+01, +1.000000E+00, -3.329586E-01, +1.571624E-01, -6.655376E-01

378-2.330710E+00, +1.000000E+00, +1.000000E+00, +5.734950E-01, +8.516300E-01

379+4.400548E+00, +2.293980E+00, +1.600000E+01, +1.000000E+00, +4.781445E-01

380-3.261134E+01, +5.961410E+00, +1.338805E+01, +4.900000E+01, +1.000000E+00

381

382

383ndim = getUnifRand(1, 7)

384ndim

385+6

386std = getUnifRand(1, 10, ndim)

387std

388+7.000000E+00, +3.000000E+00, +1.000000E+00, +8.000000E+00, +3.000000E+00, +8.000000E+00

389cor = getCovRand(1., ndim)

390cor

391+1.000000E+00, -6.865154E-01, -3.722259E-01, -4.555238E-01, +5.617024E-01, -4.812185E-01

392-6.865154E-01, +1.000000E+00, -4.175428E-02, +3.215826E-01, -5.951670E-01, +6.668763E-01

393-3.722259E-01, -4.175428E-02, +1.000000E+00, -3.474011E-01, -3.992162E-01, +4.134352E-01

394-4.555238E-01, +3.215826E-01, -3.474011E-01, +1.000000E+00, -3.370606E-01, -2.603529E-01

395+5.617024E-01, -5.951670E-01, -3.992162E-01, -3.370606E-01, +1.000000E+00, -4.563255E-01

396-4.812185E-01, +6.668763E-01, +4.134352E-01, -2.603529E-01, -4.563255E-01, +1.000000E+00

397cov = getFilled(0._TKG, ndim, ndim + 1)

398call setMatCopy(cov(1:ndim, 1:ndim), rdpack, cor, rdpack, lowDia)

399cov

400+1.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00

401-6.865154E-01, +1.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00

402-3.722259E-01, -4.175428E-02, +1.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00

403-4.555238E-01, +3.215826E-01, -3.474011E-01, +1.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00

404+5.617024E-01, -5.951670E-01, -3.992162E-01, -3.370606E-01, +1.000000E+00, +0.000000E+00, +0.000000E+00

405-4.812185E-01, +6.668763E-01, +4.134352E-01, -2.603529E-01, -4.563255E-01, +1.000000E+00, +0.000000E+00

406call 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.

407cov

408+1.000000E+00, +4.900000E+01, -1.441682E+01, -2.605582E+00, -2.550933E+01, +1.179575E+01, -2.694823E+01

409-6.865154E-01, +1.000000E+00, +9.000000E+00, -1.252628E-01, +7.717983E+00, -5.356503E+00, +1.600503E+01

410-3.722259E-01, -4.175428E-02, +1.000000E+00, +1.000000E+00, -2.779209E+00, -1.197649E+00, +3.307481E+00

411-4.555238E-01, +3.215826E-01, -3.474011E-01, +1.000000E+00, +6.400000E+01, -8.089455E+00, -1.666258E+01

412+5.617024E-01, -5.951670E-01, -3.992162E-01, -3.370606E-01, +1.000000E+00, +9.000000E+00, -1.095181E+01

413-4.812185E-01, +6.668763E-01, +4.134352E-01, -2.603529E-01, -4.563255E-01, +1.000000E+00, +6.400000E+01

414

415cov = getFilled(0._TKG, ndim, ndim + 1)

416call setMatCopy(cov(1:ndim, 2:ndim+1), rdpack, cor, rdpack, uppDia)

417cov

418+0.000000E+00, +1.000000E+00, -6.865154E-01, -3.722259E-01, -4.555238E-01, +5.617024E-01, -4.812185E-01

419+0.000000E+00, +0.000000E+00, +1.000000E+00, -4.175428E-02, +3.215826E-01, -5.951670E-01, +6.668763E-01

420+0.000000E+00, +0.000000E+00, +0.000000E+00, +1.000000E+00, -3.474011E-01, -3.992162E-01, +4.134352E-01

421+0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +1.000000E+00, -3.370606E-01, -2.603529E-01

422+0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +1.000000E+00, -4.563255E-01

423+0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +0.000000E+00, +1.000000E+00

424call 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.

425cov

426+4.900000E+01, +1.000000E+00, -6.865154E-01, -3.722259E-01, -4.555238E-01, +5.617024E-01, -4.812185E-01

427-1.441682E+01, +9.000000E+00, +1.000000E+00, -4.175428E-02, +3.215826E-01, -5.951670E-01, +6.668763E-01

428-2.605582E+00, -1.252628E-01, +1.000000E+00, +1.000000E+00, -3.474011E-01, -3.992162E-01, +4.134352E-01

429-2.550933E+01, +7.717983E+00, -2.779209E+00, +6.400000E+01, +1.000000E+00, -3.370606E-01, -2.603529E-01

430+1.179575E+01, -5.356503E+00, -1.197649E+00, -8.089455E+00, +9.000000E+00, +1.000000E+00, -4.563255E-01

431-2.694823E+01, +1.600503E+01, +3.307481E+00, -1.666258E+01, -1.095181E+01, +6.400000E+01, +1.000000E+00

432

433

434!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

435!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

436!Compute the covariance matrix of a 2-D sample.

437!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

438!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

439

440ndim = 2; nsam = 10; dim = 2

441sample = reshape(getUnifRand(1, 20, ndim * nsam), shape = [ndim, nsam], order = [2, 1])

442sample

443+15.0000000, +14.0000000, +2.00000000, +5.00000000, +5.00000000, +2.00000000, +8.00000000, +3.00000000, +17.0000000, +15.0000000

444+7.00000000, +17.0000000, +20.0000000, +5.00000000, +18.0000000, +17.0000000, +8.00000000, +6.00000000, +7.00000000, +16.0000000

445mean = getMean(sample, dim)

446mean

447+8.60000038, +12.1000004

448samShifted = getShifted(sample, dim, -mean)

449samShifted

450+6.39999962, +5.39999962, -6.60000038, -3.60000038, -3.60000038, -6.60000038, -0.600000381, -5.60000038, +8.39999962, +6.39999962

451-5.10000038, +4.89999962, +7.89999962, -7.10000038, +5.89999962, +4.89999962, -4.10000038, -6.10000038, -5.10000038, +3.89999962

452cov = getFilled(0., ndim, ndim)

453call setCov(cov, uppDia, mean, sample, dim)

454cov

455+32.6399994, -6.76000071

456+0.00000000, +31.6900005

457cov = getFilled(0., ndim, ndim)

458call setCov(cov, uppDia, samShifted, dim) ! same result as above.

459cov

460+32.6399994, -6.76000071

461+0.00000000, +31.6900005

462

463cov = getFilled(0., ndim, ndim)

464call setCov(cov, lowDia, mean, sample, dim)

465cov

466+32.6399994, +0.00000000

467-6.76000071, +31.6900005

468cov = getFilled(0., ndim, ndim)

469call setCov(cov, lowDia, samShifted, dim) ! same result as above.

470cov

471+32.6399994, +0.00000000

472-6.76000071, +31.6900005

473

474'Compute the sample covariance along the first dimension.'

475

476dim = 1

477cov = getFilled(0., ndim, ndim)

478call setCov(cov, uppDia, mean, transpose(sample), dim)

479cov

480+32.6399994, -6.76000071

481+0.00000000, +31.6900005

482cov = getFilled(0., ndim, ndim)

483call setCov(cov, uppDia, transpose(samShifted), dim) ! same result as above.

484cov

485+32.6399994, -6.76000071

486+0.00000000, +31.6900005

487

488cov = getFilled(0., ndim, ndim)

489call setCov(cov, lowDia, mean, transpose(sample), dim)

490cov

491+32.6399994, +0.00000000

492-6.76000071, +31.6900005

493cov = getFilled(0., ndim, ndim)

494call setCov(cov, lowDia, transpose(samShifted), dim) ! same result as above.

495cov

496+32.6399994, +0.00000000

497-6.76000071, +31.6900005

498

499'Compute the full sample covariance for a pair of time series.'

500

501cov = getFilled(0., ndim, ndim)

502call setCov(cov, mean, sample(1,:), sample(2,:))

503cov

504+32.6399994, -6.76000071

505-6.76000071, +31.6900005

506cov = getFilled(0., ndim, ndim)

507call setCov(cov, samShifted(1,:), samShifted(2,:)) ! same result as above.

508cov

509+32.6399994, -6.76000071

510-6.76000071, +31.6900005

511

512

513!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

514!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

515!Compute the biased covariance matrix of a weighted 2-D sample.

516!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

517!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

518

519ndim = 2; nsam = 10; dim = 2

520sample = reshape(getUnifRand(1, 20, ndim * nsam), shape = [ndim, nsam], order = [2, 1])

521sample

522+13.0000000, +17.0000000, +17.0000000, +8.00000000, +12.0000000, +12.0000000, +4.00000000, +12.0000000, +11.0000000, +4.00000000

523+6.00000000, +14.0000000, +10.0000000, +10.0000000, +12.0000000, +3.00000000, +7.00000000, +13.0000000, +7.00000000, +2.00000000

524call setResized(mean, ndim)

525iweight = getUnifRand(1, 10, nsam) ! integer-valued weights.

526iweight

527+7, +5, +4, +1, +4, +9, +10, +8, +6, +5

528call setMean(mean, sample, dim, iweight, iweisum)

529mean

530+10.6779661, +7.84745741

531iweisum

532+59

533rweight = iweight ! or real-valued weights.

534iweight

535+7, +5, +4, +1, +4, +9, +10, +8, +6, +5

536call setMean(mean, sample, dim, rweight, rweisum)

537mean

538+10.6779661, +7.84745741

539rweisum

540+59.0000000

541

542!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

543!Compute the covariance matrix integer weights.

544!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

545

546samShifted = getShifted(sample, dim, -mean)

547samShifted

548+2.32203388, +6.32203388, +6.32203388, -2.67796612, +1.32203388, +1.32203388, -6.67796612, +1.32203388, +0.322033882, -6.67796612

549-1.84745741, +6.15254259, +2.15254259, +2.15254259, +4.15254259, -4.84745741, -0.847457409, +5.15254259, -0.847457409, -5.84745741

550cov = getFilled(0., ndim, ndim)

551call setCov(cov, uppDia, mean, sample, dim, iweight, iweisum)

552cov

553+18.8284969, +8.17121601

554+0.00000000, +15.4513073

555cov = getFilled(0., ndim, ndim)

556call setCov(cov, uppDia, samShifted, dim, iweight, sum(iweight)) ! same result as above.

557cov

558+18.8284969, +8.17121601

559+0.00000000, +15.4513073

560

561cov = getFilled(0., ndim, ndim)

562call setCov(cov, lowDia, mean, sample, dim, iweight, iweisum)

563cov

564+18.8284969, +0.00000000

565+8.17121601, +15.4513073

566cov = getFilled(0., ndim, ndim)

567call setCov(cov, lowDia, samShifted, dim, iweight, sum(iweight)) ! same result as above.

568cov

569+18.8284969, +0.00000000

570+8.17121601, +15.4513073

571

572'Compute the sample covariance along the first dimension.'

573

574dim = 1

575cov = getFilled(0., ndim, ndim)

576call setCov(cov, uppDia, mean, transpose(sample), dim, iweight, iweisum)

577cov

578+18.8284969, +8.17121506

579+0.00000000, +15.4513063

580cov = getFilled(0., ndim, ndim)

581call setCov(cov, uppDia, transpose(samShifted), dim, iweight, sum(iweight)) ! same result as above.

582cov

583+18.8284969, +8.17121506

584+0.00000000, +15.4513063

585

586cov = getFilled(0., ndim, ndim)

587call setCov(cov, lowDia, mean, transpose(sample), dim, iweight, iweisum)

588cov

589+18.8284969, +0.00000000

590+8.17121506, +15.4513063

591cov = getFilled(0., ndim, ndim)

592call setCov(cov, lowDia, transpose(samShifted), dim, iweight, sum(iweight)) ! same result as above.

593cov

594+18.8284969, +0.00000000

595+8.17121506, +15.4513063

596

597'Compute the full sample covariance for a pair of time series.'

598

599cov = getFilled(0., ndim, ndim)

600call setCov(cov, mean, sample(1,:), sample(2,:), iweight, iweisum)

601cov

602+18.8284969, +8.17121601

603+8.17121601, +15.4513073

604cov = getFilled(0., ndim, ndim)

605call setCov(cov, samShifted(1,:), samShifted(2,:), iweight, sum(iweight)) ! same result as above.

606cov

607+18.8284969, +8.17121601

608+8.17121601, +15.4513073

609

610

611!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

612!Compute the covariance matrix real weights.

613!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

614

615dim = 2_IK

616samShifted = getShifted(sample, dim, -mean)

617samShifted

618+2.32203388, +6.32203388, +6.32203388, -2.67796612, +1.32203388, +1.32203388, -6.67796612, +1.32203388, +0.322033882, -6.67796612

619-1.84745741, +6.15254259, +2.15254259, +2.15254259, +4.15254259, -4.84745741, -0.847457409, +5.15254259, -0.847457409, -5.84745741

620cov = getFilled(0., ndim, ndim)

621call setCov(cov, uppDia, mean, sample, dim, rweight, rweisum)

622cov

623+18.8284969, +8.17121601

624+0.00000000, +15.4513073

625cov = getFilled(0., ndim, ndim)

626call setCov(cov, uppDia, samShifted, dim, rweight, sum(rweight)) ! same result as above.

627cov

628+18.8284969, +8.17121601

629+0.00000000, +15.4513073

630

631cov = getFilled(0., ndim, ndim)

632call setCov(cov, lowDia, mean, sample, dim, rweight, rweisum)

633cov

634+18.8284969, +0.00000000

635+8.17121601, +15.4513073

636cov = getFilled(0., ndim, ndim)

637call setCov(cov, lowDia, samShifted, dim, rweight, sum(rweight)) ! same result as above.

638cov

639+18.8284969, +0.00000000

640+8.17121601, +15.4513073

641

642'Compute the sample covariance along the first dimension.'

643

644dim = 1

645cov = getFilled(0., ndim, ndim)

646call setCov(cov, uppDia, mean, transpose(sample), dim, rweight, rweisum)

647cov

648+18.8284969, +8.17121506

649+0.00000000, +15.4513063

650cov = getFilled(0., ndim, ndim)

651call setCov(cov, uppDia, transpose(samShifted), dim, rweight, sum(rweight)) ! same result as above.

652cov

653+18.8284969, +8.17121506

654+0.00000000, +15.4513063

655

656cov = getFilled(0., ndim, ndim)

657call setCov(cov, lowDia, mean, transpose(sample), dim, rweight, rweisum)

658cov

659+18.8284969, +0.00000000

660+8.17121506, +15.4513063

661cov = getFilled(0., ndim, ndim)

662call setCov(cov, lowDia, transpose(samShifted), dim, rweight, sum(rweight)) ! same result as above.

663cov

664+18.8284969, +0.00000000

665+8.17121506, +15.4513063

666

667'Compute the full sample covariance for a pair of time series.'

668

669cov = getFilled(0., ndim, ndim)

670call setCov(cov, mean, sample(1,:), sample(2,:), rweight, rweisum)

671cov

672+18.8284969, +8.17121601

673+8.17121601, +15.4513073

674cov = getFilled(0., ndim, ndim)

675call setCov(cov, samShifted(1,:), samShifted(2,:), rweight, sum(rweight)) ! same result as above.

676cov

677+18.8284969, +8.17121601

678+8.17121601, +15.4513073

679

680

Benchmarks:

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

! Test the performance of Cholesky factorization computation using an assumed-shape interface vs. explicit-shape interface.
program benchmark
 
    use pm_kind, only: IK, LK, RKG => RKD, SK
    use pm_sampleCov, only: uppDia
    use pm_bench, only: bench_type
 
    implicit none
 
    integer(IK)                         :: itry, ntry
    integer(IK)                         :: i
    integer(IK)                         :: iarr
    integer(IK)                         :: fileUnit
    integer(IK)     , parameter         :: NARR = 18_IK
    integer(IK)     , allocatable       :: rperm(:)
    real(RKG)       , allocatable       :: samdim1(:,:)
    real(RKG)       , allocatable       :: samdim2(:,:)
    type(bench_type), allocatable       :: bench(:)
    integer(IK)     , parameter         :: nsammax = 2**NARR
    integer(IK)     , parameter         :: ndim = 5_IK
    real(RKG)                           :: mean(ndim), cov(ndim, ndim)
    integer(IK)                         :: idim, jdim, isam, nsam
    real(RKG)                           :: dumm
 
    bench = [ bench_type(name = SK_"intrinsicDIM1", exec = intrinsicDIM1, overhead = setOverhead) &
            , bench_type(name = SK_"intrinsicDIM2", exec = intrinsicDIM2, overhead = setOverhead) &
            , bench_type(name = SK_"setCovDIM1", exec = setCovDIM1, overhead = setOverhead) &
            , bench_type(name = SK_"setCovDIM2", exec = setCovDIM2, overhead = setOverhead) &
            ]
 
    write(*,"(*(g0,:,' '))")
    write(*,"(*(g0,:,' '))") "sample covariance benchmarking..."
    write(*,"(*(g0,:,' '))")
 
    open(newunit = fileUnit, file = "main.out", status = "replace")
 
        write(fileUnit, "(*(g0,:,','))") "nsam", (bench(i)%name, i = 1, size(bench))
 
        dumm = 0._RKG
        loopOverMatrixSize: do iarr = 1, NARR - 1
 
            nsam = 2**iarr
            ntry = nsammax / nsam
            allocate(samdim1(nsam, ndim), samdim2(ndim, nsam))
            write(*,"(*(g0,:,' '))") "Benchmarking setCovDIM1() vs. setCovDIM2()", nsam, ntry
 
            do i = 1, size(bench)
                bench(i)%timing = bench(i)%getTiming()
            end do
 
            write(fileUnit,"(*(g0,:,','))") nsam, (bench(i)%timing%mean / ntry, i = 1, size(bench))
            deallocate(samdim1, samdim2)
 
        end do loopOverMatrixSize
        write(*,"(*(g0,:,' '))") dumm
 
    close(fileUnit)
 
contains
 
    !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    ! procedure wrappers.
    !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 
    subroutine setOverhead()
        do itry = 1, ntry
            call random_number(mean)
            dumm = dumm - mean(1)
        end do
    end subroutine
 
    subroutine setSamDIM1()
        do isam = 1, nsam
            call random_number(samdim1(isam, 1 : ndim))
        end do
    end subroutine
 
    subroutine setSamDIM2()
        do isam = 1, nsam
            call random_number(samdim2(1 : ndim, isam))
        end do
    end subroutine
 
    subroutine setCovDIM1()
        block
            use pm_sampleCov, only: setCov
            do itry = 1, ntry
                call setSamDIM1()
                call setCov(cov, uppDia, samdim1, dim = 1_IK)
                dumm = dumm + cov(1,1)
            end do
        end block
    end subroutine
 
    subroutine setCovDIM2()
        block
            use pm_sampleCov, only: setCov
            do itry = 1, ntry
                call setSamDIM2()
                call setCov(cov, uppDia, samdim2, dim = 2_IK)
                dumm = dumm + cov(1,1)
            end do
        end block
    end subroutine
 
    subroutine intrinsicDIM1()
        do itry = 1, ntry
            call setSamDIM1()
            do jdim = 1, ndim
                do idim = 1, jdim
                    cov(idim, jdim) = dot_product(samdim1(1 : nsam, idim), samdim1(1 : nsam, jdim)) / size(samdim1, dim = 1)
                end do
            end do
            dumm = dumm + cov(1,1)
        end do
    end subroutine
 
    subroutine intrinsicDIM2()
        do itry = 1, ntry
            call setSamDIM2()
            do jdim = 1, ndim
                do idim = 1, jdim
                    cov(idim, jdim) = dot_product(samdim2(idim, 1 : nsam), samdim2(jdim, 1 : nsam)) / size(samdim2, dim = 2)
                end do
            end do
            dumm = dumm + cov(1,1)
        end do
    end subroutine
 
end program benchmark

Example Unix compile command via Intel ifort compiler ⛓

#!/usr/bin/env sh
rm main.exe
ifort -fpp -standard-semantics -O3 -Wl,-rpath,../../../lib -I../../../inc main.F90 ../../../lib/libparamonte* -o main.exe
./main.exe

Example Windows Batch compile command via Intel ifort compiler ⛓

del main.exe
set PATH=..\..\..\lib;%PATH%
ifort /fpp /standard-semantics /O3 /I:..\..\..\include main.F90 ..\..\..\lib\libparamonte*.lib /exe:main.exe
main.exe

Example Unix / MinGW compile command via GNU gfortran compiler ⛓

#!/usr/bin/env sh
rm main.exe
gfortran -cpp -ffree-line-length-none -O3 -Wl,-rpath,../../../lib -I../../../inc main.F90 ../../../lib/libparamonte* -o main.exe
./main.exe

Postprocessing of the benchmark output ⛓

#!/usr/bin/env python
 
import matplotlib.pyplot as plt
import pandas as pd
import numpy as np
 
import os
dirname = os.path.basename(os.getcwd()) 
 
fontsize = 14
 
df = pd.read_csv("main.out", delimiter = ",")
colnames = list(df.columns.values)
 
 
 
ax = plt.figure(figsize = 1.25 * np.array([6.4,4.6]), dpi = 200)
ax = plt.subplot()
 
for colname in colnames[1:]:
    plt.plot( df[colnames[0]].values
            , df[colname].values
            , linewidth = 2
            )
 
plt.xticks(fontsize = fontsize)
plt.yticks(fontsize = fontsize)
ax.set_xlabel(colnames[0], fontsize = fontsize)
ax.set_ylabel("Runtime [ seconds ]", fontsize = fontsize)
ax.set_title(" vs. ".join(colnames[1:])+"\nLower is better.", fontsize = fontsize)
ax.set_xscale("log")
ax.set_yscale("log")
plt.minorticks_on()
plt.grid(visible = True, which = "both", axis = "both", color = "0.85", linestyle = "-")
ax.tick_params(axis = "y", which = "minor")
ax.tick_params(axis = "x", which = "minor")
ax.legend   ( colnames[1:]
           #, loc='center left'
           #, bbox_to_anchor=(1, 0.5)
            , fontsize = fontsize
            )
 
plt.tight_layout()
plt.savefig("benchmark." + dirname + ".runtime.png")
 
 
 
ax = plt.figure(figsize = 1.25 * np.array([6.4,4.6]), dpi = 200)
ax = plt.subplot()
 
plt.plot( df[colnames[0]].values
        , np.ones(len(df[colnames[0]].values))
        , linestyle = "--"
       #, color = "black"
        , linewidth = 2
        )
for colname in colnames[2:]:
    plt.plot( df[colnames[0]].values
            , df[colname].values / df[colnames[1]].values
            , linewidth = 2
            )
 
plt.xticks(fontsize = fontsize)
plt.yticks(fontsize = fontsize)
ax.set_xlabel(colnames[0], fontsize = fontsize)
ax.set_ylabel("Runtime compared to {}".format(colnames[1]), fontsize = fontsize)
ax.set_title("Runtime Ratio Comparison. Lower means faster.\nLower than 1 means faster than {}().".format(colnames[1]), fontsize = fontsize)
ax.set_xscale("log")
ax.set_yscale("log")
plt.minorticks_on()
plt.grid(visible = True, which = "both", axis = "both", color = "0.85", linestyle = "-")
ax.tick_params(axis = "y", which = "minor")
ax.tick_params(axis = "x", which = "minor")
ax.legend   ( colnames[1:]
           #, bbox_to_anchor = (1, 0.5)
           #, loc = "center left"
            , fontsize = fontsize
            )
 
plt.tight_layout()
plt.savefig("benchmark." + dirname + ".runtime.ratio.png")

Visualization of the benchmark output ⛓

Benchmark moral ⛓

The procedures under the generic interface setCov can compute the covariance under two different sample axes.
Recall that C is a row-major language while Fortran is a column-major language.
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.
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.
The higher the number of data attributes, the more likely the computations along the second axis of sample will be faster.
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.

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.
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.

Copyright: Computational Data Science Lab

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:

src/fortran/main/pm_sampleCov.F90