请教关于SAS程序调试

在看一篇关于缺失值处理的文献，文献末给出proc iml程序，想改来试一试，发现跑不出来，我自己也看过一点proc iml过程的书，但无法解决这个问题，现po出其中一个模拟过程程序，希望有缘的高手和我一起尝试，看看如何运行成功。
以下程序目的：生成数据集（附件）
proc iml;
%let studypgmcode=C:\Users\SYT\Desktop\codes;
%let studysimdata=C:\Users\SYT\Desktop\dataset;
libname simdata "&studysimdata";
%include "&studypgmcode\ordinal.sas";

%let nsim=1000;
n=120;
n1=n/2;
n0=n-n1;

marginal={
0.25 0.25 0.25,
0.25 0.25 0.25,
0.25 0.25 0.25,
0.25 0.25 0.25
};
corr={
1.0 0.4 0.4,
0.4 1.0 0.68,
0.4 0.68 1.0
};

pgender=0.5;
pdaicomp=0.15;
psite={0.2 0.5 0.3};
meanp2h=16.4;
stdp2h=4.1;
ptrt=0.5;
dur_dml=0;
dur_dmu=20;
bmil=19;
bmiu=38;


do nsima=1 to ≁
kk= j(n1,1,nsima);  /*设置随机数id，共n1条观测*/
call ordinal(marginal,corr,n1,result,rc);
x7=result[,1];
x8=result[,2];
v1=result[,3];
call randseed(123);
x1=j(n1,1);call randgen(x1,"Bern",pgender);
call randseed(1234);
x4=j(n1,1);call randgen(x4,"Bern",pdaicomp);
call randseed(12345);
x5=j(n1,1);call randgen(x5,"Table",psite);
call randseed(456);
x6=j(n1,1);call randgen(x6,"normal",meanp2h,stdp2h);
call randseed(1235);
trt=j(n1,1);call randgen(trt,"Bern",0.5);
call randseed(456);
xx2=j(n1,1);call randgen(xx2,"Uniform");
call randseed(456);
xx3=j(n1,1);call randgen(xx3,"Uniform");

x2=dur_dmu*(0.5-xx2)+0.5*(dur_dml+dur_dmu);
x3=bmil*(0.5-xx3)+0.5*(bmil+bmiu);

xx=kk||x7||x8||v1||x1||x4||x5||x6||trt||x2||x3;
xxx=xx//xxx;

end;

varnames={"kk" "x7" "x8" "v1" "x1" "x4" "x5" "x6" "trt" "x2" "x3"};
create temp1 from xxx[colname=varnames];/*a matrix that contains names*/
append from xxx;
free varnames;
/*对照组*/
do nsima=1 to ≁
kk=j(n0,1,nsima);   /*设置随机数id，共n0条观测*/

call ordinal(marginal,corr,n0,result,rc);
x7=result[,1];
x8=result[,2];
v1=result[,3];
call randseed(123);
x1=j(n0,1);call randgen(x1,"Bern",pgender);
call randseed(1234);
x4=j(n0,1);call randgen(x4,"Bern",pdaicomp);
call randseed(12345);
x5=j(n0,1);call randgen(x5,"Table",psite);
call randseed(456);
x6=j(n0,1);call randgen(x6,"normal",meanp2h,stdp2h);
call randseed(1235);
trt=j(n0,1);call randgen(trt,"Bern",0.5);
call randseed(456);
xx2=j(n0,1);call randgen(xx2,"Uniform");
call randseed(456);
xx3=j(n0,1);call randgen(xx3,"Uniform");

x2=dur_dmu*(0.5-xx2)+0.5*(dur_dml+dur_dmu);
x3=bmil*(0.5-xx3)+0.5*(bmil+bmiu);

xx0=kk||x7||x8||v1||x1||x4||x5||x6||trt||x2||x3;
xxx0=xx0//xxx0;

end;

varnames={"kk" "x7" "x8" "v1" "x1" "x4" "x5" "x6" "trt" "x2" "x3"};
create temp0 from xxx0[colname=varnames];/*a matrix that contains names*/
append from xxx0;
free varnames;
quit;

/*不同缺失比例模拟10% 20% 30%*/
%macro simdata_anaset(totalmiss=10, n=120, beta=-0.05);
data temp11;
set temp1;/*处理组数据集*/
trt0=1;
run;

data temp10;
set temp0;/*对照组数据集*/
trt0=0;
run;

data allrandgen;
call streaminit(123);
set temp11 temp10;

beta0=4.85;
beta1=β/*TRT*/
beta2=0.31;/*Beseline HbA1c*/
sigma0=1;

/*Missing=10%*/
if &totalmiss=10 then do;
x7a=-1.9735;
x8a=-2.2546;
x78a=0.5603;

if trt=1 then v11a=-0.0439;
else if trt=0 then v11a=-0.1668;

ps=1/(1+1/exp(7.6175+0.4112*x1+0.0331*x2+0.0678*x3-0.1110*x4-0.1178*x5
+0.0526*x6+x7a*x7+x8a*x8+x78a*x7*x8+0.1227*10*trt+v11a*vl*v1));
obser=rand('BERNOULLI',ps);
simoutcpl=beta0+beta1*trt+beta2*x8+sigma0*rannor(123);

if obser=1 then do;
simout=simoutcpl;
end;
else if obser=0 then do;
simout=.;
end;
output;
end;

/*Missing=20%*/
else if &totalmiss=20 then do;
x7a=-2.0250;
x8a=-2.3295;
x78a=0.4631;

if trt=1 then v11a==-0.0439;
else if trt=0 then v11a=-0.1668;

ps=1/(1+1/exp(7.6175+0.4112*x1+0.0331*x2+0.0678*x3-0.1110*x4-0.1178*x5
+0.0526*x6+x7a*x7+x8a*x8+x78a*x7*x8+0.1227*10*trt+v11a*vl*v1));
obser=rand('BERNOULLI',ps);
simoutcpl=beta0+beta1*trt+beta2*x8+sigma0*rannor(123);

if obser=1 then do;
simout=simoutcpl;
end;
else if obser=0 then do;
simout=.;
end;
output;
end;

/*Missing=30%*/
else if &totalmiss=30 then do;
x7a=-2.1013;
x8a=-2.4361;
x78a=0.4139;

if trt=1 then v11a==-0.0439;
else if trt=0 then v11a=-0.1668;

ps=1/(1+1/exp(7.6175+0.4112*x1+0.0331*x2+0.0678*x3-0.1110*x4-0.1178*x5
+0.0526*x6+x7a*x7+x8a*x8+x78a*x7*x8+0.1227*10*trt+v11a*vl*v1));
obser=rand('BERNOULLI',ps);
simoutcpl=beta0+beta1*trt+beta2*x8+sigma0*rannor(123);

if obser=1 then do;
simout=simoutcpl;
end;
else if obser=0 then do;
simout=.;
end;
output;
end;

drop x7a x8a x78a beta0 beta1 beta2 sigma0;
run;

proc sort data=allrandgen;
by kk trt;
run;

data simdata.simd&n.&totalmiss;
set allrandgen;
run;

%mend simdata_anaset;
%simdata_anaset(10, 120, -0.05);

或者有其它方法实现模拟过程，也可以介绍一下，不过我主要用的软件是SAS