Hi everyone,
After reading this post: https://communities.sas.com/t5/SAS-Procedures/Proc-Multtest-on-non-parametric-data/td-p/116780 and the very helpful responses, I was wondering if anyone might be able to help me with a couple questions I had concerning the analysis of my own data.
Quick background: I am testing disease severity on detached potato leaves following 7 treatments (one of which is the ctrl). The disease severity was measured once (10 days following treatment) by estimating the infection percentage of each leaf (sampling unit). Leaves had been placed in boxes (6 leaves per box) so that environmental conditions remained very humid from the beginning. Boxes were essentially reps (trt), the experimental units. There were 3 boxes/reps per trt, 21 boxes total - subject to randomly allocated locations in a growth room. To collect enough data, the experiment was conducted 4 times, or over 4 periods. We are considering periods 1 thru 4 as blocks.
I understand that since I have one experimental factor + one blocking factor, I can use Friedman's test (Shah and Madden, 2004).
My question is:
Can I analyze this for multiple comparisons as SteveDenham described in the link above, such that:
My SAS code:
data CA9d10; /*proc mixed*/
input score trt rep leaf block;
cards;
(...)
/*Friedman test for a randomized complete block design with one experimental
factor + one blocking factor*/
proc freq data = CA9d10;
tables block*trt*score / cmh2 scores = rank noprint ;
title 'Friedman''s Chi-Square';
run ;
/* The test statistic and P-value are given in the */
/* "Row Mean Scores Differ" row of the output. */
/**************************************************************/
proc means median;
class trt;
var score;
run;
/*Non-parametric multiple comparison procedure*/
proc rank data=CA9d10 out=rankeddata;
var score;
ranks rank_score;
run;
/* Change to PROC MIXED to get LSMESTIMATEs, which can be adjusted for multiple comparisons, and to accommodate the comparisons which are not just one group to another. Multiple comparisons are the stepdown Bonferroni (Holm) adjustment */
proc mixed data=rankeddata ;
class trt;
model rank_score=trt;
lsmestimate trt 'trt 1 vs 2' 1 -1 0 0 0 0 0/adjust=bonferroni stepdown;
lsmestimate trt 'trt 1 vs 3' 1 0 -1 0 0 0 0/adjust=bonferroni stepdown;
lsmestimate trt 'trt 1 vs 4' 1 0 0 -1 0 0 0/adjust=bonferroni stepdown;
lsmestimate trt 'trt 1 vs 5' 1 0 0 0 -1 0 0/adjust=bonferroni stepdown;
lsmestimate trt 'trt 1 vs 6' 1 0 0 0 0 -1 0/adjust=bonferroni stepdown;
lsmestimate trt 'trt 1 vs 7' 1 0 0 0 0 0 -1/adjust=bonferroni stepdown;
lsmestimate trt 'trt 2 vs 3' 0 1 -1 0 0 0 0/adjust=bonferroni stepdown;
lsmestimate trt 'trt 2 vs 4' 0 1 0 -1 0 0 0/adjust=bonferroni stepdown;
lsmestimate trt 'trt 2 vs 5' 0 1 0 0 -1 0 0/adjust=bonferroni stepdown;
lsmestimate trt 'trt 2 vs 6' 0 1 0 0 0 -1 0/adjust=bonferroni stepdown;
lsmestimate trt 'trt 2 vs 7' 0 1 0 0 0 0 -1/adjust=bonferroni stepdown;
lsmestimate trt 'trt 3 vs 4' 0 0 1 -1 0 0 0/adjust=bonferroni stepdown;
lsmestimate trt 'trt 3 vs 5' 0 0 1 0 -1 0 0/adjust=bonferroni stepdown;
lsmestimate trt 'trt 3 vs 6' 0 0 1 0 0 -1 0/adjust=bonferroni stepdown;
lsmestimate trt 'trt 3 vs 7' 0 0 1 0 0 0 -1/adjust=bonferroni stepdown;
lsmestimate trt 'trt 4 vs 5' 0 0 0 1 -1 0 0/adjust=bonferroni stepdown;
lsmestimate trt 'trt 4 vs 6' 0 0 0 1 0 -1 0/adjust=bonferroni stepdown;
lsmestimate trt 'trt 4 vs 7' 0 0 0 1 0 0 -1/adjust=bonferroni stepdown;
lsmestimate trt 'trt 5 vs 6' 0 0 0 0 1 -1 0/adjust=bonferroni stepdown;
lsmestimate trt 'trt 5 vs 7' 0 0 0 0 1 0 -1/adjust=bonferroni stepdown;
lsmestimate trt 'trt 6 vs 7' 0 0 0 0 0 1 -1/adjust=bonferroni stepdown;
ods output lsmestimates=lsmestimates;
run;
/* Hochberg adjustment. See Example 61.5 Inputting Raw p-Values in the MULTTEST Procedure documentation. In particular, read the last paragraph in this example, describing some general assumptions regarding the methods of adjustment */
/* First, pare down the lsmestimates dataset to just the labels and the p values */
data lsmestimates_for_multtest;
set lsmestimates;
raw_p=probt;
keep label raw_p;
run;
/*Now running the adjustments. This gives adjustments using the stepdown Bonferroni (Holm, same as in PROC MIXED, hopefully), Hochberg and false discovery rate methods*/
proc multtest inpvalues=lsmestimates_for_multtest holm hoc fdr;
run;
Questions:
-I didn't include it in the SAS code yet, but shouldn't I be accounting for rep (trt) as a random effect?
-Is it permisible to test for potential significance of the interaction trt*block by including it in the model? I'm assuming if the interaction is significant I will therefore be forced to analyze the data block by block separately, correct? In this case kruskal-wallis could be used I guess...
-Which adjustment (stepdown bonferroni, holm, hochberg...) should I be using for p values? I'm kind of confused about this for non-parametric multiple comparisons...
I also found this suggestion online: http://www.lexjansen.com/sugi/sugi12/sd/001-12.pdf by David Ipe
which claims to prove mathematically that the EXACT Friedman test and an associated multiple comparison procedure can be reduced to (in SAS):
proc rank data = CA9d10 out = rdata;
by block;
var score;
ranks r;
run;
proc glm data = rdata;
class block trt;
model r = block trt;
lsmeans trt/pdiff;
run;
Questions:
-Again, can I add rep(trt) into the model statement and add a random statement:
proc glm data = rdata;
class block trt rep;
model r = block trt rep(trt) block*trt;
random rep(trt);
lsmeans trt/pdiff;
run;
-Should I be adding an adjustment here like in the first method? If so, which one?
I get similar results using both methods... do you see any differences in the two approaches and what would your suggestion for me be?
Thank you so much!!
-Nick
Click image to register for webinar
