求教蒙特卡洛随机数据生成

感谢关注！

我有2个时间序列，Y(t)和X(t)，通过VAR(1)检验如下关系。
Y(t)=a+b*X(t-1)+u;
X(t)=c+d*X(t-1)+v;

回归结果为a=0.02; b=0.288；c=0.04; d=0.562;
冲击u的方差为0.007；冲击v的方差为0.0014；冲击u和冲击v的协方差为-0.001；冲击u和v的相关系数-0.319。

先在希望根据冲击u和冲击v的相关系数为限制条件：冲击u和冲击v符合正态分布
1）冲击v的方差为0.0014；

2）冲击u的方差为0.007；

3）冲击u和v的相关系数-0.319；

4）b=0.288和d=0.562;

5)  X(0)初始值为0.038；


6)  Y(0)初始值为0.012；

生成随机数据Y(t)和X(t).


随机生成序列长度为100+T （T=216），重复50000次，
然后剔除前100个观察值，用T=216进行线性回归，Y(t)=a+b*X(t-1)+u; 求出随机生成序列的b的平均值。


请高手指教，如何code.

u,v服从正态分布？
如果是的话你的问题就是生成一组二维的正态分布随机数
版上过去有人问过这个问题，你搜搜看先

proc iml;
a=0.02; b=0.288;c=0.04; d=0.562;
T=216;
N=100+T;
loop=5000;       /*循环5000次*/
temp_b=j(loop,1,0);
do i=1 to loop;
cov={0.007 -0.001,
     -0.001 0.0014};
mu={0 0};
uv=randnormal(N,mu,cov);
y=j(N,1,0);
x=y;
y[1]=0.012;
x[1]=0.038;
do j=1 to N-1;
y[j+1]=a+b*x[j]+uv[j,1];
x[j+1]=c+d*x[j]+uv[j,2];
end;
x=x[100:N-1];
y=y[101:N];
temp_b=t(x-x[:])*(y-y[:])/(t(x-x[:])*(x-x[:]));
end;
para_b=temp_b[:];
print para_b;
quit;