options formdlim = '-';

data asparagus;
  infile 'asparagus.txt' firstobs = 2;
  input trt year block yield;
  lyield =  log(yield);

/* repeated measures as a split plot, using proc mixed */
/*  outp = dsname creates a data set with residuals and pred values */

proc mixed  method = type3;
  class block year trt;
  model yield =  block trt year year*trt 
     /ddfm = kr outp = resids;
  random block*trt;
  estimate 'trt 1 - trt 2' trt 1 -1;
  estimate 'year 2 - year 3' year 1 -1;
    
  title 'MIXED rep. meas. as split plot';
run;
     
proc plot;
  plot resid*pred = year;
run;
       
/* another way to fit the split plot model, */
/*   by specifying the split plot correlation structure */
/* when you use a repeated statement, can not use type3 method */

/*  adding the rcorr and/or r options to the repeated statement */
/*    are often useful */
/*    RCORR gives you the within subject correl. matrix */
/*    R gives you the within subject covariance matrix */
/*  that are fit by the specified model */

/* a couple of small differences when you use a repeated statement */
/*   method = type3 does not work.  You need to use REML (the default) */
/*   and you really should use the kenward-rogers method (ddfm = kr) */
/*     to approximate the error d.f. */
/*   for a few correlation matrices (e.g. splitplot with balanced data), */
/*    kr is same as satterthwaite.  */
/*   when you use AR(1) or other structures, KR works much better */

proc mixed;
  class block year trt;
  model yield =  block trt year year*trt /ddfm = kr;

  repeated year /type = cs subject = block*trt r rcorr;
  
  estimate 'trt 1 - trt 2' trt 1 -1;
  estimate 'year 2 - year 3' year 1 -1;
    
  title 'MIXED rep. meas. as split plot using repeated';
run;


/* repeated measures with ar(1) correlated errors */

proc mixed;
  class block year trt;
  model yield =  block trt year year*trt /ddfm = kr;
  repeated year /type = ar(1) subject = block*trt r rcorr;
  estimate 'trt 1 - trt 2' trt 1 -1;
  estimate 'year 2 - year 3' year 1 -1;
    
  title 'MIXED rep. meas. with AR(1) temporal correlation';
run;


/* repeated measures with a  AR(1)+RE model */
/*  needs both a repeated and a random statement */
/*  the repeated specifies the AR(1) part, */
/*  the random adds the additional random effect */
/*  sometimes these models are hard to fit */
/*   if so, probably don't have enough data */
/*   or the model isn't appropriate */
/*  check the log window.  if SAS complains about: */
/*    not converging or */
/*    matrix non-positive definite */
/*  model is probably not appropriate */
/*  if the biology says it is, see me for help */

proc mixed;
  class block year trt;
  model yield =  block trt year year*trt /ddfm = kr;

  random block*trt;  
  repeated year /type = ar(1) subject = block*trt r rcorr;

  estimate 'trt 1 - trt 2' trt 1 -1;
  estimate 'year 2 - year 3' year 1 -1;
    
  title 'MIXED rep. meas. with AR(1)+RE temporal correlation';
run;

/* other correlation structures fit with different type = */
/* repeated measures with EXP(year) model */
/*   this is like AR(1) but for unequally spaced obs */
/*   correl =  exp(-(separation in time ) / param ) */

/*  repeated year /type = sp(exp)(year) subject = block*trt r rcorr;

/* repeated measures with a AR(1) and heterogenous variances */
/*   repeated year /type = arh(1) subject = block*trt r rcorr;

/* repeated measures with a  ARH(1)+RE  */
/*  random block*trt;
/*  repeated year /type = arh(1) subject = block*trt r rcorr;

/* repeated measures with a  ANTE(1) model */
/*  repeated year /type = ante(1) subject = block*trt r rcorr;

/* repeated measures with an unstructured correlation matrix */
/*  R and RCORR are especially useful here because */
/*  the UN model is so general, the matrices are */
/*  the correlations (or covariances) in the data */
/*  repeated year /type = un subject = block*trt rcorr r;

/* these are the common (CS, AR(1), UN) and some */
/*  useful not-so-common correlation models */
/* there are many more.  The SAS help files list them all */