1-way ANOVA

proc import datafile="F:\publish\data analysis using R\data\plantgrowth.csv" out=plantgrowth dbms=csv replace;
run;

proc sql;
  select group, mean(weight), var(weight)
  from plantgrowth
  group by group;
quit;

proc glm data=plantgrowth;
  class group;
  model weight = group;
run; quit;

proc mixed data=plantgrowth;
  class group;
  model weight = group;
run;

proc mixed data=plantgrowth;
  class group;
  model weight = group / ddfm=Satterthwaite;
  repeated / group=group;
run;

paired t-test

proc import datafile="F:\publish\data analysis using R\data\CBT.csv" out=cbt dbms=csv replace;
run;

proc sql;
  select t(diff) as t, prt(diff) as pvalue
  from cbt
;quit;

proc sql;
  create table cbt_long as
  select monotonic() as id, 1 as period, prewt as weight from cbt union
  select monotonic() as id, 2 as period, postwt as weigt from cbt
;quit;

proc mixed data=cbt_long;
  class id period;
  model weight = period;
  random id;
run;

independent t-test

proc import datafile="F:\publish\data analysis using R\data\indept_ttest.csv" out=dental dbms=csv replace;
run; 

proc ttest data=dental;
  class treatment;
  var resp;
run;

proc mixed data=dental;
  class  treatment;
  model resp=treatment;
run;

proc mixed data=dental;
  class  treatment;
  model resp=treatment / ddfm=Satterthwaite;
  repeated / group=treatment;
run;

/* ML: df is not right */
proc nlmixed data=dental;
  parms mu=100, tau=10, s2e=2;
  if (treatment=1) then eta=mu; 
                   else eta=mu+tau; 
  model resp ~ normal(eta, s2e);
run;

Yij = mu + bi + error

proc import datafile="F:\publish\data analysis using R\data\rail.csv" out=rail dbms=csv replace;
run;

proc mixed data=rail;
  class rail;
  model travel= / s;
  random rail;
run;

proc mixed data=rail method=ml;
  class rail;
  model travel= / s;
  random rail;
run;

* GEE;
proc genmod data=rail;
  class rail;
  model travel = ;
  repeated subject = Rail / type=cs;
run;

proc nlmixed data=rail;
  parms mu=100, s2b=10, s2e=10;
  model travel ~ normal(mu+bi, s2e);
  random bi ~ normal(0, s2b) subject=rail;
run;

y = a + b*x + error

proc import datafile="F:\publish\data analysis using R\data\cars.csv" out=cars dbms=csv replace;
run;
 
proc mixed data=cars;
  model dist = speed / s;
run;

proc mixed data=cars method=ml;
  model dist = speed / s;
run;

proc genmod data=cars;
  model dist = speed;
run;

proc glimmix data=cars;
  model dist = speed / s;
run;

proc nlmixed data=cars;
  parms alpha=10 beta=10 s2e=100;
  model dist ~ normal(alpha+beta*speed, s2e);
run;

y = mu + e

data test;
  input y @@;
  datalines;
59.9 91.8 82.4 59.1 82.9 70.6 77.2 69.5 78.5 85.6
run;

proc sql;
  select mean(y) as mean,  css(y)/count(y) as MLE, css(y)/(count(y)-1) as REMLE
  from test
;quit;

proc mixed data=test;
  model y= / s;
run;

proc mixed data=test method=ml;
  model y= / s;
run;

proc genmod data=test;
  model y=;
run;

/* REML: default */
proc glimmix data=test method=rspl;
  model y= / s;
run;

/* ML */
proc glimmix data=test method=mspl;
  model y= / s;
run;

proc nlmixed data=test;
  parms mu=100 s2e=100;
  model y~normal(mu,s2e);
run;

Modeling McNemar's test

proc import datafile='F:\consulting\dich.csv' out=dich dbms=csv replace; run;

/* McNemar test */
proc freq data=dich;
by group;
table before*after / agree;
run;

proc sql;
create table dich2 as
select group, id, before="TRUE" as before, after="TRUE" as after
from dich
;quit;

/* Logistic Regression */
proc sql;
create table logist as
select group, after-before as diff
from dich2
;quit;

proc logistic data=logist (where=(diff ne 0));
by group;
model diff=;
run;

proc logistic data=logist (where=(diff ne 0));
class group;
model diff=group;
run;

/* Generalized Mixed Models */
proc sql;
create table dich3 as
select group, id, 0 as period, before as y from dich2 union
select group, id, 1 as period, after as y from dich2
;quit;

proc glimmix data=dich3;
by group;
class group id period;
model y = period / link=logit dist=binomial;
random id;
run;

proc glimmix data=dich3;
class group id period;
model y = period | group / link=logit dist=binomial;
random id;
run;

cross over design (binary output)

data crossover;
input group $ seq1 seq2 count;
seq='(' || put(seq1,1.) || ',' || put(seq2,1.) || ')';
diff=seq2-seq1;
trt=diff*(group='1(AB)') - diff*(group='2(BA)');
ig=-trt;
datalines;
1(AB) 0 0 5
1(AB) 0 1 14
1(AB) 1 0 3
1(AB) 1 1 6
2(BA) 0 0 10
2(BA) 0 1 7
2(BA) 1 0 18
2(BA) 1 1 10
;

proc freq;
weight count;
table group*seq1*seq2;
run;

/* Odds Ratio */
* carry-over effect;
proc freq data=crossover (where=(diff=0));
weight count;
table group*seq / chisq;
run;
proc freq data=crossover (where=(diff ne 0));
weight count;
table group*seq / chisq;
run;

proc logistic data=crossover (where=(diff ne 0));
weight count;
class group;
model diff = group;
run;

/* Prescott's test */
proc freq data=crossover;
weight count;
table group*diff / chisq;
run;

/*** McNemar's test ***/
* 1(AB);
proc freq data=crossover (where=(group="1(AB)"));
weight count;
table seq1*seq2 / agree;
run;
proc logistic data=crossover (where=(group="1(AB)" & diff ne 0));
weight count;
model diff =;
run;

* 2(BA);
proc freq data=crossover (where=(group="2(BA)"));
weight count;
table seq1*seq2 / agree;
run;
proc logistic data=crossover (where=(group="2(BA)" & diff ne 0));
weight count;
model diff =;
run;

* group differ ;
proc freq data=crossover (where=(diff ne 0));
weight count;
table ig*trt / agree;
run;