options ls=80 FORMDLIM="#"; title 'FRUITFLY DATA ANALYSIS -- with multiple comparisons'; /* Read in data from webpage (url) */ /* data fruitfly; infile "http://tofu.byu.edu/stat230/fruitfly3.txt" url; input fly companions type lifespan thorax sleeppercent group $; run; */ /* OR */ /* Read in data by pasting data file into SAS code using datalines */ data fruitfly; input fly companions type lifespan thorax sleeppercent group $; datalines; 1 8 0 35 0.64 22 preg8 2 8 0 37 0.68 9 preg8 3 8 0 49 0.68 49 preg8 4 8 0 46 0.72 1 preg8 5 8 0 63 0.72 23 preg8 6 8 0 39 0.76 83 preg8 7 8 0 46 0.76 23 preg8 8 8 0 56 0.76 15 preg8 9 8 0 63 0.76 9 preg8 10 8 0 65 0.76 81 preg8 11 8 0 56 0.80 12 preg8 12 8 0 65 0.80 15 preg8 13 8 0 70 0.80 37 preg8 14 8 0 63 0.84 24 preg8 15 8 0 65 0.84 26 preg8 16 8 0 70 0.84 17 preg8 17 8 0 77 0.84 14 preg8 18 8 0 81 0.84 14 preg8 19 8 0 86 0.84 6 preg8 20 8 0 70 0.88 25 preg8 21 8 0 70 0.88 18 preg8 22 8 0 77 0.92 26 preg8 23 8 0 77 0.92 24 preg8 24 8 0 81 0.92 29 preg8 25 8 0 77 0.94 27 preg8 1 0 9 40 0.64 18 nocomp 2 0 9 37 0.70 6 nocomp 3 0 9 44 0.72 19 nocomp 4 0 9 47 0.72 7 nocomp 5 0 9 47 0.72 16 nocomp 6 0 9 47 0.76 13 nocomp 7 0 9 68 0.78 35 nocomp 8 0 9 47 0.80 2 nocomp 9 0 9 54 0.84 35 nocomp 10 0 9 61 0.84 6 nocomp 11 0 9 71 0.84 15 nocomp 12 0 9 75 0.84 14 nocomp 13 0 9 89 0.84 18 nocomp 14 0 9 58 0.88 50 nocomp 15 0 9 59 0.88 25 nocomp 16 0 9 62 0.88 10 nocomp 17 0 9 79 0.88 33 nocomp 18 0 9 96 0.88 43 nocomp 19 0 9 58 0.92 35 nocomp 20 0 9 62 0.92 17 nocomp 21 0 9 70 0.92 27 nocomp 22 0 9 72 0.92 22 nocomp 23 0 9 75 0.92 16 nocomp 24 0 9 96 0.92 20 nocomp 25 0 9 75 0.94 37 nocomp 1 1 0 46 0.64 23 preg1 2 1 0 42 0.68 4 preg1 3 1 0 65 0.72 20 preg1 4 1 0 46 0.76 42 preg1 5 1 0 58 0.76 9 preg1 6 1 0 42 0.80 32 preg1 7 1 0 48 0.80 66 preg1 8 1 0 58 0.80 28 preg1 9 1 0 50 0.82 10 preg1 10 1 0 80 0.82 4 preg1 11 1 0 63 0.84 12 preg1 12 1 0 65 0.84 17 preg1 13 1 0 70 0.84 12 preg1 14 1 0 70 0.84 23 preg1 15 1 0 72 0.84 40 preg1 16 1 0 97 0.84 18 preg1 17 1 0 46 0.88 10 preg1 18 1 0 56 0.88 38 preg1 19 1 0 70 0.88 7 preg1 20 1 0 70 0.88 23 preg1 21 1 0 72 0.88 36 preg1 22 1 0 76 0.88 9 preg1 23 1 0 90 0.88 21 preg1 24 1 0 76 0.92 62 preg1 25 1 0 92 0.92 36 preg1 1 1 1 21 0.68 23 virgin1 2 1 1 40 0.68 62 virgin1 3 1 1 44 0.72 28 virgin1 4 1 1 54 0.76 18 virgin1 5 1 1 36 0.78 10 virgin1 6 1 1 40 0.80 28 virgin1 7 1 1 56 0.80 22 virgin1 8 1 1 60 0.80 29 virgin1 9 1 1 48 0.84 15 virgin1 10 1 1 53 0.84 73 virgin1 11 1 1 60 0.84 10 virgin1 12 1 1 60 0.84 5 virgin1 13 1 1 65 0.84 13 virgin1 14 1 1 68 0.84 27 virgin1 15 1 1 60 0.88 20 virgin1 16 1 1 81 0.88 21 virgin1 17 1 1 81 0.88 12 virgin1 18 1 1 48 0.90 49 virgin1 19 1 1 48 0.90 17 virgin1 20 1 1 56 0.90 22 virgin1 21 1 1 68 0.90 71 virgin1 22 1 1 75 0.90 17 virgin1 23 1 1 81 0.90 10 virgin1 24 1 1 48 0.92 24 virgin1 25 1 1 68 0.92 18 virgin1 1 8 1 16 0.64 34 virgin8 2 8 1 19 0.64 6 virgin8 3 8 1 19 0.68 4 virgin8 4 8 1 32 0.72 22 virgin8 5 8 1 33 0.72 28 virgin8 6 8 1 33 0.74 31 virgin8 7 8 1 30 0.76 16 virgin8 8 8 1 42 0.76 27 virgin8 9 8 1 42 0.76 8 virgin8 10 8 1 33 0.78 32 virgin8 11 8 1 26 0.80 20 virgin8 12 8 1 30 0.80 35 virgin8 13 8 1 40 0.82 12 virgin8 14 8 1 54 0.82 14 virgin8 15 8 1 34 0.84 17 virgin8 16 8 1 34 0.84 29 virgin8 17 8 1 47 0.84 31 virgin8 18 8 1 47 0.84 6 virgin8 19 8 1 42 0.88 30 virgin8 20 8 1 47 0.88 27 virgin8 21 8 1 54 0.88 40 virgin8 22 8 1 54 0.88 19 virgin8 23 8 1 56 0.88 8 virgin8 24 8 1 60 0.88 8 virgin8 25 8 1 44 0.92 15 virgin8 ; run; proc print data=fruitfly; run; proc glm data=fruitfly; class group; model lifespan = group; /* The next 4 lines show which groups are sig diff, but do not give C.I.'s */ means group / t ; means group / bon ; means group / tukey ; means group / scheffe ; /* The next 4 lines show all pairwise C.I.'s */ *means group / t cldiff; *means group / bon cldiff; *means group / tukey cldiff; *means group / scheffe cldiff; run; title 'Fruitfly analysis -- with CONTRASTS'; proc glm data=fruitfly; class group; model lifespan = group / clparm; /* Contrast statements give F tests */ contrast "nocomp - comps" group 1 -.25 -.25 -.25 -.25 / e; /* The '/ e' option tells SAS to print the complete contrast being used for the comparison */ contrast "pregnant - virgin" group 0 .5 .5 -.5 -.5; contrast "1comp - 8comp" group 0 .5 -.5 .5 -.5; contrast "interaction" group 0 .5 -.5 -.5 .5; /* Estimate statements give t tests (and potentially C.I.s, if you use clparm option in model statement) */ *estimate "nocomp - comps" group 1 -.25 -.25 -.25 -.25 / e; *estimate "pregnant - virgin" group 0 .5 .5 -.5 -.5; *estimate "1comp - 8comp" group 0 .5 -.5 .5 -.5; *estimate "interaction" group 0 .5 -.5 -.5 .5; run; /* How SAS handles the fruitfly data as an unbalanced BF[2] -- Type III analysis below agrees with analysis above*/ proc glm data=fruitfly; class companions type; model lifespan = companions|type; run; title 'FABRIC WEAR DATA ANALYSIS'; /* You can read in the file as below from your local drive... */ /* data wear; infile 'My Documents\stat230\wear.csv' delimiter="," firstobs=2; INPUT prop filler wear1 wear2 wear3; run; */ /* ... OR You can read in the file from the web (url) ... */ /* filename wearfile url "http://tofu.byu.edu/stat230/wear.csv"; *first create a nickname 'wearfile' for the data file; data wear; infile wearfile delimiter="," firstobs=2; * now I read in wearfile with some delimiter and firstobs options; INPUT prop filler wear1 wear2 wear3; run; */ /* OR, you can read in data by pasting the data file into SAS code using datalines */ data wear; INPUT prop filler wear1 wear2 wear3; datalines; 1 1 155 169 151 1 1 173 152 141 1 2 137 82 77 1 2 160 82 83 2 1 198 187 176 2 1 177 196 167 2 2 129 94 78 2 2 98 89 48 3 1 235 225 166 3 1 229 270 183 3 2 155 76 92 3 2 132 105 67 ; run; proc print; run; proc glm data=wear; class filler prop; /* MAKE SURE to enter your class variables in the order of your interaction variable 'filler*prop' */ model wear1 = filler prop filler*prop / clparm; /* For this example, I do not adjust the alpha for multiple comparisons because I am just illustrating different contrasts. In an actual analysis, I would adjust alpha if I were evaluating the impact of a group of non-orthogonal contrasts. For k non-orthogonal contrasts, this is generally done using Bonferroni, where we can simply say that a p-value is significant if the p-value is less than alpha/k. */ *means prop / tukey; *means filler / tukey; *means filler*prop ; estimate 'prop75-25' prop -1 0 1 /e; estimate 'prop50-25' prop -1 1 0 /e; estimate 'prop75-50' prop 0 -1 1 /e; estimate 'filler1-filler2' filler 1 -1 / e; estimate 'filler1-filler2 ALT' filler 1 -1 filler*prop .33333333 .33333333 .33333333 -.33333333 -.33333333 -.33333333 / e; /* Gives same result as previous contrast */ estimate 'f1prop75 - f2prop75' filler 1 -1 filler*prop 0 0 1 0 0 -1 / e; estimate 'f1p2 - f2' filler 1 -1 prop -.333333333 .666666666 -.3333333333 filler*prop 0 1 0 -.333333333 -.333333333 -.333333333 / e; estimate 'f1p2 - f2 : same thing' filler 3 -3 prop -1 2 -1 filler*prop 0 3 0 -1 -1 -1 / e divisor=3; run;