/* Gaussian Copula */
new; cls;
rho=0.9;
n=1000;
U=Cgaus(rho,n);

library pgraph;
graphset;
_plctrl=-1;
xy(U[.,1],U[.,2]);




proc Cgaus(rho,n);
   local Sig,A,dim,Z,X,U;
   Sig=1~rho|
       rho~1;
   A=chol(Sig);
   dim=cols(Sig);
   Z=rndn(n,dim);
   Z=(Z-meanc(Z)')*inv(chol(vcx(Z))); /* This moment match is optional. You could erase. */
   X=Z*A;
   U=cdfn(X);
   retp(U);
endp;