1

这是我第一次问关于 R 的问题,如果我遗漏了必要的细节,请原谅我。我会尽量做到彻底但简洁。

我的数据有 4 列:分数(值 -1、0 或 1);物种(至少 10 种不同的);栖息地(三个不同的);和站点(至少 8 个不同的站点)。所以每一行都是分数、物种、栖息地和地点的独特组合。请参阅下面的顶行:

score species habitat  site
0    NOCA     NAG  NAG1
0    BRTH     NAG  NAG1
0    BARS     NAG  NAG1
1    COYE     NAG  NAG1
0    HOWR     NAG  NAG1
0    SAVS     NAG  NAG1
0    CEDW     NAG  NAG1
1    CHSP     NAG  NAG1
0    EAKI     NAG  NAG1
1    MODO     NAG  NAG1
0    NOCA     NAG NAG16
0    BRTH     NAG NAG16

我运行了一个 lme 来评估物种和栖息地对分数的影响,地点是随机效应。代码如下:

anova(model<-lme(score~species*habitat, random=~1|site))

然后,我想进行成对比较,以确定显着差异在哪里。TukeyHSD 不适用于 lme,所以我使用了 multcomp 包,如下所示:

summary(glht(model, linfct=mcp(habitat="Tukey")))

据我所知,这将根据栖息地类型对分数进行成对比较。我还尝试用物种代替栖息地,以根据物种进行比较。在这两种情况下,我都收到以下错误:

contrMat(table(mf[[nm]]), type = types[pm]) 中的错误:少于两组

我看不出我怎么能在任何地方拥有少于两个组,因为我对任何列的最小组数是三个。也许我误解了 R 所指的“组”?关于可能出现的问题以及如何解决问题的任何建议?

编辑根据建议(谢谢),这里有一些关于我的数据的更多信息,可能会使其可重现。首先,我完整使用的代码。

library(nlme)
library(multcomp)
null<-read.csv("baci_null_red.csv", head=TRUE)
attach(null)
anova(model<-lme(score~species*habitat, random=~1|site))
              numDF   denDF  F-value  p-value

(Intercept)         1   1392   1.929021  0.1651

species            29   1392   2.691207    <.0001

habitat             2    48    3.412485  0.0411

species:habitat    58  1392    1.239267  0.1099

摘要(glht(模型, linfct=mcp(habitat="Tukey")))

Error in contrMat(table(mf[[nm]]), type = types[pm]) : 
  less than two groups

以下是一些关于可能有助于重现它的数据的摘要……它有点大,我不确定如何使它最小化,但仍然会出现错误。

str(null)

'data.frame':   1530 obs. of  4 variables:
 $ score  : int  0 0 -1 0 0 0 0 0 0 0 ...
 $ species: Factor w/ 30 levels "AMGO","AMRO",..: 27 5 7 2 18 12 22 28 6 13 ...
 $ habitat: Factor w/ 3 levels "CUM","IAG","NAG": 3 3 3 3 3 3 3 3 3 3 ...
 $ site   : Factor w/ 51 levels " CUM6","CUM1",..: 37 37 37 37 37 37 37 37 37 37 ..

summary(null)

 score             species     habitat        site     

 Min.   :-1.00000   AMGO   :  51   CUM:210    CUM6  :  30  
 1st Qu.: 0.00000   AMRO   :  51   IAG:660   CUM1   :  30  
 Median : 0.00000   BARS   :  51   NAG:660   CUM10  :  30  
 Mean   :-0.01569   BCCH   :  51             CUM12  :  30  
 3rd Qu.: 0.00000   BHCO   :  51             CUM2   :  30  
 Max.   : 1.00000   BLJA   :  51             CUM3   :  30  
                (Other):1224             (Other):1350  

dput(null[somerows,c("species","habitat","site")])


structure(list(species = structure(c(27L, 5L, 7L, 2L, 18L, 12L, 
22L, 28L, 6L, 13L, 15L, 23L, 4L, 10L, 14L, 11L, 20L, 30L, 9L, 
19L, 26L, 8L, 16L, 25L, 3L, 17L, 21L, 29L, 24L, 1L, 27L, 5L, 
7L, 2L, 18L, 12L, 22L, 28L, 6L, 13L, 15L, 23L, 4L, 10L, 14L, 
11L, 20L, 30L, 9L, 19L, 26L, 8L, 16L, 25L, 3L, 17L, 21L, 29L, 
24L, 1L, 27L, 5L, 7L, 2L, 18L, 12L, 22L, 28L, 6L, 13L, 15L, 23L, 
4L, 10L, 14L, 11L, 20L, 30L, 9L, 19L, 26L, 8L, 16L, 25L, 3L, 
17L, 21L, 29L, 24L, 1L, 27L, 5L, 7L, 2L, 18L, 12L, 22L, 28L, 
6L, 13L), .Label = c("AMGO", "AMRO", "BARS", "BCCH", "BHCO", 
"BLJA", "BOBO", "BRTH", "CEDW", "CSWA", "EABL", "EAKI", "EAWP", 
"FISP", "GCFL", "HOLA", "HOWR", "INBU", "KILL", "LEFL", "MODO", 
"NOCA", "RBGR", "RWBL", "SAVS", "SOSP", "VESP", "WIFL", "YSFL", 
"YWAR"), class = "factor"), habitat = structure(c(3L, 3L, 3L, 
3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 
3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 
3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 
3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 
3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 
3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 
3L), .Label = c("CUM", "IAG", "NAG"), class = "factor"), site = structure(c(37L, 
37L, 37L, 37L, 37L, 37L, 37L, 37L, 37L, 37L, 37L, 37L, 37L, 37L, 
37L, 37L, 37L, 37L, 37L, 37L, 37L, 37L, 37L, 37L, 37L, 37L, 37L, 
37L, 37L, 37L, 46L, 46L, 46L, 46L, 46L, 46L, 46L, 46L, 46L, 46L, 
46L, 46L, 46L, 46L, 46L, 46L, 46L, 46L, 46L, 46L, 46L, 46L, 46L, 
46L, 46L, 46L, 46L, 46L, 46L, 46L, 47L, 47L, 47L, 47L, 47L, 47L, 
47L, 47L, 47L, 47L, 47L, 47L, 47L, 47L, 47L, 47L, 47L, 47L, 47L, 
47L, 47L, 47L, 47L, 47L, 47L, 47L, 47L, 47L, 47L, 47L, 48L, 48L, 
48L, 48L, 48L, 48L, 48L, 48L, 48L, 48L), .Label = c(" CUM6", 
"CUM1", "CUM10", "CUM12", "CUM2", "CUM3", "CUM8", "IAG1", "IAG10", 
"IAG13", "IAG14", "IAG15", "IAG16", "IAG18", "IAG19", "IAG21", 
"IAG22", "IAG23", "IAG24", "IAG25", "IAG26", "IAG27", "IAG28", 
"IAG3", "IAG4", "IAG5", "IAG6", "IAG8", "IAG9", "NAG10", "NAG11", 
"NAG13", "NAG14", "NAG15", "NAG18", "NAG19", "NAG2", "NAG21", 
"NAG22", "NAG23", "NAG24", "NAG25", "NAG26", "NAG27", "NAG28", 
"NAG3", "NAG4", "NAG5", "NAG6", "NAG7", "NAG8"), class = "factor")), .Names =
  c("species", "habitat", "site"), row.names = c(NA, 100L), class = "data.frame")

你有它:s 感谢您的任何指导。

4

1 回答 1

0

我刚刚通过避免附加数据框解决了这个问题。我在网上读到,当尝试将拟合模型传递给函数时,附加数据框会产生问题。

我建议使用这样的代码而不是使用attach()

summary(model<-lme(fix~treat*week, random=~1|rep/week, data=dataframename))
summary(glht(model, linfct=mcp(treat="Tukey")))
于 2016-03-12T00:45:34.923 回答