(*:Mathematica:: V2.2 *) (*:Contex: "exfitz`" *) (*Title: binary cross-over data Kenward and Jones complete data *) (*References: model related to (but not identical with) Fitzmaurice, G.M. and Laird, N. M. 1993. A likelihood based method for analysing longitudinal binary responses, Biometrika, 1993, 80, 143-151. Zhao, L. P. and Prentice, R. P. 1990, Correlated binary regression using a quadratic exponential model, Biometrika, 77, 642-8 Zeger, S.L. and Liang,K. 1992. An overview of methods for the analysis of longitudinal data, Statistics in Medicine, 11, 1825-1839. Kenward M.G, and Jones, B. 1987. A log-linear model for binary cross-over data. Applied Statistics, 36 192-204. *) (*Keywords: *) (*Requirements: "clm`" *) (*History: 1992 Stuart G. Baker *) Clear[exfitz] (*the order is abnormal, normal*) (*the logit is Log[pr(abormal)/pr(normal)]*) (*the table in Zhao and Prentice has a 5 instead of a 4, in the abnormal/normal cell of group 2*) y=Transpose[{{6,0,6,22,9,4,2,18}}]; yadj=Transpose[{{6,.000001,6,22,9,4,2,18}}]; mod1="cts2"; mod2="ct2"; mod3="c2"; mod4="cts1"; mod5="ct1"; mod6="c1"; mod={mod1,mod2,mod3,mod4,mod5,mod6}; Clear[exfitz,createx1,createparametername] exfitz[y_,modelname_]:= Module[{q1,w1,g1,h1,z1,x1, func1,m1,components1,nv1,ncore1,arg1, nv,ncore,n,q,w,g,h,z,x, components,mstepfunc,msteparg, parametername,rationame, c,model,ratio,name}, (*model 1*) (*order of parameter*) (*p1,q1,pi1,p2,q2,pi2*) (*p1=1+;q1=++,pi=11*) q1=Log; w1=J[2,1] ~Dir~ {{0},{0},{0},{1}}; g1= Id[2] ~Dir~ {{0,0,1}, {1,0,-1}, {0,1,-1}, {-1,-1,1}}; h1= Exp[#2]/(1+Exp[#2])&; z1=J[6,1]; x1=createx1[modelname]; func1=NGLR; m1= Id[2] ~Dir~ {{1,1,0,0}, {0,0,1,1}, {1,0,1,0}, {0,1,0,1}, {1,0,0,0}, {0,1,1,1}}; nv1=J[8,1]; ncore1={{1}}; components1={{nv1,ncore1},{q1},{w1},{g1},{h1},{z1},{x1}}; arg1={m1,components1}; (*combination*) ncore=Dir[Id[2],J[1,4]].y; nv=ncore ~Dir~ J[4,1]; n={nv,ncore}; q={q1}; w={w1}; g={g1}; h={h1}; z={z1}; x={x1}; c=Id[8]; components={n,q,w,g,h,z,x}; mstepfunc={func1}; msteparg={arg1}; model={components,mstepfunc,msteparg}; ratio="None"; parametername=createparametername[modelname]; rationame=Automatic; name={parametername,rationame}; Return[{c,model,ratio,name}]]/;MemberQ[mod,modelname] (*starting values*) t1=Transpose[{{1.3,-.627,-.489,-.72,-1.24}}]; t={t1}; createx1[modelname_]:= Module[{constant,treatment,sequence,int2,int1 }, constant=J[6,1]; treatment={{0},{1},{0},{1},{0},{0}}; sequence={{0},{0},{0},{1},{1},{0}}; int2={{0,0},{0,0},{1,0},{0,0},{0,0},{0,1}}; int1={{0},{0},{1},{0},{0},{1}}; x1=Switch[modelname, "cts2",Hcat[constant,treatment,sequence,int2], "ct2", Hcat[constant,treatment,int2], "c2", Hcat[constant,int2], "cts1",Hcat[constant,treatment,sequence,int1], "ct1", Hcat[constant,treatment,int1], "c1", Hcat[constant,int1]]; Return[x1]] createparametername[modelname_]:= Switch[modelname, "cts2", {"constant","treatment","sequence", "theta1","theta2"}, "ct2", {"constant","treatment", "theta1","theta2"}, "c2", {"constant", "theta1","theta2"}, "cts1", {"constant","treatment","sequence", "theta"}, "ct1", {"constant","treatment","theta"}, "c1", {"constant","theta"}]