proc glm data=ma2612.asphalt plots (only)=diagnostics; class site type; model wear=type site; means type/cldiff lsd bon tukey; output out=asphalt_out p=predicted; run; quit; * The next bit of code constructs the Tukey 1 df for nonadditivity test * *Get mean of wear from output: it equals 53.0 *Rescale the predicted values; title 'Tukey Test for Interaction'; data asphalt_out1; set asphalt_out; rescale=predicted**2/(2*53); run; proc glm data=asphalt_out1; class site type; model wear=site type rescale; contrast 'Tukey eta' rescale 1; ods select contrasts; run; quit; /* The interaction term is not significant, so we can analyze the fit of the additive model.*/ title 'Interaction Plot 1'; proc glm data=ma2612.asphalt; class site type; model wear=site type site*type; ods select intplot; run; quit; title 'Interaction Plot 2'; proc glm data=ma2612.asphalt; class type site; model wear=type site type*site; ods select intplot; run; quit; title ' '; /* There are a couple of ways to obtain the RCBD model parameter estimates. Perhaps the easiest is to use the following macro. This macro will get the treatment and block effects for all but the last level of each. Since the effects add to 0, the last can be obtained by subtraction. */ %macro rcbd_est(data,response,treatment,block); proc transreg data=&data design; model class(&treatment &block/EFFECT); output out=xdesign; run; data xdesign; merge &data xdesign; run; proc reg data=xdesign; model &response= &_TrgInd; ods select ParameterEstimates; run; quit; %mend rcbd_est; title 'Parameter Estimates Using proc transreg'; * This call to the macro will get the results for the asphalt example; %rcbd_est(ma2612.asphalt,wear,type,site); /* The code below for the asphalt data shows how to estimate RCBD model parameters using proc glm. The model is y_ij=mu+tau_i+beta_j+epsilon_ij, where i=1,2,3,4 and j=1,...,9. If the factor has k levels, the estimate command has form estimate 'identifier' factor_name sequence_of_k_integers/divisor=k; So to estimate type2, o the identifier I used was type 2 (you can use what you like) o the factor_name is type o the sequence of 4 integers has the number of levels minus 1 (4-1=3) in position 2 (since we're estimating level 2 of type) and -1 in the other 3 places o the divisor equals 4. */ title 'Parameter Estimates Using proc glm'; proc glm data=ma2612.asphalt plots (only)=none; class site type; model wear=type site; estimate 'mu' intercept 1; estimate 'type 1' type 3 -1 -1 -1/divisor=4; estimate 'type 2' type -1 3 -1 -1/divisor=4; estimate 'type 3' type -1 -1 3 -1/divisor=4; estimate 'type 4' type -1 -1 -1 3/divisor=4; estimate 'site 1' site 8 -1 -1 -1 -1 -1 -1 -1 -1/divisor=9; estimate 'site 2' site -1 8 -1 -1 -1 -1 -1 -1 -1/divisor=9; estimate 'site 3' site -1 -1 8 -1 -1 -1 -1 -1 -1/divisor=9; estimate 'site 4' site -1 -1 -1 8 -1 -1 -1 -1 -1/divisor=9; estimate 'site 5' site -1 -1 -1 -1 8 -1 -1 -1 -1/divisor=9; estimate 'site 6' site -1 -1 -1 -1 -1 8 -1 -1 -1/divisor=9; estimate 'site 7' site -1 -1 -1 -1 -1 -1 8 -1 -1/divisor=9; estimate 'site 8' site -1 -1 -1 -1 -1 -1 -1 8 -1/divisor=9; estimate 'site 9' site -1 -1 -1 -1 -1 -1 -1 -1 8/divisor=9; ods select estimates; run; quit;