I had a code that worked fine so far. I want to test things with gls
, lme
and gamm
(from packages nlme and mgcv), and I compared different models with anova()
. However, I needed another package, that did not work with my R version (which was almost one year old). Thus, I updated R (via the updater
package) and RStudio.
The issue now is, that anova()
does not give any output after running or only "Denom. DF: 91" and nothing else.
Now I tried different things and searched a lot, but I found no current threat dealing with such a problem, while looking at the help files just says, it should work that way I use it. Thus, I am suspecting that I miss something essential (probably even obvious), but I don't get it. I hope you can tell me where I do something wrong.
Here is some data to play with (copied from txt-file):
"treat" "x" "time" "nest"
"1" "1" 49.37 1 "K1"
"2" "1" 48.68 1 "K2"
"3" "2" 44.7 1 "T7"
"4" "2" 49.3 1 "T8"
"5" "1" 48.78 1 "K3"
"6" "2" 42.37 1 "T10"
"7" "1" 39.26 1 "K4"
"8" "2" 46.36 1 "T11"
"9" "1" 40.36 1 "K5"
"10" "2" 47.14 1 "T9"
"11" "1" 48.81 1 "K6"
"12" "1" 40.4 1 "K10"
"13" "2" 53.42 1 "T4"
"14" "2" 46.85 1 "T5"
"15" "2" 44.58 1 "T2"
"16" "2" 47.51 1 "T6"
"17" "1" 51.7 1 "K8"
"18" "1" 48.16 1 "K7"
"19" "2" 48.86 1 "T3"
"20" "1" 44.6 1 "K11"
"21" "1" 49.71 1 "K9"
"22" "2" 44.54 1 "T1"
"23" "2" 41.55 2 "T3"
"24" "1" 32.55 2 "K3"
"25" "1" 42.15 2 "K1"
"26" "2" 51.06 2 "T1"
"27" "1" 38.43 2 "K11"
"28" "2" 39.91 2 "T11"
"29" "1" 36.73 2 "K7"
"30" "2" 50.19 2 "T4"
"31" "1" 42.26 2 "K8"
"32" "1" 43.02 2 "K6"
"33" "2" 37.6 2 "T10"
"34" "1" 33.42 2 "K4"
"35" "2" 39.64 2 "T5"
"36" "2" 43.56 2 "T2"
"37" "2" 35.31 2 "T7"
"38" "2" 37 2 "T8"
"39" "2" 40.87 2 "T6"
"40" "1" 35.29 2 "K9"
"41" "2" 41.83 2 "T9"
"42" "1" 37.88 2 "K10"
"43" "1" 36.5 2 "K5"
"44" "1" 34.21 3 "K4"
"45" "1" 38.04 3 "K6"
"46" "1" 35.14 3 "K3"
"47" "2" 38.18 3 "T10"
"48" "1" 40.26 3 "K11"
"49" "2" 37.09 3 "T3"
"50" "2" 43.1 3 "T11"
"51" "2" 34.26 3 "T7"
"52" "1" 36.58 3 "K9"
"53" "1" 35.81 3 "K2"
"54" "1" 39.83 3 "K10"
"55" "2" 37.65 3 "T6"
"56" "1" 39.8 3 "K7"
"57" "1" 36.41 3 "K8"
"58" "1" 35.22 3 "K5"
"59" "2" 39.68 3 "T8"
"60" "2" 41.12 3 "T1"
"61" "2" 36.93 3 "T9"
"62" "1" 35.66 3 "K1"
"63" "2" 36.91 3 "T4"
"64" "2" 38.84 3 "T5"
"65" "2" 34.31 3 "T2"
"66" "1" 32.71 4 "K9"
"67" "2" 37.84 4 "T11"
"68" "1" 28.01 4 "K10"
"69" "2" 39.69 5 "T11"
"70" "2" 35.08 4 "T10"
"71" "2" 34.43 4 "T9"
"72" "1" 32.12 4 "T8"
"73" "2" 30.41 4 "T7"
"74" "1" 31.81 4 "K7"
"75" "2" 36.41 4 "T6"
"76" "1" 29.17 5 "K6"
"77" "1" 28.59 4 "K6"
"78" "2" 33.99 4 "T5"
"79" "1" 30.41 4 "K5"
"80" "1" 29.8 4 "K4"
"81" "2" 34.72 4 "T4"
"82" "2" 34.38 4 "T3"
"83" "1" 28.12 4 "K3"
"84" "2" 34.62 4 "T2"
"85" "1" 31.88 4 "K2"
"86" "1" 29.35 4 "K1"
"87" "2" 37.95 4 "T1"
"88" "2" 40.85 5 "T4"
"89" "2" 35.07 5 "T5"
"90" "2" 36.15 5 "T8"
"91" "2" 36.48 5 "T10"
"92" "1" 33.73 4 "K8"
"93" "1" 28.17 5 "K9"
"94" "1" 32.81 5 "K10"
"95" "1" 32.17 4 "K11"
And this is basically one of the models I try to run:
test <- read.table(file="C:/Users/marvi_000/Desktop/testdata.txt")
str(test)
test$treat <- as.factor(test$treat)
test$nest <- as.factor(test$nest)
library(nlme)
m.test <- gls(x ~ treat * time,
correlation = corAR1(form =~ time | nest),
test, na.action = na.omit)
anova(m.test)
the output is:
Denom. DF: 91
When comparing models with anova(m1, m2)
nothing happens at all.
The same is true when I run a gamm
from package mgcv
and using anova(m$lme)
or anova(m1$lme, m2$lme)
.
I would appreciate any help or hint, pointing me towards the right direction. Thanks a lot!
EDIT:
After some discussion, I found out, that it is a problem with the scripts. I'm using RStudio and RMarkdown. However, when I run the code (with cntrl+enter, line by line) within the markdown script, the anova(lmemodel)
command does not work as supposed to. However, if I just copy this single command into a plane r script (still using the current environment), the command is executed properly showing the desired output.
I have no clue what is happening there. If anybody has an idea where the problem is, or how to solve it, I would still be happy to hear it.