/*First order factor model*/
proc nlmixed data=kaz;
Parms b1=0 b2=0 b3=0 b4=0 b5=1 b6=1 b7=1 b8=1;
Eta = b1*d1 + b2*d2 + b3*d3 + b4*d4 + (b5*d1 + b6*d2 + b7*d3 + b8*d4)*fwm;
expecta = exp(eta);
p = expecta/(1 + expecta);
model response ~ binary(p);
random fwm ~ normal(0,1) subject=ID;
run;
PROC PRINTTO PRINT=PRINT;
run;
quit;

/*Second order factor model*/
proc nlmixed data=final;
Parms b1=0 b2=0 b3=0 b4=0 b5=1 b6=1 b7=1 b8=1 cov=0.5;
Eta = b1*d1 + b2*d2 + b3*d3 + b4*d4 + (b5*d1 + b6*d2)*fw + (b7*d3 + b8*d4)*fm;
expecta = exp(eta);
p = expecta/(1 + expecta);
model response ~ binary(p);
random fw  fm ~ normal([0,0],[1,cov,1]) subject=teacher_ID;
run;