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

Generate and return the (Pearson) correlation coefficient or matrix of a pair of (weighted) time series x(1:nsam) and y(1:nsam) or of an input (weighted) array of shape (ndim, nsam) or (nsam, ndim) where ndim is the number of data dimensions (the number of data attributes) and nsam is the number of data points.
More...

Detailed Description

Generate and return the (Pearson) correlation coefficient or matrix of a pair of (weighted) time series x(1:nsam) and y(1:nsam) or of an input (weighted) array of shape (ndim, nsam) or (nsam, ndim) 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 one of the following computational tasks:

  1. compute the correlation matrix corresponding to an input covariance matrix and vector of standard deviations in arbitrary ndim dimensions.
  2. compute the sample correlation matrix of a random sample of nsam observations in arbitrary ndim dimensional space.
  3. compute the sample correlation matrix of a pair of nsam observations in two separated data vectors x and y.

See the documentation of the parent module pm_sampleCor for algorithmic details and sample correlation matrix definition.

Note
The sample correlation matrix can be also readily computed from the sample covariance matrix as in the following example,
cor(:,:) = getCor(getCov(sample, dim))
cor(:,:) = getCor(getCov(sample, dim, weight)) ! weighted sample.
!
Parameters
[in]cov: The input positive semi-definite square matrix of shape (1:ndim, 1:ndim) of the same type and kind as the output cor, representing the correlation matrix based upon which the output correlation matrix cor is to be computed.
(optional. It must be present if and only if the input arguments sample and x and y are missing.)
[in]subsetv: 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 covariance matrix must be used.
  2. The constant uppDia, implying that only the upper-diagonal subset of the input covariance matrix must be used.
  3. The constant low, implying that only the lower subset of the input covariance matrix must be used.
    This option is available if and only if the input argument stdinv is present (otherwise, how do we know the variances?).
  4. The constant upp, implying that only the upper subset of the input covariance matrix must be used.
    This option is available if and only if the input argument stdinv is present (otherwise, how do we know the variances?).
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 argument cov is present.)
[in]stdinv: The input positive vector of shape (1 : ndim) of type real of the same kind as the input cor, containing the inverse of the standard deviations of the ndim data attributes based upon which the output correlation matrix cor is to be computed.

\begin{equation} \ms{stdinv}(i) = \frac{1} {\sqrt{\ms{cov}(i,i)}} ~, \end{equation}

(optional. It must be present if and only if the input argument subsetv is set to either low or upp.)
[in]x: The input contiguous vector of shape (nsam) of the same type and kind as the output cor, 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 sample is missing.)
[in]y: The input contiguous vector of shape (nsam) of the same type and kind as the output cor, 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 sample is missing.)
[in]sample: The input contiguous array of shape (ndim, nsam) or (nsam, ndim) of the same type and kind as the output cor, containing the sample comprised of nsam observations each with ndim attributes.
If sample is a matrix, then the direction along which the correlation 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 correlation 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 cor,
containing the corresponding weights of individual nsam observations in sample.
(optional. default = getFilled(1, nsam). It can be present only if the input arguments cov, subsetv, and stdinv are missing.)
Returns
cor : 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 its Pearson correlation matrix.
When the input arguments x and y are present, then ndim == 2 holds by definition and the output cor is a scalar of value \(r_{xy}\),

\begin{equation} \ms{cor} = \begin{bmatrix} 1 && r_{xy} \\ r_{yx} && 1 ~. \end{bmatrix} \end{equation}


Possible calling interfaces

use pm_sampleCor, only: getCor
! correlation matrix to correlation matrix.
cor(1:ndim, 1:ndim) = getCor(cov(1:ndim, 1:ndim), subsetv) ! subsetv = uppDia, lowDia
cor(1:ndim, 1:ndim) = getCor(cov(1:ndim, 1:ndim), subsetv, stdinv(1:ndim)) ! subsetv = upp, low
! XY time series correlation matrix.
cor = getCor(x(1:nsam), y(1:nsam)) ! correlation coefficient.
cor = getCor(x(1:nsam), y(1:nsam), weight(1:nsam)) ! correlation coefficient.
! sample correlation matrix.
cor(:,:) = getCor(sample(:,:), dim) ! full correlation matrix.
cor(:,:) = getCor(sample(:,:), dim, weight(1:nsam)) ! upper/lower correlation matrix.
!
Generate and return the (Pearson) correlation coefficient or matrix of a pair of (weighted) time seri...
This module contains classes and procedures for computing properties related to the correlation matri...
Warning
This generic interface is merely a flexible wrapper around the generic subroutine interface setCor.
As such, all conditions and warnings associated with setCor equally hold for this generic interface.
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.
See also
getCor
setCor
getRho
setRho
getCov
setCov
setECDF
getMean
setMean
getShifted
setShifted
getVar
setVar


Example usage

1program example
2
3 use pm_kind, only: SK, IK, LK
4 use pm_io, only: display_type
5 use pm_distUnif, only: getUnifRand
6 use pm_distCov, only: getCovRand
7 use pm_arrayRange, only: getRange
8 use pm_sampleCor, only: getCor
9 use pm_sampleCor, only: uppDia
10 use pm_sampleCor, only: lowDia
11 use pm_sampleCor, only: upp
12 use pm_sampleCor, only: low
13 use pm_sampleCov, only: getCov
14 use pm_sampleMean, only: getMean
15 use pm_sampleMean, only: setMean
16 use pm_sampleShift, only: getShifted
17 use pm_arraySpace, only: getLinSpace
18 use pm_arrayResize, only: setResized
19 use pm_arrayFill, only: getFilled
20 use pm_io, only: getFormat
21
22 implicit none
23
24 type(display_type) :: disp
25 integer(IK) :: itry, ntry = 10
26 character(:), allocatable :: format
27 disp = display_type(file = "main.out.F90")
28
29 call disp%skip()
30 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
31 call disp%show("!Compute the correlation matrix from covariance matrix.")
32 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
33 call disp%skip()
34
35 block
36 use pm_kind, only: TKG => RKS ! All other real types are also supported.
37 real(TKG), allocatable :: cov(:,:), cor(:,:), std(:)
38 integer(IK) :: ndim
39 do itry = 1, 10
40 call disp%skip()
41 call disp%show("ndim = getUnifRand(0, 7)")
42 ndim = getUnifRand(0, 7)
43 call disp%show("ndim")
44 call disp%show( ndim )
45 call disp%show("std = getUnifRand(1, ndim, ndim)")
46 std = getUnifRand(1, ndim, ndim)
47 call disp%show("std")
48 call disp%show( std )
49 call disp%show("cov = getCovRand(1._TKG, std)")
50 cov = getCovRand(1._TKG, std)
51 call disp%show("cov")
52 call disp%show( cov )
53 call disp%show("cor = getCor(cov, uppDia)")
54 cor = getCor(cov, uppDia)
55 call disp%show("cor")
56 call disp%show( cor )
57 call disp%show("cor = getCor(cov, upp, stdinv = 1 / std)")
58 cor = getCor(cov, upp, stdinv = 1 / std)
59 call disp%show("cor")
60 call disp%show( cor )
61 call disp%show("cor = getCor(cov, lowDia)")
62 cor = getCor(cov, lowDia)
63 call disp%show("cor")
64 call disp%show( cor )
65 call disp%show("cor = getCor(cov, low, stdinv = 1 / std)")
66 cor = getCor(cov, low, stdinv = 1 / std)
67 call disp%show("cor")
68 call disp%show( cor )
69 call disp%show("getCov(getCor(cov, lowDia), lowDia, std) ! reconstruct the original covariance matrix.")
70 call disp%show( getCov(getCor(cov, lowDia), lowDia, std) , format = format )
71 call disp%show("getCov(getCor(cov, uppDia), uppDia, std) ! reconstruct the original covariance matrix.")
72 call disp%show( getCov(getCor(cov, uppDia), uppDia, std) , format = format )
73 call disp%skip()
74 end do
75 end block
76
77 call disp%skip()
78 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
79 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
80 call disp%show("!Compute the Pearson correlation matrix for a pair of time series.")
81 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
82 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
83 call disp%skip()
84
85 block
86 use pm_kind, only: TKG => RKS ! All other real types are also supported.
87 integer(IK) :: ndim, nsam, dim
88 real(TKG), allocatable :: sample(:,:), cor(:,:), mean(:)
89 format = getFormat(mold = [0._TKG], ed = SK_"es", signed = .true._LK)
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 , format = format )
96 call disp%show("mean = getMean(sample, dim)")
97 mean = getMean(sample, dim)
98 call disp%show("mean")
99 call disp%show( mean , format = format )
100 call disp%show("cor = getCor(sample, dim) ! same result as above.")
101 cor = getCor(sample, dim)
102 call disp%show("cor")
103 call disp%show( cor , format = format )
104 call disp%skip()
105 call disp%show("Compute the sample correlation along the first dimension.", deliml = SK_'''')
106 call disp%skip()
107 call disp%show("dim = 1")
108 dim = 1
109 call disp%show("cor = getCor(transpose(sample), dim) ! same result as above.")
110 cor = getCor(transpose(sample), dim)
111 call disp%show("cor")
112 call disp%show( cor , format = format )
113 call disp%skip()
114 call disp%show("Compute the full sample correlation for a pair of time series.", deliml = SK_'''')
115 call disp%skip()
116 call disp%show("cor(1,1) = getCor(sample(1,:), sample(2,:)) ! same result as above.")
117 cor(1,1) = getCor(sample(1,:), sample(2,:))
118 call disp%show("cor(1,1)")
119 call disp%show( cor(1,1) , format = format )
120 call disp%skip()
121 end block
122
123 call disp%skip()
124 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
125 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
126 call disp%show("!Compute the Pearson correlation matrix for a weighted pair of time series.")
127 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
128 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
129 call disp%skip()
130
131 block
132 use pm_arrayVerbose, only: getVerbose
133 use pm_kind, only: TKG => RKS ! All other real types are also supported.
134 integer(IK), allocatable :: iweight(:)
135 real(TKG), allocatable :: rweight(:)
136 real(TKG) :: rweisum
137 integer(IK) :: iweisum
138 integer(IK) :: ndim, nsam, dim
139 real(TKG), allocatable :: sample(:,:), cor(:,:), mean(:)
140 format = getFormat(mold = [0._TKG], ed = SK_"es", signed = .true._LK)
141 call disp%show("ndim = 2; nsam = 10; dim = 2")
142 ndim = 2; nsam = 10; dim = 2
143 call disp%show("sample = reshape(getUnifRand(1, 20, ndim * nsam), shape = [ndim, nsam], order = [2, 1])")
144 sample = reshape(getUnifRand(1, 20, ndim * nsam), shape = [ndim, nsam], order = [2, 1])
145 call disp%show("sample")
146 call disp%show( sample , format = format )
147 call disp%show("call setResized(mean, ndim)")
148 call setResized(mean, ndim)
149 call disp%show("iweight = getUnifRand(1, 10, nsam) ! integer-valued weights.")
150 iweight = getUnifRand(1, 10, nsam) ! integer-valued weights.
151 call disp%show("iweight")
152 call disp%show( iweight )
153 call disp%show("call setMean(mean, sample, dim, iweight, iweisum)")
154 call setMean(mean, sample, dim, iweight, iweisum)
155 call disp%show("mean")
156 call disp%show( mean )
157 call disp%show("iweisum")
158 call disp%show( iweisum )
159 call disp%show("rweight = iweight ! or real-valued weights.")
160 rweight = iweight ! or real-valued weights.
161 call disp%show("iweight")
162 call disp%show( iweight )
163 call disp%show("call setMean(mean, sample, dim, rweight, rweisum)")
164 call setMean(mean, sample, dim, rweight, rweisum)
165 call disp%show("mean")
166 call disp%show( mean , format = format )
167 call disp%show("rweisum")
168 call disp%show( rweisum , format = format )
169
170 call disp%skip()
171 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
172 call disp%show("!Compute the correlation matrix with integer weights.")
173 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
174 call disp%skip()
175
176 call disp%show("cor = getCor(sample, dim, iweight) ! same result as above.")
177 cor = getCor(sample, dim, iweight)
178 call disp%show("cor")
179 call disp%show( cor , format = format )
180 call disp%skip()
181 call disp%show("Compute the sample correlation along the first dimension.", deliml = SK_'''')
182 call disp%skip()
183 call disp%show("dim = 1")
184 dim = 1
185 call disp%show("cor = getCor(transpose(sample), dim, iweight) ! same result as above.")
186 cor = getCor(transpose(sample), dim, iweight)
187 call disp%show("cor")
188 call disp%show( cor , format = format )
189 call disp%skip()
190 call disp%show("Compute the full sample correlation for a pair of time series.", deliml = SK_'''')
191 call disp%skip()
192 call disp%show("cor(1,1) = getCor(sample(1,:), sample(2,:), iweight) ! same result as above.")
193 cor(1,1) = getCor(sample(1,:), sample(2,:), iweight)
194 call disp%show("cor(1,1)")
195 call disp%show( cor(1,1) , format = format )
196 call disp%skip()
197
198 call disp%skip()
199 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
200 call disp%show("!Compute the correlation matrix with real weights.")
201 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
202 call disp%skip()
203
204 call disp%show("dim = 2")
205 dim = 2
206 call disp%show("cor = getCor(sample, dim, rweight) ! same result as above.")
207 cor = getCor(sample, dim, rweight)
208 call disp%show("cor")
209 call disp%show( cor , format = format )
210 call disp%skip()
211 call disp%show("Compute the sample correlation along the first dimension.", deliml = SK_'''')
212 call disp%skip()
213 call disp%show("dim = 1")
214 dim = 1
215 call disp%show("cor = getCor(transpose(sample), dim, rweight) ! same result as above.")
216 cor = getCor(transpose(sample), dim, rweight)
217 call disp%show("cor")
218 call disp%show( cor , format = format )
219 call disp%skip()
220 call disp%show("Compute the full sample correlation for a pair of time series.", deliml = SK_'''')
221 call disp%skip()
222 call disp%show("cor(1,1) = getCor(sample(1,:), sample(2,:), rweight) ! same result as above.")
223 cor(1,1) = getCor(sample(1,:), sample(2,:), rweight)
224 call disp%show("cor(1,1)")
225 call disp%show( cor(1,1) , format = format )
226 call disp%skip()
227 end block
228
229end program example
Generate and return an array of the specified rank and shape of arbitrary intrinsic type and kind wit...
Generate minimally-spaced character, integer, real sequences or sequences at fixed intervals of size ...
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 an equally-weighted (verbose or flattened) array of the input weighted array of rank 1 or 2.
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
Generate and return the (optionally unbiased) covariance matrix of a pair of (potentially weighted) t...
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 generating ranges of discrete character,...
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 procedures and generic interfaces for flattening (duplicating the elements of) a...
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 classes and procedures for computing the properties related to the covariance ma...
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!Compute the correlation matrix from covariance matrix.
4!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
5
6
7ndim = getUnifRand(0, 7)
8ndim
9+7
10std = getUnifRand(1, ndim, ndim)
11std
12+4.00000000, +1.00000000, +7.00000000, +1.00000000, +2.00000000, +6.00000000, +1.00000000
13cov = getCovRand(1._TKG, std)
14cov
15+15.9999981, -3.62055397, +21.8242397, -0.647497892, +2.41851354, -8.25101376, -1.86213255
16-3.62055397, +1.00000000, -6.40514708, +0.185297295, -0.339333564, +0.220663548E-1, +0.274384558
17+21.8242397, -6.40514708, +48.9999962, +0.403621674, -1.22789717, +8.26337337, -1.21427906
18-0.647497892, +0.185297295, +0.403621674, +0.999999821, +0.265933573, -0.912627339, +0.351888090
19+2.41851354, -0.339333564, -1.22789717, +0.265933573, +4.00000000, -9.75441742, -0.300967097
20-8.25101376, +0.220663548E-1, +8.26337337, -0.912627339, -9.75441742, +36.0000000, +2.71584511
21-1.86213255, +0.274384558, -1.21427906, +0.351888090, -0.300967097, +2.71584511, +1.00000012
22cor = getCor(cov, uppDia)
23cor
24+1.00000000, -0.905138612, +0.779437244, -0.161874518, +0.302314222, -0.343792289, -0.465533197
25-0.905138612, +1.00000000, -0.915021062, +0.185297310, -0.169666782, +0.367772579E-2, +0.274384558
26+0.779437244, -0.915021062, +1.00000000, +0.576602481E-1, -0.877069458E-1, +0.196747005, -0.173468441
27-0.161874518, +0.185297310, +0.576602481E-1, +1.00000000, +0.132966802, -0.152104571, +0.351888120
28+0.302314222, -0.169666782, -0.877069458E-1, +0.132966802, +1.00000000, -0.812868118, -0.150483549
29-0.343792289, +0.367772579E-2, +0.196747005, -0.152104571, -0.812868118, +1.00000000, +0.452640861
30-0.465533197, +0.274384558, -0.173468441, +0.351888120, -0.150483549, +0.452640861, +1.00000000
31cor = getCor(cov, upp, stdinv = 1 / std)
32cor
33+1.00000000, -0.905138493, +0.779437184, -0.161874473, +0.302314192, -0.343792260, -0.465533137
34-0.905138493, +1.00000000, -0.915021062, +0.185297295, -0.169666782, +0.367772579E-2, +0.274384558
35+0.779437184, -0.915021062, +1.00000000, +0.576602407E-1, -0.877069458E-1, +0.196747005, -0.173468441
36-0.161874473, +0.185297295, +0.576602407E-1, +1.00000000, +0.132966787, -0.152104557, +0.351888090
37+0.302314192, -0.169666782, -0.877069458E-1, +0.132966787, +1.00000000, -0.812868118, -0.150483549
38-0.343792260, +0.367772579E-2, +0.196747005, -0.152104557, -0.812868118, +1.00000000, +0.452640861
39-0.465533137, +0.274384558, -0.173468441, +0.351888090, -0.150483549, +0.452640861, +1.00000000
40cor = getCor(cov, lowDia)
41cor
42+1.00000000, -0.905138612, +0.779437244, -0.161874518, +0.302314222, -0.343792289, -0.465533197
43-0.905138612, +1.00000000, -0.915021062, +0.185297310, -0.169666782, +0.367772579E-2, +0.274384558
44+0.779437244, -0.915021062, +1.00000000, +0.576602481E-1, -0.877069458E-1, +0.196747005, -0.173468441
45-0.161874518, +0.185297310, +0.576602481E-1, +1.00000000, +0.132966802, -0.152104571, +0.351888120
46+0.302314222, -0.169666782, -0.877069458E-1, +0.132966802, +1.00000000, -0.812868118, -0.150483549
47-0.343792289, +0.367772579E-2, +0.196747005, -0.152104571, -0.812868118, +1.00000000, +0.452640861
48-0.465533197, +0.274384558, -0.173468441, +0.351888120, -0.150483549, +0.452640861, +1.00000000
49cor = getCor(cov, low, stdinv = 1 / std)
50cor
51+1.00000000, -0.905138493, +0.779437184, -0.161874473, +0.302314192, -0.343792260, -0.465533137
52-0.905138493, +1.00000000, -0.915021062, +0.185297295, -0.169666782, +0.367772579E-2, +0.274384558
53+0.779437184, -0.915021062, +1.00000000, +0.576602407E-1, -0.877069458E-1, +0.196747005, -0.173468441
54-0.161874473, +0.185297295, +0.576602407E-1, +1.00000000, +0.132966787, -0.152104557, +0.351888090
55+0.302314192, -0.169666782, -0.877069458E-1, +0.132966787, +1.00000000, -0.812868118, -0.150483549
56-0.343792260, +0.367772579E-2, +0.196747005, -0.152104557, -0.812868118, +1.00000000, +0.452640861
57-0.465533137, +0.274384558, -0.173468441, +0.351888090, -0.150483549, +0.452640861, +1.00000000
58getCov(getCor(cov, lowDia), lowDia, std) ! reconstruct the original covariance matrix.
59+16.0000000, -3.62055445, +21.8242435, -0.647498071, +2.41851377, -8.25101471, -1.86213279
60-3.62055445, +1.00000000, -6.40514755, +0.185297310, -0.339333564, +0.220663548E-1, +0.274384558
61+21.8242435, -6.40514755, +49.0000000, +0.403621733, -1.22789729, +8.26337433, -1.21427906
62-0.647498071, +0.185297310, +0.403621733, +1.00000000, +0.265933603, -0.912627459, +0.351888120
63+2.41851377, -0.339333564, -1.22789729, +0.265933603, +4.00000000, -9.75441742, -0.300967097
64-8.25101471, +0.220663548E-1, +8.26337433, -0.912627459, -9.75441742, +36.0000000, +2.71584511
65-1.86213279, +0.274384558, -1.21427906, +0.351888120, -0.300967097, +2.71584511, +1.00000000
66getCov(getCor(cov, uppDia), uppDia, std) ! reconstruct the original covariance matrix.
67+16.0000000, -3.62055445, +21.8242435, -0.647498071, +2.41851377, -8.25101471, -1.86213279
68-3.62055445, +1.00000000, -6.40514755, +0.185297310, -0.339333564, +0.220663548E-1, +0.274384558
69+21.8242435, -6.40514755, +49.0000000, +0.403621733, -1.22789729, +8.26337433, -1.21427906
70-0.647498071, +0.185297310, +0.403621733, +1.00000000, +0.265933603, -0.912627459, +0.351888120
71+2.41851377, -0.339333564, -1.22789729, +0.265933603, +4.00000000, -9.75441742, -0.300967097
72-8.25101471, +0.220663548E-1, +8.26337433, -0.912627459, -9.75441742, +36.0000000, +2.71584511
73-1.86213279, +0.274384558, -1.21427906, +0.351888120, -0.300967097, +2.71584511, +1.00000000
74
75
76ndim = getUnifRand(0, 7)
77ndim
78+2
79std = getUnifRand(1, ndim, ndim)
80std
81+2.00000000, +2.00000000
82cov = getCovRand(1._TKG, std)
83cov
84+4.00000000, +3.99492455
85+3.99492455, +4.00000048
86cor = getCor(cov, uppDia)
87cor
88+1.00000000, +0.998731136
89+0.998731136, +1.00000000
90cor = getCor(cov, upp, stdinv = 1 / std)
91cor
92+1.00000000, +0.998731136
93+0.998731136, +1.00000000
94cor = getCor(cov, lowDia)
95cor
96+1.00000000, +0.998731136
97+0.998731136, +1.00000000
98cor = getCor(cov, low, stdinv = 1 / std)
99cor
100+1.00000000, +0.998731136
101+0.998731136, +1.00000000
102getCov(getCor(cov, lowDia), lowDia, std) ! reconstruct the original covariance matrix.
103+4.00000000, +3.99492455
104+3.99492455, +4.00000000
105getCov(getCor(cov, uppDia), uppDia, std) ! reconstruct the original covariance matrix.
106+4.00000000, +3.99492455
107+3.99492455, +4.00000000
108
109
110ndim = getUnifRand(0, 7)
111ndim
112+6
113std = getUnifRand(1, ndim, ndim)
114std
115+2.00000000, +1.00000000, +3.00000000, +4.00000000, +2.00000000, +1.00000000
116cov = getCovRand(1._TKG, std)
117cov
118+3.99999952, +1.18026459, +4.16013670, +0.893203139, +1.67780912, -0.483220547
119+1.18026459, +1.00000000, -0.513817668, -1.36749721, -0.634094179E-1, +0.392314851
120+4.16013670, -0.513817668, +8.99999809, +5.69412661, +3.14331579, -2.00817108
121+0.893203139, -1.36749721, +5.69412661, +16.0000000, -2.30509329, -2.83024764
122+1.67780912, -0.634094179E-1, +3.14331579, -2.30509329, +4.00000000, -0.859302998
123-0.483220547, +0.392314851, -2.00817108, -2.83024764, -0.859302998, +0.999999940
124cor = getCor(cov, uppDia)
125cor
126+1.00000000, +0.590132356, +0.693356276, +0.111650407, +0.419452339, -0.241610333
127+0.590132356, +1.00000000, -0.171272576, -0.341874301, -0.317047089E-1, +0.392314911
128+0.693356276, -0.171272576, +1.00000000, +0.474510610, +0.523886025, -0.669390500
129+0.111650407, -0.341874301, +0.474510610, +1.00000000, -0.288136661, -0.707561970
130+0.419452339, -0.317047089E-1, +0.523886025, -0.288136661, +1.00000000, -0.429651558
131-0.241610333, +0.392314911, -0.669390500, -0.707561970, -0.429651558, +1.00000000
132cor = getCor(cov, upp, stdinv = 1 / std)
133cor
134+1.00000000, +0.590132296, +0.693356156, +0.111650392, +0.419452280, -0.241610274
135+0.590132296, +1.00000000, -0.171272561, -0.341874301, -0.317047089E-1, +0.392314851
136+0.693356156, -0.171272561, +1.00000000, +0.474510550, +0.523885965, -0.669390380
137+0.111650392, -0.341874301, +0.474510550, +1.00000000, -0.288136661, -0.707561910
138+0.419452280, -0.317047089E-1, +0.523885965, -0.288136661, +1.00000000, -0.429651499
139-0.241610274, +0.392314851, -0.669390380, -0.707561910, -0.429651499, +1.00000000
140cor = getCor(cov, lowDia)
141cor
142+1.00000000, +0.590132356, +0.693356276, +0.111650407, +0.419452339, -0.241610333
143+0.590132356, +1.00000000, -0.171272576, -0.341874301, -0.317047089E-1, +0.392314911
144+0.693356276, -0.171272576, +1.00000000, +0.474510610, +0.523886025, -0.669390500
145+0.111650407, -0.341874301, +0.474510610, +1.00000000, -0.288136661, -0.707561970
146+0.419452339, -0.317047089E-1, +0.523886025, -0.288136661, +1.00000000, -0.429651558
147-0.241610333, +0.392314911, -0.669390500, -0.707561970, -0.429651558, +1.00000000
148cor = getCor(cov, low, stdinv = 1 / std)
149cor
150+1.00000000, +0.590132296, +0.693356156, +0.111650392, +0.419452280, -0.241610274
151+0.590132296, +1.00000000, -0.171272561, -0.341874301, -0.317047089E-1, +0.392314851
152+0.693356156, -0.171272561, +1.00000000, +0.474510550, +0.523885965, -0.669390380
153+0.111650392, -0.341874301, +0.474510550, +1.00000000, -0.288136661, -0.707561910
154+0.419452280, -0.317047089E-1, +0.523885965, -0.288136661, +1.00000000, -0.429651499
155-0.241610274, +0.392314851, -0.669390380, -0.707561910, -0.429651499, +1.00000000
156getCov(getCor(cov, lowDia), lowDia, std) ! reconstruct the original covariance matrix.
157+4.00000000, +1.18026471, +4.16013765, +0.893203259, +1.67780936, -0.483220667
158+1.18026471, +1.00000000, -0.513817728, -1.36749721, -0.634094179E-1, +0.392314911
159+4.16013765, -0.513817728, +9.00000000, +5.69412708, +3.14331627, -2.00817156
160+0.893203259, -1.36749721, +5.69412708, +16.0000000, -2.30509329, -2.83024788
161+1.67780936, -0.634094179E-1, +3.14331627, -2.30509329, +4.00000000, -0.859303117
162-0.483220667, +0.392314911, -2.00817156, -2.83024788, -0.859303117, +1.00000000
163getCov(getCor(cov, uppDia), uppDia, std) ! reconstruct the original covariance matrix.
164+4.00000000, +1.18026471, +4.16013765, +0.893203259, +1.67780936, -0.483220667
165+1.18026471, +1.00000000, -0.513817728, -1.36749721, -0.634094179E-1, +0.392314911
166+4.16013765, -0.513817728, +9.00000000, +5.69412708, +3.14331627, -2.00817156
167+0.893203259, -1.36749721, +5.69412708, +16.0000000, -2.30509329, -2.83024788
168+1.67780936, -0.634094179E-1, +3.14331627, -2.30509329, +4.00000000, -0.859303117
169-0.483220667, +0.392314911, -2.00817156, -2.83024788, -0.859303117, +1.00000000
170
171
172ndim = getUnifRand(0, 7)
173ndim
174+4
175std = getUnifRand(1, ndim, ndim)
176std
177+4.00000000, +2.00000000, +2.00000000, +1.00000000
178cov = getCovRand(1._TKG, std)
179cov
180+16.0000000, -5.58781815, +2.32739043, +0.626870215
181-5.58781815, +4.00000000, -1.44074559, +0.461182684
182+2.32739043, -1.44074559, +3.99999952, +1.27536035
183+0.626870215, +0.461182684, +1.27536035, +1.00000000
184cor = getCor(cov, uppDia)
185cor
186+1.00000000, -0.698477268, +0.290923834, +0.156717554
187-0.698477268, +1.00000000, -0.360186428, +0.230591342
188+0.290923834, -0.360186428, +1.00000000, +0.637680233
189+0.156717554, +0.230591342, +0.637680233, +1.00000000
190cor = getCor(cov, upp, stdinv = 1 / std)
191cor
192+1.00000000, -0.698477268, +0.290923804, +0.156717554
193-0.698477268, +1.00000000, -0.360186398, +0.230591342
194+0.290923804, -0.360186398, +1.00000000, +0.637680173
195+0.156717554, +0.230591342, +0.637680173, +1.00000000
196cor = getCor(cov, lowDia)
197cor
198+1.00000000, -0.698477268, +0.290923834, +0.156717554
199-0.698477268, +1.00000000, -0.360186428, +0.230591342
200+0.290923834, -0.360186428, +1.00000000, +0.637680233
201+0.156717554, +0.230591342, +0.637680233, +1.00000000
202cor = getCor(cov, low, stdinv = 1 / std)
203cor
204+1.00000000, -0.698477268, +0.290923804, +0.156717554
205-0.698477268, +1.00000000, -0.360186398, +0.230591342
206+0.290923804, -0.360186398, +1.00000000, +0.637680173
207+0.156717554, +0.230591342, +0.637680173, +1.00000000
208getCov(getCor(cov, lowDia), lowDia, std) ! reconstruct the original covariance matrix.
209+16.0000000, -5.58781815, +2.32739067, +0.626870215
210-5.58781815, +4.00000000, -1.44074571, +0.461182684
211+2.32739067, -1.44074571, +4.00000000, +1.27536047
212+0.626870215, +0.461182684, +1.27536047, +1.00000000
213getCov(getCor(cov, uppDia), uppDia, std) ! reconstruct the original covariance matrix.
214+16.0000000, -5.58781815, +2.32739067, +0.626870215
215-5.58781815, +4.00000000, -1.44074571, +0.461182684
216+2.32739067, -1.44074571, +4.00000000, +1.27536047
217+0.626870215, +0.461182684, +1.27536047, +1.00000000
218
219
220ndim = getUnifRand(0, 7)
221ndim
222+5
223std = getUnifRand(1, ndim, ndim)
224std
225+4.00000000, +1.00000000, +3.00000000, +4.00000000, +2.00000000
226cov = getCovRand(1._TKG, std)
227cov
228+16.0000000, +3.02767801, -9.61742592, -12.0814743, -5.63244200
229+3.02767801, +1.00000000, -1.70646954, -3.79589796, -0.428531528
230-9.61742592, -1.70646954, +9.00000191, +5.00026083, +2.13427949
231-12.0814743, -3.79589796, +5.00026083, +16.0000019, +3.10046911
232-5.63244200, -0.428531528, +2.13427949, +3.10046911, +4.00000000
233cor = getCor(cov, uppDia)
234cor
235+1.00000000, +0.756919503, -0.801452100, -0.755092144, -0.704055250
236+0.756919503, +1.00000000, -0.568823159, -0.948974490, -0.214265764
237-0.801452100, -0.568823159, +1.00000000, +0.416688383, +0.355713218
238-0.755092144, -0.948974490, +0.416688383, +1.00000000, +0.387558639
239-0.704055250, -0.214265764, +0.355713218, +0.387558639, +1.00000000
240cor = getCor(cov, upp, stdinv = 1 / std)
241cor
242+1.00000000, +0.756919503, -0.801452160, -0.755092144, -0.704055250
243+0.756919503, +1.00000000, -0.568823218, -0.948974490, -0.214265764
244-0.801452160, -0.568823218, +1.00000000, +0.416688412, +0.355713248
245-0.755092144, -0.948974490, +0.416688412, +1.00000000, +0.387558639
246-0.704055250, -0.214265764, +0.355713248, +0.387558639, +1.00000000
247cor = getCor(cov, lowDia)
248cor
249+1.00000000, +0.756919503, -0.801452100, -0.755092144, -0.704055250
250+0.756919503, +1.00000000, -0.568823159, -0.948974490, -0.214265764
251-0.801452100, -0.568823159, +1.00000000, +0.416688383, +0.355713218
252-0.755092144, -0.948974490, +0.416688383, +1.00000000, +0.387558639
253-0.704055250, -0.214265764, +0.355713218, +0.387558639, +1.00000000
254cor = getCor(cov, low, stdinv = 1 / std)
255cor
256+1.00000000, +0.756919503, -0.801452160, -0.755092144, -0.704055250
257+0.756919503, +1.00000000, -0.568823218, -0.948974490, -0.214265764
258-0.801452160, -0.568823218, +1.00000000, +0.416688412, +0.355713248
259-0.755092144, -0.948974490, +0.416688412, +1.00000000, +0.387558639
260-0.704055250, -0.214265764, +0.355713248, +0.387558639, +1.00000000
261getCov(getCor(cov, lowDia), lowDia, std) ! reconstruct the original covariance matrix.
262+16.0000000, +3.02767801, -9.61742496, -12.0814743, -5.63244200
263+3.02767801, +1.00000000, -1.70646954, -3.79589796, -0.428531528
264-9.61742496, -1.70646954, +9.00000000, +5.00026035, +2.13427925
265-12.0814743, -3.79589796, +5.00026035, +16.0000000, +3.10046911
266-5.63244200, -0.428531528, +2.13427925, +3.10046911, +4.00000000
267getCov(getCor(cov, uppDia), uppDia, std) ! reconstruct the original covariance matrix.
268+16.0000000, +3.02767801, -9.61742496, -12.0814743, -5.63244200
269+3.02767801, +1.00000000, -1.70646954, -3.79589796, -0.428531528
270-9.61742496, -1.70646954, +9.00000000, +5.00026035, +2.13427925
271-12.0814743, -3.79589796, +5.00026035, +16.0000000, +3.10046911
272-5.63244200, -0.428531528, +2.13427925, +3.10046911, +4.00000000
273
274
275ndim = getUnifRand(0, 7)
276ndim
277+5
278std = getUnifRand(1, ndim, ndim)
279std
280+2.00000000, +3.00000000, +2.00000000, +2.00000000, +1.00000000
281cov = getCovRand(1._TKG, std)
282cov
283+4.00000000, +2.49345994, +3.08140135, +2.95499659, -0.922998667
284+2.49345994, +9.00000000, -1.55313754, +2.68888330, -1.09782410
285+3.08140135, -1.55313754, +4.00000000, +1.89834201, -0.490174472
286+2.95499659, +2.68888330, +1.89834201, +4.00000000, -1.57213426
287-0.922998667, -1.09782410, -0.490174472, -1.57213426, +1.00000000
288cor = getCor(cov, uppDia)
289cor
290+1.00000000, +0.415576667, +0.770350337, +0.738749146, -0.461499333
291+0.415576667, +1.00000000, -0.258856267, +0.448147237, -0.365941375
292+0.770350337, -0.258856267, +1.00000000, +0.474585503, -0.245087236
293+0.738749146, +0.448147237, +0.474585503, +1.00000000, -0.786067128
294-0.461499333, -0.365941375, -0.245087236, -0.786067128, +1.00000000
295cor = getCor(cov, upp, stdinv = 1 / std)
296cor
297+1.00000000, +0.415576667, +0.770350337, +0.738749146, -0.461499333
298+0.415576667, +1.00000000, -0.258856267, +0.448147237, -0.365941375
299+0.770350337, -0.258856267, +1.00000000, +0.474585503, -0.245087236
300+0.738749146, +0.448147237, +0.474585503, +1.00000000, -0.786067128
301-0.461499333, -0.365941375, -0.245087236, -0.786067128, +1.00000000
302cor = getCor(cov, lowDia)
303cor
304+1.00000000, +0.415576667, +0.770350337, +0.738749146, -0.461499333
305+0.415576667, +1.00000000, -0.258856267, +0.448147237, -0.365941375
306+0.770350337, -0.258856267, +1.00000000, +0.474585503, -0.245087236
307+0.738749146, +0.448147237, +0.474585503, +1.00000000, -0.786067128
308-0.461499333, -0.365941375, -0.245087236, -0.786067128, +1.00000000
309cor = getCor(cov, low, stdinv = 1 / std)
310cor
311+1.00000000, +0.415576667, +0.770350337, +0.738749146, -0.461499333
312+0.415576667, +1.00000000, -0.258856267, +0.448147237, -0.365941375
313+0.770350337, -0.258856267, +1.00000000, +0.474585503, -0.245087236
314+0.738749146, +0.448147237, +0.474585503, +1.00000000, -0.786067128
315-0.461499333, -0.365941375, -0.245087236, -0.786067128, +1.00000000
316getCov(getCor(cov, lowDia), lowDia, std) ! reconstruct the original covariance matrix.
317+4.00000000, +2.49345994, +3.08140135, +2.95499659, -0.922998667
318+2.49345994, +9.00000000, -1.55313754, +2.68888330, -1.09782410
319+3.08140135, -1.55313754, +4.00000000, +1.89834201, -0.490174472
320+2.95499659, +2.68888330, +1.89834201, +4.00000000, -1.57213426
321-0.922998667, -1.09782410, -0.490174472, -1.57213426, +1.00000000
322getCov(getCor(cov, uppDia), uppDia, std) ! reconstruct the original covariance matrix.
323+4.00000000, +2.49345994, +3.08140135, +2.95499659, -0.922998667
324+2.49345994, +9.00000000, -1.55313754, +2.68888330, -1.09782410
325+3.08140135, -1.55313754, +4.00000000, +1.89834201, -0.490174472
326+2.95499659, +2.68888330, +1.89834201, +4.00000000, -1.57213426
327-0.922998667, -1.09782410, -0.490174472, -1.57213426, +1.00000000
328
329
330ndim = getUnifRand(0, 7)
331ndim
332+2
333std = getUnifRand(1, ndim, ndim)
334std
335+2.00000000, +1.00000000
336cov = getCovRand(1._TKG, std)
337cov
338+4.00000000, -0.182149988E-1
339-0.182149988E-1, +1.00000000
340cor = getCor(cov, uppDia)
341cor
342+1.00000000, -0.910749938E-2
343-0.910749938E-2, +1.00000000
344cor = getCor(cov, upp, stdinv = 1 / std)
345cor
346+1.00000000, -0.910749938E-2
347-0.910749938E-2, +1.00000000
348cor = getCor(cov, lowDia)
349cor
350+1.00000000, -0.910749938E-2
351-0.910749938E-2, +1.00000000
352cor = getCor(cov, low, stdinv = 1 / std)
353cor
354+1.00000000, -0.910749938E-2
355-0.910749938E-2, +1.00000000
356getCov(getCor(cov, lowDia), lowDia, std) ! reconstruct the original covariance matrix.
357+4.00000000, -0.182149988E-1
358-0.182149988E-1, +1.00000000
359getCov(getCor(cov, uppDia), uppDia, std) ! reconstruct the original covariance matrix.
360+4.00000000, -0.182149988E-1
361-0.182149988E-1, +1.00000000
362
363
364ndim = getUnifRand(0, 7)
365ndim
366+5
367std = getUnifRand(1, ndim, ndim)
368std
369+3.00000000, +4.00000000, +4.00000000, +1.00000000, +1.00000000
370cov = getCovRand(1._TKG, std)
371cov
372+9.00000000, +2.55421615, +8.16738319, -1.29766917, -0.384563580E-1
373+2.55421615, +16.0000019, +4.42552280, -0.704014182, -1.86062086
374+8.16738319, +4.42552280, +15.9999990, -0.147946119, -1.76953816
375-1.29766917, -0.704014182, -0.147946119, +1.00000024, -0.721210301
376-0.384563580E-1, -1.86062086, -1.76953816, -0.721210301, +0.999999940
377cor = getCor(cov, uppDia)
378cor
379+1.00000000, +0.212851346, +0.680615366, -0.432556361, -0.128187872E-1
380+0.212851346, +1.00000000, +0.276595205, -0.176003531, -0.465155274
381+0.680615366, +0.276595205, +1.00000000, -0.369865298E-1, -0.442384660
382-0.432556361, -0.176003531, -0.369865298E-1, +1.00000000, -0.721210301
383-0.128187872E-1, -0.465155274, -0.442384660, -0.721210301, +1.00000000
384cor = getCor(cov, upp, stdinv = 1 / std)
385cor
386+1.00000000, +0.212851346, +0.680615306, -0.432556391, -0.128187863E-1
387+0.212851346, +1.00000000, +0.276595175, -0.176003546, -0.465155214
388+0.680615306, +0.276595175, +1.00000000, -0.369865298E-1, -0.442384541
389-0.432556391, -0.176003546, -0.369865298E-1, +1.00000000, -0.721210301
390-0.128187863E-1, -0.465155214, -0.442384541, -0.721210301, +1.00000000
391cor = getCor(cov, lowDia)
392cor
393+1.00000000, +0.212851346, +0.680615366, -0.432556361, -0.128187872E-1
394+0.212851346, +1.00000000, +0.276595205, -0.176003531, -0.465155274
395+0.680615366, +0.276595205, +1.00000000, -0.369865298E-1, -0.442384660
396-0.432556361, -0.176003531, -0.369865298E-1, +1.00000000, -0.721210301
397-0.128187872E-1, -0.465155274, -0.442384660, -0.721210301, +1.00000000
398cor = getCor(cov, low, stdinv = 1 / std)
399cor
400+1.00000000, +0.212851346, +0.680615306, -0.432556391, -0.128187863E-1
401+0.212851346, +1.00000000, +0.276595175, -0.176003546, -0.465155214
402+0.680615306, +0.276595175, +1.00000000, -0.369865298E-1, -0.442384541
403-0.432556391, -0.176003546, -0.369865298E-1, +1.00000000, -0.721210301
404-0.128187863E-1, -0.465155214, -0.442384541, -0.721210301, +1.00000000
405getCov(getCor(cov, lowDia), lowDia, std) ! reconstruct the original covariance matrix.
406+9.00000000, +2.55421615, +8.16738415, -1.29766905, -0.384563617E-1
407+2.55421615, +16.0000000, +4.42552328, -0.704014122, -1.86062109
408+8.16738415, +4.42552328, +16.0000000, -0.147946119, -1.76953864
409-1.29766905, -0.704014122, -0.147946119, +1.00000000, -0.721210301
410-0.384563617E-1, -1.86062109, -1.76953864, -0.721210301, +1.00000000
411getCov(getCor(cov, uppDia), uppDia, std) ! reconstruct the original covariance matrix.
412+9.00000000, +2.55421615, +8.16738415, -1.29766905, -0.384563617E-1
413+2.55421615, +16.0000000, +4.42552328, -0.704014122, -1.86062109
414+8.16738415, +4.42552328, +16.0000000, -0.147946119, -1.76953864
415-1.29766905, -0.704014122, -0.147946119, +1.00000000, -0.721210301
416-0.384563617E-1, -1.86062109, -1.76953864, -0.721210301, +1.00000000
417
418
419ndim = getUnifRand(0, 7)
420ndim
421+7
422std = getUnifRand(1, ndim, ndim)
423std
424+1.00000000, +3.00000000, +2.00000000, +4.00000000, +2.00000000, +3.00000000, +2.00000000
425cov = getCovRand(1._TKG, std)
426cov
427+1.00000000, -2.99995232, +0.928574502, -2.10579324, +1.14359212, -1.54214978, -0.607439354E-1
428-2.99995232, +9.00000000, -2.80412865, +6.28827095, -3.41697836, +4.60675716, +0.203933954
429+0.928574502, -2.80412865, +4.00000048, -0.193841815, +0.973897874, +1.65067041, -1.39805067
430-2.10579324, +6.28827095, -0.193841815, +15.9999990, -0.858241320, +4.37173319, -1.73380888
431+1.14359212, -3.41697836, +0.973897874, -0.858241320, +3.99999976, -1.01856732, +1.48593700
432-1.54214978, +4.60675716, +1.65067041, +4.37173319, -1.01856732, +8.99999905, -1.60488367
433-0.607439354E-1, +0.203933954, -1.39805067, -1.73380888, +1.48593700, -1.60488367, +4.00000000
434cor = getCor(cov, uppDia)
435cor
436+1.00000000, -0.999984145, +0.464287251, -0.526448369, +0.571796119, -0.514050007, -0.303719677E-1
437-0.999984145, +1.00000000, -0.467354774, +0.524022639, -0.569496453, +0.511861980, +0.339889936E-1
438+0.464287251, -0.467354774, +1.00000000, -0.242302306E-1, +0.243474498, +0.275111765, -0.349512666
439-0.526448369, +0.524022639, -0.242302306E-1, +1.00000000, -0.107280187, +0.364311188, -0.216726139
440+0.571796119, -0.569496453, +0.243474498, -0.107280187, +1.00000000, -0.169761255, +0.371484280
441-0.514050007, +0.511861980, +0.275111765, +0.364311188, -0.169761255, +1.00000000, -0.267480642
442-0.303719677E-1, +0.339889936E-1, -0.349512666, -0.216726139, +0.371484280, -0.267480642, +1.00000000
443cor = getCor(cov, upp, stdinv = 1 / std)
444cor
445+1.00000000, -0.999984145, +0.464287251, -0.526448309, +0.571796060, -0.514049947, -0.303719677E-1
446-0.999984145, +1.00000000, -0.467354774, +0.524022579, -0.569496393, +0.511861920, +0.339889936E-1
447+0.464287251, -0.467354774, +1.00000000, -0.242302269E-1, +0.243474469, +0.275111735, -0.349512666
448-0.526448309, +0.524022579, -0.242302269E-1, +1.00000000, -0.107280165, +0.364311099, -0.216726109
449+0.571796060, -0.569496393, +0.243474469, -0.107280165, +1.00000000, -0.169761226, +0.371484250
450-0.514049947, +0.511861920, +0.275111735, +0.364311099, -0.169761226, +1.00000000, -0.267480612
451-0.303719677E-1, +0.339889936E-1, -0.349512666, -0.216726109, +0.371484250, -0.267480612, +1.00000000
452cor = getCor(cov, lowDia)
453cor
454+1.00000000, -0.999984145, +0.464287251, -0.526448369, +0.571796119, -0.514050007, -0.303719677E-1
455-0.999984145, +1.00000000, -0.467354774, +0.524022639, -0.569496453, +0.511861980, +0.339889936E-1
456+0.464287251, -0.467354774, +1.00000000, -0.242302306E-1, +0.243474498, +0.275111765, -0.349512666
457-0.526448369, +0.524022639, -0.242302306E-1, +1.00000000, -0.107280187, +0.364311188, -0.216726139
458+0.571796119, -0.569496453, +0.243474498, -0.107280187, +1.00000000, -0.169761255, +0.371484280
459-0.514050007, +0.511861980, +0.275111765, +0.364311188, -0.169761255, +1.00000000, -0.267480642
460-0.303719677E-1, +0.339889936E-1, -0.349512666, -0.216726139, +0.371484280, -0.267480642, +1.00000000
461cor = getCor(cov, low, stdinv = 1 / std)
462cor
463+1.00000000, -0.999984145, +0.464287251, -0.526448309, +0.571796060, -0.514049947, -0.303719677E-1
464-0.999984145, +1.00000000, -0.467354774, +0.524022579, -0.569496393, +0.511861920, +0.339889936E-1
465+0.464287251, -0.467354774, +1.00000000, -0.242302269E-1, +0.243474469, +0.275111735, -0.349512666
466-0.526448309, +0.524022579, -0.242302269E-1, +1.00000000, -0.107280165, +0.364311099, -0.216726109
467+0.571796060, -0.569496393, +0.243474469, -0.107280165, +1.00000000, -0.169761226, +0.371484250
468-0.514049947, +0.511861920, +0.275111735, +0.364311099, -0.169761226, +1.00000000, -0.267480612
469-0.303719677E-1, +0.339889936E-1, -0.349512666, -0.216726109, +0.371484250, -0.267480612, +1.00000000
470getCov(getCor(cov, lowDia), lowDia, std) ! reconstruct the original covariance matrix.
471+1.00000000, -2.99995232, +0.928574502, -2.10579348, +1.14359224, -1.54215002, -0.607439354E-1
472-2.99995232, +9.00000000, -2.80412865, +6.28827190, -3.41697884, +4.60675764, +0.203933954
473+0.928574502, -2.80412865, +4.00000000, -0.193841845, +0.973897994, +1.65067053, -1.39805067
474-2.10579348, +6.28827190, -0.193841845, +16.0000000, -0.858241498, +4.37173414, -1.73380911
475+1.14359224, -3.41697884, +0.973897994, -0.858241498, +4.00000000, -1.01856756, +1.48593712
476-1.54215002, +4.60675764, +1.65067053, +4.37173414, -1.01856756, +9.00000000, -1.60488391
477-0.607439354E-1, +0.203933954, -1.39805067, -1.73380911, +1.48593712, -1.60488391, +4.00000000
478getCov(getCor(cov, uppDia), uppDia, std) ! reconstruct the original covariance matrix.
479+1.00000000, -2.99995232, +0.928574502, -2.10579348, +1.14359224, -1.54215002, -0.607439354E-1
480-2.99995232, +9.00000000, -2.80412865, +6.28827190, -3.41697884, +4.60675764, +0.203933954
481+0.928574502, -2.80412865, +4.00000000, -0.193841845, +0.973897994, +1.65067053, -1.39805067
482-2.10579348, +6.28827190, -0.193841845, +16.0000000, -0.858241498, +4.37173414, -1.73380911
483+1.14359224, -3.41697884, +0.973897994, -0.858241498, +4.00000000, -1.01856756, +1.48593712
484-1.54215002, +4.60675764, +1.65067053, +4.37173414, -1.01856756, +9.00000000, -1.60488391
485-0.607439354E-1, +0.203933954, -1.39805067, -1.73380911, +1.48593712, -1.60488391, +4.00000000
486
487
488ndim = getUnifRand(0, 7)
489ndim
490+7
491std = getUnifRand(1, ndim, ndim)
492std
493+3.00000000, +2.00000000, +5.00000000, +7.00000000, +7.00000000, +7.00000000, +3.00000000
494cov = getCovRand(1._TKG, std)
495cov
496+9.00000000, +5.69375610, +12.5114079, -11.6907473, -11.6906328, +7.42474079, -3.75195885
497+5.69375610, +3.99999952, +9.63698769, -10.2755671, -4.53542519, +1.25932574, -2.87770605
498+12.5114079, +9.63698769, +25.0000019, -29.9349041, -3.78695035, -3.89857817, -7.74310684
499-11.6907473, -10.2755671, -29.9349041, +48.9999886, -2.16257381, +6.43660259, +13.1313171
500-11.6906328, -4.53542519, -3.78695035, -2.16257381, +48.9999886, -37.4038734, +1.04546726
501+7.42474079, +1.25932574, -3.89857817, +6.43660259, -37.4038734, +49.0000076, +1.09809399
502-3.75195885, -2.87770605, -7.74310684, +13.1313171, +1.04546726, +1.09809399, +9.00000000
503cor = getCor(cov, uppDia)
504cor
505+1.00000000, +0.948959470, +0.834093928, -0.556702375, -0.556696892, +0.353559077, -0.416884333
506+0.948959470, +1.00000000, +0.963698924, -0.733969271, -0.323959023, +0.899518430E-1, -0.479617745
507+0.834093928, +0.963698924, +1.00000000, -0.855283082, -0.108198598, -0.111387938, -0.516207159
508-0.556702375, -0.733969271, -0.855283082, +1.00000000, -0.441341698E-1, +0.131359249, +0.625300944
509-0.556696892, -0.323959023, -0.108198598, -0.441341698E-1, +1.00000000, -0.763344407, +0.497841649E-1
510+0.353559077, +0.899518430E-1, -0.111387938, +0.131359249, -0.763344407, +1.00000000, +0.522901863E-1
511-0.416884333, -0.479617745, -0.516207159, +0.625300944, +0.497841649E-1, +0.522901863E-1, +1.00000000
512cor = getCor(cov, upp, stdinv = 1 / std)
513cor
514+1.00000000, +0.948959351, +0.834093928, -0.556702316, -0.556696832, +0.353559107, -0.416884333
515+0.948959351, +1.00000000, +0.963698804, -0.733969092, -0.323958963, +0.899518430E-1, -0.479617685
516+0.834093928, +0.963698804, +1.00000000, -0.855283022, -0.108198591, -0.111387953, -0.516207159
517-0.556702316, -0.733969092, -0.855283022, +1.00000000, -0.441341624E-1, +0.131359249, +0.625300884
518-0.556696832, -0.323958963, -0.108198591, -0.441341624E-1, +1.00000000, -0.763344407, +0.497841612E-1
519+0.353559107, +0.899518430E-1, -0.111387953, +0.131359249, -0.763344407, +1.00000000, +0.522901937E-1
520-0.416884333, -0.479617685, -0.516207159, +0.625300884, +0.497841612E-1, +0.522901937E-1, +1.00000000
521cor = getCor(cov, lowDia)
522cor
523+1.00000000, +0.948959470, +0.834093928, -0.556702375, -0.556696892, +0.353559077, -0.416884333
524+0.948959470, +1.00000000, +0.963698924, -0.733969271, -0.323959023, +0.899518430E-1, -0.479617745
525+0.834093928, +0.963698924, +1.00000000, -0.855283082, -0.108198598, -0.111387938, -0.516207159
526-0.556702375, -0.733969271, -0.855283082, +1.00000000, -0.441341698E-1, +0.131359249, +0.625300944
527-0.556696892, -0.323959023, -0.108198598, -0.441341698E-1, +1.00000000, -0.763344407, +0.497841649E-1
528+0.353559077, +0.899518430E-1, -0.111387938, +0.131359249, -0.763344407, +1.00000000, +0.522901863E-1
529-0.416884333, -0.479617745, -0.516207159, +0.625300944, +0.497841649E-1, +0.522901863E-1, +1.00000000
530cor = getCor(cov, low, stdinv = 1 / std)
531cor
532+1.00000000, +0.948959351, +0.834093928, -0.556702316, -0.556696832, +0.353559107, -0.416884333
533+0.948959351, +1.00000000, +0.963698804, -0.733969092, -0.323958963, +0.899518430E-1, -0.479617685
534+0.834093928, +0.963698804, +1.00000000, -0.855283022, -0.108198591, -0.111387953, -0.516207159
535-0.556702316, -0.733969092, -0.855283022, +1.00000000, -0.441341624E-1, +0.131359249, +0.625300884
536-0.556696832, -0.323958963, -0.108198591, -0.441341624E-1, +1.00000000, -0.763344407, +0.497841612E-1
537+0.353559107, +0.899518430E-1, -0.111387953, +0.131359249, -0.763344407, +1.00000000, +0.522901937E-1
538-0.416884333, -0.479617685, -0.516207159, +0.625300884, +0.497841612E-1, +0.522901937E-1, +1.00000000
539getCov(getCor(cov, lowDia), lowDia, std) ! reconstruct the original covariance matrix.
540+9.00000000, +5.69375706, +12.5114088, -11.6907501, -11.6906347, +7.42474079, -3.75195909
541+5.69375706, +4.00000000, +9.63698959, -10.2755699, -4.53542614, +1.25932574, -2.87770653
542+12.5114088, +9.63698959, +25.0000000, -29.9349079, -3.78695083, -3.89857793, -7.74310732
543-11.6907501, -10.2755699, -29.9349079, +49.0000000, -2.16257429, +6.43660307, +13.1313200
544-11.6906347, -4.53542614, -3.78695083, -2.16257429, +49.0000000, -37.4038773, +1.04546750
545+7.42474079, +1.25932574, -3.89857793, +6.43660307, -37.4038773, +49.0000000, +1.09809387
546-3.75195909, -2.87770653, -7.74310732, +13.1313200, +1.04546750, +1.09809387, +9.00000000
547getCov(getCor(cov, uppDia), uppDia, std) ! reconstruct the original covariance matrix.
548+9.00000000, +5.69375706, +12.5114088, -11.6907501, -11.6906347, +7.42474079, -3.75195909
549+5.69375706, +4.00000000, +9.63698959, -10.2755699, -4.53542614, +1.25932574, -2.87770653
550+12.5114088, +9.63698959, +25.0000000, -29.9349079, -3.78695083, -3.89857793, -7.74310732
551-11.6907501, -10.2755699, -29.9349079, +49.0000000, -2.16257429, +6.43660307, +13.1313200
552-11.6906347, -4.53542614, -3.78695083, -2.16257429, +49.0000000, -37.4038773, +1.04546750
553+7.42474079, +1.25932574, -3.89857793, +6.43660307, -37.4038773, +49.0000000, +1.09809387
554-3.75195909, -2.87770653, -7.74310732, +13.1313200, +1.04546750, +1.09809387, +9.00000000
555
556
557!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
558!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
559!Compute the Pearson correlation matrix for a pair of time series.
560!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
561!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
562
563ndim = 2; nsam = 10; dim = 2
564sample = reshape(getUnifRand(1, 20, ndim * nsam), shape = [ndim, nsam], order = [2, 1])
565sample
566+1.100000E+01, +1.200000E+01, +5.000000E+00, +3.000000E+00, +1.100000E+01, +1.500000E+01, +2.000000E+01, +8.000000E+00, +1.600000E+01, +9.000000E+00
567+2.000000E+00, +1.000000E+01, +1.000000E+00, +1.000000E+01, +3.000000E+00, +2.000000E+01, +3.000000E+00, +1.800000E+01, +1.000000E+00, +4.000000E+00
568mean = getMean(sample, dim)
569mean
570+1.100000E+01, +7.200000E+00
571cor = getCor(sample, dim) ! same result as above.
572cor
573+1.000000E+00, -8.017607E-02
574-8.017607E-02, +1.000000E+00
575
576'Compute the sample correlation along the first dimension.'
577
578dim = 1
579cor = getCor(transpose(sample), dim) ! same result as above.
580cor
581+1.000000E+00, -8.017607E-02
582-8.017607E-02, +1.000000E+00
583
584'Compute the full sample correlation for a pair of time series.'
585
586cor(1,1) = getCor(sample(1,:), sample(2,:)) ! same result as above.
587cor(1,1)
588-8.017609E-02
589
590
591!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
592!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
593!Compute the Pearson correlation matrix for a weighted pair of time series.
594!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
595!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
596
597ndim = 2; nsam = 10; dim = 2
598sample = reshape(getUnifRand(1, 20, ndim * nsam), shape = [ndim, nsam], order = [2, 1])
599sample
600+1.400000E+01, +1.700000E+01, +2.000000E+01, +1.700000E+01, +9.000000E+00, +1.900000E+01, +1.600000E+01, +9.000000E+00, +2.000000E+01, +3.000000E+00
601+1.300000E+01, +1.000000E+00, +4.000000E+00, +3.000000E+00, +1.300000E+01, +7.000000E+00, +1.500000E+01, +1.000000E+01, +1.400000E+01, +1.600000E+01
602call setResized(mean, ndim)
603iweight = getUnifRand(1, 10, nsam) ! integer-valued weights.
604iweight
605+5, +1, +2, +8, +8, +2, +7, +4, +6, +7
606call setMean(mean, sample, dim, iweight, iweisum)
607mean
608+13.2399998, +11.1399994
609iweisum
610+50
611rweight = iweight ! or real-valued weights.
612iweight
613+5, +1, +2, +8, +8, +2, +7, +4, +6, +7
614call setMean(mean, sample, dim, rweight, rweisum)
615mean
616+1.324000E+01, +1.114000E+01
617rweisum
618+5.000000E+01
619
620!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
621!Compute the correlation matrix with integer weights.
622!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
623
624cor = getCor(sample, dim, iweight) ! same result as above.
625cor
626+1.000000E+00, -4.658200E-01
627-4.658200E-01, +1.000000E+00
628
629'Compute the sample correlation along the first dimension.'
630
631dim = 1
632cor = getCor(transpose(sample), dim, iweight) ! same result as above.
633cor
634+1.000000E+00, -4.658199E-01
635-4.658199E-01, +1.000000E+00
636
637'Compute the full sample correlation for a pair of time series.'
638
639cor(1,1) = getCor(sample(1,:), sample(2,:), iweight) ! same result as above.
640cor(1,1)
641-4.658200E-01
642
643
644!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
645!Compute the correlation matrix with real weights.
646!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
647
648dim = 2
649cor = getCor(sample, dim, rweight) ! same result as above.
650cor
651+1.000000E+00, -4.658200E-01
652-4.658200E-01, +1.000000E+00
653
654'Compute the sample correlation along the first dimension.'
655
656dim = 1
657cor = getCor(transpose(sample), dim, rweight) ! same result as above.
658cor
659+1.000000E+00, -4.658199E-01
660-4.658199E-01, +1.000000E+00
661
662'Compute the full sample correlation for a pair of time series.'
663
664cor(1,1) = getCor(sample(1,:), sample(2,:), rweight) ! same result as above.
665cor(1,1)
666-4.658200E-01
667
668
Test:
test_pm_sampleCor


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

Definition at line 765 of file pm_sampleCor.F90.


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