1
$\begingroup$

thanks in advance for the help. I've been getting different values for the mean squares (and as a result, F-values) using the base R lmversus lmer. Here's the data (sorry, but I believe all of it needs to be there to reproduce the problem):

tmean <- structure(list(site = structure(c(1L, 1L, 1L, 1L, 1L, 1L, 1L,1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 2L,2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 3L, 3L, 3L, 3L, 3L, 3L, 3L,3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 4L, 4L, 4L, 4L, 4L, 4L, 4L,4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 5L, 5L, 5L, 5L, 5L, 5L, 5L,5L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 6L, 6L, 6L, 6L, 6L, 6L, 6L,6L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 7L, 7L, 7L, 7L, 7L, 7L, 7L,7L, 7L, 7L, 7L, 7L, 7L, 7L, 7L, 7L, 8L, 8L, 8L, 8L, 8L, 8L, 8L,8L, 8L, 8L, 8L, 8L, 8L, 8L, 8L, 8L), .Label = c("Batiquitos","Cardiff", "Imperial Beach", "Moonlight", "North Carlsbad", "Oceanside","Solana Beach", "South Carlsbad"), class = "factor"),
     TSD = c(1,1, 1, 1, 4, 4, 4, 4, 12, 12, 12, 12, 15, 15, 15, 15, 1, 1, 1,1, 4, 4, 4, 4, 12, 12, 12, 12, 15, 15, 15, 15, 1, 1, 1, 1, 4,4, 4, 4, 12, 12, 12, 12, 15, 15, 15, 15, 1, 1, 1, 1, 4, 4, 4,4, 12, 12, 12, 12, 15, 15, 15, 15, 1, 1, 1, 1, 4, 4, 4, 4, 12,12, 12, 12, 15, 15, 15, 15, 1, 1, 1, 1, 4, 4, 4, 4, 12, 12, 12,12, 15, 15, 15, 15, 1, 1, 1, 1, 4, 4, 4, 4, 12, 12, 12, 12, 15,15, 15, 15, 1, 1, 1, 1, 4, 4, 4, 4, 12, 12, 12, 12, 15, 15, 15,15),
     sand = structure(c(1L, 1L, 2L, 2L, 1L, 1L, 2L, 2L, 1L, 1L,2L, 2L, 1L, 1L, 2L, 2L, 1L, 1L, 2L, 2L, 1L, 1L, 2L, 2L, 1L, 1L,2L, 2L, 1L, 1L, 2L, 2L, 1L, 1L, 2L, 2L, 1L, 1L, 2L, 2L, 1L, 1L,2L, 2L, 1L, 1L, 2L, 2L, 1L, 1L, 2L, 2L, 1L, 1L, 2L, 2L, 1L, 1L,2L, 2L, 1L, 1L, 2L, 2L, 1L, 1L, 2L, 2L, 1L, 1L, 2L, 2L, 1L, 1L,2L, 2L, 1L, 1L, 2L, 2L, 1L, 1L, 2L, 2L, 1L, 1L, 2L,
     2L, 1L, 1L,2L, 2L, 1L, 1L, 2L, 2L, 1L, 1L, 2L, 2L, 1L, 1L, 2L, 2L, 1L, 1L,2L, 2L, 1L, 1L, 2L, 2L, 1L, 1L, 2L, 2L, 1L, 1L, 2L, 2L, 1L, 1L,2L, 2L, 1L, 1L, 2L, 2L), .Label = c("control", "nourished"), class = "factor"),Total = c(0.722222222222222, 0.5625, 0.222222222222222, 0.222222222222222,0.0526315789473684, 0.55, 1.53846153846154, 2.53846153846154,0.923076923076923,
     1.85714285714286, 1.41176470588235, 1.23529411764706,1, 0.8125, 0.705882352941177, 0.647058823529412, 1.17391304347826,1.875, 0.285714285714286, 0.1, 0.277777777777778, 0.266666666666667,0.0714285714285714, 0.133333333333333, 0.846153846153846,1.15384615384615, 0.952380952380952, 1.85714285714286, 2.85,2.5, 0.625, 1.26923076923077, 1.04545454545455,
     1.23809523809524,0.1, 0.0526315789473684, 0.5, 2.375, 0, 0.0769230769230769,2.05263157894737, 2, 0.294117647058824, 0.529411764705882,3.78947368421053, 3.42105263157895, 0.466666666666667, 0.611111111111111,1.44444444444444, 1, 0.523809523809524, 0.0526315789473684,1.6875, 0.6875, 0.294117647058824, 0.0625, 2.65, 2.5, 1.55,1.2, 3.16666666666667, 3.68421052631579,
     1.53333333333333,0.625, 0.1875, 0.6875, 0.2, 0.2, 1.05882352941176, 3.94117647058824,5.28571428571429, 5.26315789473684, 0.3125, 0.0625, 0.176470588235294,0.235294117647059, 0.15, 0.0555555555555556, 0.210526315789474,0.1875, 1.46666666666667, 1.78571428571429, 0.142857142857143,0.0714285714285714, 0.571428571428571, 0.428571428571429,0.153846153846154,
     0.142857142857143, 0.357142857142857,0.5, 0.636363636363636, 0.384615384615385, 0.0526315789473684,0.684210526315789, 0.611111111111111, 0.388888888888889,0.933333333333333, 3.85714285714286, 0, 0, 3.92307692307692,4.75, 6, 5, 3.75, 2.25, 2.52941176470588, 2.88235294117647,0.944444444444444, 2, 1.69230769230769, 0.769230769230769,0.333333333333333, 0.466666666666667,
     0, 0, 1.08333333333333,0.384615384615385, 0.916666666666667, 3.91666666666667, 1.23076923076923,0.785714285714286, 0.928571428571429, 0.538461538461538,0.470588235294118, 0.235294117647059, 1, 1.04761904761905)), .Names = c("site", "TSD", "sand", "Total"), row.names = c(NA,128L), class = "data.frame")

The packages I've been using are:

require(lme4)
require(lmerTest)

And I'm interested in the effects of three factors: Sand (fixed, categorical), TSD (fixed, quantitative), and Site (random, categorical). I want to know the mean squares for all three factors and their interactions, so I use anova(lm(Total~TSD*sand+site,data=tmean)) to obtain this:

Analysis of Variance Table

Response: Total
           Df  Sum Sq Mean Sq F value    Pr(>F)    
TSD         1   1.090  1.0901  0.7514 0.3877938    
sand        1   5.723  5.7232  3.9450 0.0493478 *  
site        7  42.427  6.0610  4.1779 0.0003856 ***
TSD:sand    1   0.324  0.3242  0.2235 0.6372812    
Residuals 117 169.735  1.4507                      
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Now that I know all the MS values, I can create the actual model I want (appropriately treating "site" as a random factor) with lmer and see if the results make sense:

# Model using lme4

> n<-lmer(Total~TSD*sand+(1|site),data=tmean)
> anova(n)
Analysis of Variance Table of type III  with  Satterthwaite 
approximation for degrees of freedom
          Sum Sq Mean Sq NumDF DenDF F.value Pr(>F)
TSD      1.09014 1.09014     1   117 0.75144 0.3878
sand     0.85496 0.85496     1   117 0.58933 0.4442
TSD:sand 0.32421 0.32421     1   117 0.22348 0.6373

And that's where things get confusing. The Mean Sq for TSD and TSD:sand match exactly what the earlier lmoutput showed, but now the Mean Sq for sand is 0.85496, versus the 5.7232 returned earlier. Clearly, this significantly reduces the resulting F value when I divide by the residuals 1.4507, resulting in a very high p-value (oh no!).

Can anyone shed light on why this may be happening? I would assume that regardless of the model or package used, the Mean squares calculated for a given factor should always be the same (while F-values vary depending on what's random and p-values are subject to method used). Or am I wrong in this?

Thank you!

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Browse other questions tagged or ask your own question.