发帖
雷达卡
Here is the SAS pgm for the simulation data + estimation based on simulation data for a conditional logit model/a discreta choice model.
The estimation uses the procedures both SAS canned procedure(proc MDC) and me written codes with SAS proc nlmixed.
The results are much each other from both approaches. The results are concured with simulation parameters.
**************************************************************************************************************************;
/* cdf of extreme value type 1 is F(x)= exp(-exp(-x))
an inverse transform on the cdf yields
x = -log(-log(F(x)) */
%let size=500;
%let Eulerconstant=0.5772156649015328606;
data err_extreme;
seed=7813430;
do i=1 to &size;
e1=-log(-log(ranuni(seed)) ) - &Eulerconstant;
e2=-log(-log(ranuni(seed)) ) - &Eulerconstant;
e3=-log(-log(ranuni(seed)) ) - &Eulerconstant;
e4=-log(-log(ranuni(seed)) ) - &Eulerconstant;
output;
end;
keep e:;
run;
*verify the sample data;
proc corr data=err_extreme;
var e: ;
run;
data simudata;
set err_extreme;
***4 choices**;
seed= 6791534;
x1= rannor(seed); ***dependent on choice, cost/benefit of particular choice;
x2= rannor(seed); ***dependent on choice, cost/benefit of particular choice;
x3= rannor(seed); ***dependent on choice, cost/benefit of particular choice;
x4= rannor(seed); ***dependent on choice, cost/benefit of particular choice;
z1= ranuni(seed)-1.5; ***dependent on choice, cost/benefit of particular choice;
z2= ranuni(seed)-0.2; ***dependent on choice, cost/benefit of particular choice;
z3= ranuni(seed)-0.3; ***dependent on choice, cost/benefit of particular choice;
z4= ranuni(seed)-1.5; ***dependent on choice, cost/benefit of particular choice;
***choice related fixed cost/benefit only 3 can be identified;
*simulation parameters ;
int1=1;int2=2;int3=3;int4=0;
u1= int1 + 0.4 * x1 - 1*z1 + 1* e1;
u2= int2 + 0.4 * x2 - 1*z2 + 1* e2;
u3= int3 + 0.4 * x3 - 1*z3 + 1* e3;
u4= int4 + 0.4 * x4 - 1*z4 + 1* e4;
max=max(u1,u2,u3,u4);
if u1=max then do;
choice=1;
z=z1;
x=x1;
end;
else if u2=max then do;
choice=2;
z=z2;
x=x2;
end;
else if u3=max then do;
choice=3;
z=z3;
x=x3;
end;
else do;
choice=4;
z=z4;
x=x4;
end;
keep z: x: choice;
run;
proc print data=simudata(obs=10);
run;
proc freq data=simudata;
table choice;
run;
****shape the data for sas procedure;
data simudata2;
set simudata;
array xv(4) x1-x4;
array zv(4) z1-z4;
retain pid 0;
pid + 1;
do i = 1 to 4;
mode = i;
decision = ( choice = i );
dmy1=(mode=1);dmy2=(mode=2);dmy3=(mode=3);dmy4=(mode=4);
xx =xv(i);
zz =zv(i);
output;
end;
keep dmy1 dmy2 dmy3 xx zz decision pid mode;
run;
proc print data=simudata2(obs=10);
run;
***write your own ML estimation of discrete choice model***;
proc nlmixed data=simudata tech=DBLDOG ;
parms a1=0 a2=0 a3=0 b1=0.1 c1=-0.1 ;
a4=0.000;
v1= min(1e50, a1+ b1 * x1 + c1* z1) ;
v2= min(1e50, a2+ b1 * x2+ c1* z2) ;
v3= min(1e50, a3+ b1 * x3 + c1* z3) ;
v4= min(1e50, a4+ b1 * x4 + c1* z4) ;
pr1 =min(1e50, exp(v1));
pr2 =min(1e50, exp(v2));
pr3 =min(1e50, exp(v3));
pr4 =min(1e50, exp(v4));
tpr=pr1+pr2+pr3+pr4;
if choice=1 then p=pr1/tpr;
else if choice=2 then p=pr2/tpr;
else if choice=3 then p=pr3/tpr;
else if choice=4 then p=pr4/tpr;
p=max(1e-12,p);
log_like =log(p);
model choice ~ general(log_like);
run;
***SAS ***;
proc mdc data=simudata2( keep=dmy1 dmy2 dmy3 xx zz pid mode decision
rename=(dmy1=a1 dmy2=a2 dmy3=a3 xx=x zz = z)) ;
model decision = a1 a2 a3 x z / type=clogit choice=(mode 1 2 3 4)
covest=hess;
id pid;
run;
京公网安备 11010802022788号
