0

I am having trouble finding a post hoc test to decipher at what "Session" or time I have a treatment within session affect.

This is my data:

TR  SESSION MAC FISHD      ID
1   1      1    3.285714286 1
2   1      2    0.571428571 2
2   1      3    3.571428571 3
1   1      4    4           4
1   2      1    4           5
2   2      2    6.571428571 6
2   2      3    3.142857143 7
1   2      4    8.857142857 8
1   3      1    0.714285714 9
2   3      2    1.714285714 10
2   3      3    4.428571429 11
1   3      4    0.714285714 12

This is how I got repeated measures:

model.b = lme(FISHD ~ TR + SESSION + TR*SESSION, 
            random = ~1|MAC,
            data=TTDall2)

> ACF(model.b)
  lag        ACF
1   0  1.0000000
2   1 -0.7547232
3   2  0.4852727

> model2 = lme(FISHD ~ TR + SESSION + TR*SESSION, 
+ random = ~1|MAC,
+ correlation = corAR1(form = ~ SESSION | MAC,
+                                                                                                                                               value = -0.7547232),
+   data=TTDall2,
+  method="REML")
> Anova(model2)
Analysis of Deviance Table (Type II tests)

Response: FISHD
                     Chisq      Df Pr(>Chisq)    
TR               0.2014     1     0.6536    
SESSION          25.0418    1   5.61e-07 ***
TR:SESSION       103.9113   1  < 2.2e-16 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

I am trying to decipher when there is a treatment affect. Any ideas?

IRTFM
  • 258,963
  • 21
  • 364
  • 487
Berin
  • 11
  • 1
  • 5

1 Answers1

0

At the moment you are modeling session as a continuous variable and so the reported coefficient for TR:SESSION is sort of a "test of trend" across the range 1:3 for that value of SESSION in TRT=2 .... is that what you were hoping to discover? Let's do it the right way, with factors and in particular create an ordered factor:

model2 = lme(FISHD ~ factor(TR)*ordered(SESSION), # same as TRT+SESSION+TRT*ORDERED
 random = ~1|MAC,
 correlation = corAR1(form = ~ SESSION | MAC,
 value = -0.7547232),
   data=dat,
  method="REML")
 anova(model2)
 #------
                            numDF denDF  F-value p-value
(Intercept)                     1     4 37.48966  0.0036
factor(TR)                      1     2  0.20124  0.6976
ordered(SESSION)                2     4 13.94379  0.0157
factor(TR):ordered(SESSION)     2     4 52.02864  0.0014

So this is saying that an analysis that takes into account the interaction of TRT with sequence of SESSIONs finds a significant difference for the effect of TRT. If you do not model the interaction (using formula: FISHD ~ factor(TR)+ordered(SESSION), then neither component is significant:

model1 = lme(FISHD ~  factor(TR)+ordered(SESSION), 
 random = ~1|MAC,
 correlation = corAR1(form = ~ SESSION | MAC, value = -0.7547232),
   data=dat,
  method="REML")
 anova(model1)
 #-------------
                 numDF denDF  F-value p-value
(Intercept)          1     6 44.53272  0.0005
factor(TR)           1     2  0.14143  0.7430
ordered(SESSION)     2     6  2.55154  0.1578

And:

>  anova(model1, model2)
       Model df      AIC      BIC    logLik   Test  L.Ratio p-value
model1     1  7 54.74813 55.30422 -20.37406                        
model2     2  9 44.04435 42.17018 -13.02217 1 vs 2 14.70378   6e-04
Warning message:
In anova.lme(model1, model2) :
  fitted objects with different fixed effects. REML comparisons are not meaningful.

That is a warning but I do not think it applies here since the varaialbes in the model are the same and only the structure is different.

IRTFM
  • 258,963
  • 21
  • 364
  • 487
  • I need some clarification. What is the warning against? What you have done is just for the repeated measures in general correct? We were not confident in our approach to begin with so this is helpful. So you are saying that when the interaction of Treatment on Time, or "SESSION" is accounted for, we get the significant values for time, and treatment within time, which means time and treatment at a certain session in time is significantly affecting the Fish Density. Correct? So we are trying to find a post hoc test to determine at which point in time treatment is significant. – Berin Jul 08 '17 at 19:42
  • That's not how I see it. The tests are all comparisons across times or treatments, so there no basis for saying one SESSION is responsible. If you do a plot, you see that the FISHD rises and then falls, so the "quadratic" effect over time is what I suspect is the basis for the significance of the model comparison. – IRTFM Jul 08 '17 at 20:14
  • Sorry, one more question. So you are saying there's no way to determine the affect of treatment within time? I have 3 distinct sessions where we sampled, one, two and three. On the third session, the treated macrocosms had higher fish densities so we expect the third session is where treatment has an affect. There's no post hoc test for testing after a repeated measures anova in R? I believe we are looking for an equivalent to "slice" it in SASS. Thank you so much for your help. – Berin Jul 08 '17 at 20:21
  • I didn't say that. I said the was no one session that could be pointed to with a statement that the effect was "there". (Looking at the data, I;m not "seeing" the pattern you describe.) This is not the right place to find people who have joint experience with SAS and R and such questions are off-topic, anyway, since they are more about theory and less about coding. The mixed effects gurus in R seem to think that SAS results and claims are not always reliable, but I don't have enough (zero in fact) experience with "slicing" whatever that might be to offer anything other than "do a search". – IRTFM Jul 08 '17 at 20:26
  • I am not familiar with SAS either, that is just what my advisor is comparing it to and I was hoping someone may know an equivalent because I don't know R or SAS, so any information is helpful. What I am trying to say is that there is a significant TR:SESSION effect (p=0.0014), so we are trying to figure out what post hoc test in R would be able to determine the TR:SESSION effect? – Berin Jul 08 '17 at 20:40
  • My perspective is that talking about a "TR:SESSION effect" (separate from a TR effect) is not really sensible. (So I've been contaminated by R experts and would get into arguments with your SAS-believing advisor.) You can talk about predictions for different scenarios and you can talk about a coefficient for TR:SESSION but predictions always require specification of particular TR and SESSION values and need to calculate the sum of the "main effect" and the "interaction effect". I'm not going to respond further. You might find discussions on this topic on CrossValidated.com – IRTFM Jul 08 '17 at 21:28