/* Lynch1.sas */
title 'Hovland and Sears Lynching data';

data lynch;
     infile '/folders/myfolders/Lynch.data.txt' firstobs=7; 
                                  /* Skipping the header */
     input Year Ayres Cotton1 Cotton2 Blynch Totlynch;
     label ayres = 'Ayres'' composite economic index'
           cotton1 = 'Per-acre value of cotton'
           cotton2 = 'Farm value of cotton'
           blynch  = 'Negro lynchings'
           totlynch = 'Total lynchings';
     olynch = totlynch-blynch;
     label olynch = 'Non-Black lynchings';
     one=1; /* Used below -- one "subject." */

proc sgplot;
     title2 'Lynchings and Cotton Prices over Time';
     series x=year y=cotton2;
     series x=year y=totlynch / y2axis;

proc means;
     var ayres -- olynch;
proc corr;
     var year -- olynch;


proc reg plots=none;
     title2 'Naive regression, but request Durbin-Watson statistic';
     title3 'Lower Durbin-Watson critical value = 1.5';
     model totlynch = year cotton2 / dw;
     output out = Dixie  residual = epsilonhat; /* Save residuals */

proc arima;
     title2 'Autocorrelations and partial autocorrelations';
     identify var = epsilonhat;

/* That looks like an AR(1) */

/* Apparently SAS university Edition does not include proc autoreg, 
   so this does not work. 
   
*********************************************************
proc autoreg; 
     title2 'Fit AR(1) with proc autoreg';
     model totlynch = year cotton2 /  nlag=1 method=ml;
*********************************************************/

proc mixed;
     title2 'Fit AR(1) with proc mixed';
     model  totlynch = year cotton2 / solution;
     repeated / type=ar(1) subject=one;

/* Test for difference between AR(1) and ARMA(1,1), using
   a likelihood ratio chi-squared test. AR(1) is the
   reduced model. Fit the models with maximum likelihood. */

proc mixed method=ml;
     title2 'Reduced model is AR(1)';
     model  totlynch = year cotton2 / solution;
     repeated / type=ar(1) subject=one;

proc mixed method=ml;
     title2 'Full model is Autoregressive Moving Average: ARMA(1,1)';
     model  totlynch = year cotton2 / solution;
     repeated / type=arma(1,1) subject=one;

proc iml;
     title2 'Difference between -2 Log Likelihood values is chi-squared';
     Chisquare = 443.7-441.5; df=1; pvalue = 1-probchi(Chisquare,df);
     print Chisquare df pvalue;