data dinput; input diagnose treat r1 r2 r3 count @@; datalines; 0 0 1 1 1 16 0 0 1 1 0 13 0 0 1 0 1 9 0 0 1 0 0 3 0 0 0 1 1 14 0 0 0 1 0 4 0 0 0 0 1 15 0 0 0 0 0 6 0 1 1 1 1 31 0 1 1 1 0 0 0 1 1 0 1 6 0 1 1 0 0 0 0 1 0 1 1 22 0 1 0 1 0 2 0 1 0 0 1 9 0 1 0 0 0 0 1 0 1 1 1 2 1 0 1 1 0 2 1 0 1 0 1 8 1 0 1 0 0 9 1 0 0 1 1 9 1 0 0 1 0 15 1 0 0 0 1 27 1 0 0 0 0 28 1 1 1 1 1 7 1 1 1 1 0 2 1 1 1 0 1 5 1 1 1 0 0 2 1 1 0 1 1 31 1 1 0 1 0 5 1 1 0 0 1 32 1 1 0 0 0 6 run; data depress; set dinput; do i=1 to count; case+1; outcome=r1; time=1; q1=1; q2=0; output; outcome=r2; time=2; q1=0; q2=1; output; outcome=r3; time=4; q1=0; q2=0; output; end; run; *Marginal model with GEE; proc genmod data=depress descending; class case time; model outcome = diagnose treat time treat*time / dist=bin link=logit type3; repeated subject=case / type=exch corrw; run; *Random effects model; proc nlmixed data=depress qpoints=200; parms a=1 b1=-1 b2=2 b3=-1 b4=-0.5 b5=-2 b6=-1 sig=.1; eta = a+b1*diagnose+b2*treat+b3*q1+b4*q2+b5*q1*treat+b6*q2*treat+u; p = exp(eta)/(1+exp(eta)); model outcome ~ binary(p); random u ~ normal(0, sig*sig) subject=case; run; *Model with no random effects; proc nlmixed data=depress; parms a=1 b1=-1 b2=1 b3=-1.5 b4=-1 b5=-0.5 b6=-0.5; eta = a+b1*diag+b2*treat+b3*q1+b4*q2+b5*q1*diag+b6*q2*diag; p = exp(eta)/(1+exp(eta)); model outcome ~ binary(p); run;