This example uses method Random.NextGaussianCopula
to generate a multivariate sequence GCdevt[k=0..nseq-1][j=0..nvar-1]
whose marginal distributions are user-defined and imprinted with a user-specified correlation matrix CorrMtrxIn[i=0..nvar-1][j=0..nvar-1]
and then uses method Random.CanonicalCorrelation
to extract from this multivariate random sequence a canonical correlation matrix CorrMtrx[i=0..nvar-1][j=0..nvar-1]
.
This example illustrates two useful copula related procedures. The first procedure generates a random multivariate sequence with arbitrary user-defined marginal deviates whose dependence is specified by a user-defined correlation matrix. The second procedue is the inverse of the first: an arbitary multivariate deviate input sequence is first mapped to a corresponding sequence of empirically derived variates, i.e. cumulative distribution function values representing the probability that each random variable has a value less than or equal to the input deviate. The variates are then inverted, using the inverse Normal(0,1) function, to N(0,1) deviates; and finally, a canonical covariance matrix is extracted from the multivariate N(0,1) sequence using the standard sum of products.
This example demonstrates that the nextGaussianCopula
method correctly imbeds the user-defined correlation information into an arbitrary marginal distribution sequence by extracting the canonical correlation from these sequences and showing that they differ from the original correlation matrix by a small relative error, which generally decreases as the number of multivariate sequence vectors increases.
using System; using Imsl.Math; using Imsl.Stat; public class RandomEx2 { internal static Imsl.Stat.Random IMSLRandom() { Imsl.Stat.Random r = new Imsl.Stat.Random(123457); r.Multiplier = 16807; return r; } public static void Main(string[] args) { double[,] CorrMtrxIn = new double[,]{ {1.0, - 0.9486832980505138, 0.8164965809277261}, {- 0.9486832980505138, 1.0, -0.6454972243679028}, {0.8164965809277261, -0.6454972243679028, 1.0}}; int nvar = 3; Console.WriteLine("Random Example 2:"); Console.WriteLine(); for (int i = 0; i < nvar; i++) { for (int j = 0; j < i; j++) { Console.WriteLine("CorrMtrxIn[" + i + "," + j + "] = " + CorrMtrxIn[i, j]); } } PrintMatrixFormat pmf = new PrintMatrixFormat(); pmf.NumberFormat = "0.000000000"; new PrintMatrix("Input Correlation Matrix: ").Print(pmf, CorrMtrxIn); Console.WriteLine("Correlation Matrices calculated from"); Console.WriteLine(" Gaussian Copula imprinted multivariate sequence:"); Console.WriteLine(); // Compute the Cholesky factorization of CorrMtrxIn Cholesky CholMtrx = new Cholesky(CorrMtrxIn); for (int kmax = 500; kmax < 1000000; kmax *= 10) { Console.WriteLine("# vectors in multivariate sequence: " + kmax); double[][] GCvart = new double[kmax][]; double[][] GCdevt = new double[kmax][]; for (int i2 = 0; i2 < kmax; i2++) { GCdevt[i2] = new double[nvar]; } Imsl.Stat.Random r = IMSLRandom(); for (int k = 0; k < kmax; k++) { GCvart[k] = r.NextGaussianCopula(CholMtrx); //probs for (int j = 0; j < nvar; j++) { /* * invert Gaussian Copula probabilities to deviates using * variable-specific inversions: j = 0: Chi Square; 1: F; * 2: Normal(0,1); will end up with deviate sequences ready * for mapping to canonical correlation matrix: */ if (j == 0) { //convert probs into ChiSquare(df=10) deviates: GCdevt[k][j] = InvCdf.Chi(GCvart[k][j], 10.0); } else if (j == 1) { //convert probs into F(dfn=15,dfd=10) deviates: GCdevt[k][j] = InvCdf.F(GCvart[k][j], 15.0, 10.0); } else { //convert probs into Normal(mean=0,variance=1) deviates: GCdevt[k][j] = InvCdf.Normal(GCvart[k][j]); } } } /* * extract Canonical Correlation matrix from arbitrarily distributed * deviate sequences GCdevt[k=0..kmax-1][j=0..nvar-1] which have been * imprinted with CorrMtrxIn[i=1..nvar][j=1..nvar] above: */ double[][] CorrMtrx = r.CanonicalCorrelation(GCdevt); double relerr; for (int i = 0; i < nvar; i++) { for (int j = 0; j < i; j++) { relerr = Math.Abs(1.0 - (CorrMtrx[i][j] / CorrMtrxIn[i, j])); Console.WriteLine("CorrMtrx[" + i + "][" + j + "] = " + CorrMtrx[i][j] + "; relerr = " + relerr); } } new PrintMatrix("Correlation Matrix: ").Print(pmf, CorrMtrx); } } }
Random Example 2: CorrMtrxIn[1,0] = -0.948683298050514 CorrMtrxIn[2,0] = 0.816496580927726 CorrMtrxIn[2,1] = -0.645497224367903 Input Correlation Matrix: 0 1 2 0 1.000000000 -0.948683298 0.816496581 1 -0.948683298 1.000000000 -0.645497224 2 0.816496581 -0.645497224 1.000000000 Correlation Matrices calculated from Gaussian Copula imprinted multivariate sequence: # vectors in multivariate sequence: 500 CorrMtrx[1][0] = -0.95029565568146; relerr = 0.00169957417218125 CorrMtrx[2][0] = 0.805260514673251; relerr = 0.0137613145197848 CorrMtrx[2][1] = -0.640202740166643; relerr = 0.00820217965529457 Correlation Matrix: 0 1 2 0 1.000000000 -0.950295656 0.805260515 1 -0.950295656 1.000000000 -0.640202740 2 0.805260515 -0.640202740 1.000000000 # vectors in multivariate sequence: 5000 CorrMtrx[1][0] = -0.948611734364931; relerr = 7.54347480657058E-05 CorrMtrx[2][0] = 0.815532446740145; relerr = 0.00118081840157325 CorrMtrx[2][1] = -0.646255361981671; relerr = 0.00117450174090306 Correlation Matrix: 0 1 2 0 1.000000000 -0.948611734 0.815532447 1 -0.948611734 1.000000000 -0.646255362 2 0.815532447 -0.646255362 1.000000000 # vectors in multivariate sequence: 50000 CorrMtrx[1][0] = -0.948314953115713; relerr = 0.00038826965285188 CorrMtrx[2][0] = 0.817827046711005; relerr = 0.00162948114463246 CorrMtrx[2][1] = -0.646669093465958; relerr = 0.00181545180028175 Correlation Matrix: 0 1 2 0 1.000000000 -0.948314953 0.817827047 1 -0.948314953 1.000000000 -0.646669093 2 0.817827047 -0.646669093 1.000000000 # vectors in multivariate sequence: 500000 CorrMtrx[1][0] = -0.94873295551488; relerr = 5.23435634081082E-05 CorrMtrx[2][0] = 0.81728795761655; relerr = 0.000969234540976194 CorrMtrx[2][1] = -0.646929884520875; relerr = 0.00221946756529401 Correlation Matrix: 0 1 2 0 1.000000000 -0.948732956 0.817287958 1 -0.948732956 1.000000000 -0.646929885 2 0.817287958 -0.646929885 1.000000000Link to C# source.