/******************** RepeatedNoise.sas ***********************
* Participants listened to political discussions under different  *
* levels of background noise. "Discrimination score" is a measure *
* of how well they could tell what was being said.                *
******************************************************************/

title 'Noise data';

proc format;      value sexfmt    0 = 'Male'  1 = 'Female' ;

data loud;
     infile '/folders/myfolders/noise.data.txt';  /* Univariate data read */
     input ident  interest  sex  age  noise  time  discrim ;
     format sex sexfmt.;
     label interest = 'Interest in topic (politics)'
           time     = 'Order of presenting noise level';

proc freq;
     tables sex*age*noise / norow nocol nopercent;

proc glm;
     title2 'Classical mixed model repeated measures';
     class ident age sex noise;
     model discrim = age|sex|noise ident(age*sex);
     random ident(age*sex) / test; 
     means noise age;

/* Expected mean squares show that not all the hypotheses 
   are testable with the classical mixed model approach. 
   SAS is buying testability of the hypotheses by assuming 
   that you're only interested in an effect if all the 
   higher-order interactions involving that effect are 
   absent.
   
   Also, we can't use the covariate and it's not clear how 
   to do multiple comparisons.  */

proc mixed;
     title2 'Covariance structure (cs) with proc mixed';
     title3 'Replicate the classical tests'; 
     class age sex noise;
     model discrim = age|sex|noise;
     repeated / type=cs subject=ident;


   
proc mixed;
     title2 'Include covariate, multiple comparisons, contrasts';
     class age sex noise;
     model discrim = interest age|sex|noise / solution;
     /* solution shows the beta-hat values (interest) */
     repeated / type=cs subject=ident r;
     lsmeans age / pdiff tdiff adjust = bon;
     lsmeans noise / pdiff tdiff adjust = bon;
     contrast 'Noise Linear' noise  1 -2  1  0  0,
                             noise  0  1 -2  1  0,
                             noise  0  0  1 -2  1;
/* In the test for departure from linearity,
   mu1-mu2 = mu2-mu3 => 1 -2  1  0  0
   mu2-mu3 = mu3-mu4 => 0  1 -2  1  0
   mu3-mu4 = mu4-mu5 => 0  0  1 -2  1     */

/* Compare results with unknown covariance structure */

proc mixed cl;
     title2 'Unknown covariance structure';
     class age sex noise;
     model discrim = interest age|sex|noise;
     repeated / type=un subject=ident r;
     lsmeans age / pdiff tdiff adjust = bon;
     lsmeans noise / pdiff tdiff adjust = bon;
     contrast 'Noise Linear' noise  1 -2  1  0  0,
                             noise  0  1 -2  1  0,
                             noise  0  0  1 -2  1;


/* Test for difference between covariance structures
   using a full versus reduced model approach. 
   Difference between Likelihood Ratio chisquares is 
   chi-squared, with df = difference of degrees of freedom.
   H0 is compound symmetry, alternative is unstructured. 
   Any model for the covariance matrix is reduced compared to 
   the unknown structure. */

proc iml;
     title2 'Test for difference between covariance structures';
     title3 'Compound symmetry versus unknown';
     ChiSquared = 51.44-42.90; df = 14-1;
     pvalue = 1-probchi(ChiSquared,df);
     print ChiSquared df pvalue;