random normal generation with any covariance structure

是 否 +2 论坛币

k人 参与回答

经管之家送您一份

应届毕业生专属福利!

求职就业群

赵安豆老师微信：zhaoandou666

经管之家联合CDA

送您一个全额奖学金名额~ !

立即领取

感谢您参与论坛问题回答

经管之家送您两个论坛币！

+2 论坛币

Any covariance structure(symmetric positive definite COV)  can be decomposed into COV=L*L`.
L is lower triangle matrix. We can use Cholesky decomposition routine to solve for L.
We generate a group of standard normal random variable X. The X*L` will have COV as its covariance.

We can do the same job as vnormal routine in IML within a data step.

Here is an example.

COV=( 1  0.8 2,
           0.8 3  1,
           2   1  5)

The program results as


                Covariance Matrix, DF = 99999

                     x1                x2                x3

   x1       1.003938580       0.805219557       2.008134939
   x2       0.805219557       3.022280567       1.006021545
   x3       2.008134939       1.006021545       5.018633530



proc fcmp outlib=sasuser.funcs.temp;;
     subroutine Choleskydecom (mat[*,*],lmat[*,*] );
      outargs lmat;
      call chol(mat, lmat);
      call transpose(lmat, lmat);
      endsub;
run;
quit;

options cmplib=sasuser.funcs;

%let n=100000;
data random;
  array a(3,3) _temporary_   ( 1  0.8 2,
                               0.8 3  1,
                               2   1  5

                             );
  array b(3,3) _temporary_ ;
  call Choleskydecom (a, b);

  array x(3);
  array norm(3);
  seed=190;
  do j=1 to &n;
    do i=1 to dim(x);
   norm=rannor(seed);
x=.;
end;
do i=1 to dim(x);
do k=1 to i;
x=sum(x,b[k,i]*norm[k]);
end;

end;
output;
end;
run;

proc corr data=random cov ;
var x1 x2 x3;
run;

扫码加我 拉你入群

请注明：姓名-公司-职位

以便审核进群资格，未注明则拒绝

分享0 收藏1 回帖

关键词：Generation covariance Structure variance struct generation structure positive generate standard

norm and x should have index. I don't know why it become like that. It looks good in editing window.

norm-->norm[i]

x--->x[i]

norm[i]=rannor(seed);
x[i]=.;
end;
do i=1 to dim(x);
do k=1 to i;
x[i]=sum(x[i],b[k,i]*norm[k]);