R MuMIn model.avg() - 相对重要性不适用于粘贴的公式

Question

我正在尝试将 MuMIn 的 model.avg 函数与粘贴的模型公式一起使用，并使用索引而不是直接输入，例如：

m1<-gls(as.formula(paste(response,"~",paste(combns[,j], collapse="+"))), data=dat)

'combns' 是由 combn() 创建的包含预测变量组合的二维数组。这产生的模型平均系数和 AICc 值与 gls 函数直接包含公式时产生的相同，例如：

m1<-gls(median_Ta ~ day_of_season + hour_of_day + pct_grey_cover + 
    foliage_height_diversity + tree_shannon_diversity + median_patch_size, data=dat)

然而，相对变量的重要性不是计算，我相信这与使用 for 循环或使用变量访问列表的索引有关，模型存储在列表中以某种方式导致组件模型术语不正确 'read'（参见模型的术语代码）：

Component models: 
         df  logLik   AICc delta weight
1234567b  7 -233.08 481.43  0.00   0.59
1234567f  3 -237.97 482.21  0.78   0.40
1234567e  4 -241.32 491.08  9.65   0.00
1234567a  9 -241.15 502.39 20.96   0.00
1234567c  6 -248.37 509.68 28.25   0.00
1234567d  5 -250.22 511.11 29.68   0.00

Term codes: 
           day_of_season foliage_height_diversity              hour_of_day 
                       1                        2                        3 
       median_patch_size           pct_grey_cover   tree_shannon_diversity 
                       4                        5                        6 
 urban_boundary_distance 
                       7

这导致相对变量重要性被给出为：

Relative variable importance: 
                     day_of_season foliage_height_diversity hour_of_day
Importance:          1             1                        1          
N containing models: 6             6                        6          
                     median_patch_size pct_grey_cover     tree_shannon_diversity
Importance:          1                 1              1                     
N containing models: 6                 6              6                     
                     urban_boundary_distance
Importance:          1                      
N containing models: 6

而如果我在相同模型上使用 model.avg 并单独输入公式，我会得到以下正确输出：

Component models: 
        df  logLik   AICc delta weight
23456    7 -233.08 481.43  0.00   0.59
1        3 -237.97 482.21  0.78   0.40
57       4 -241.32 491.08  9.65   0.00
1234567  9 -241.15 502.39 20.96   0.00
1467     6 -248.37 509.68 28.25   0.00
147      5 -250.22 511.11 29.68   0.00

Relative variable importance: 
                     pct_grey_cover median_patch_size tree_shannon_diversity
Importance:            0.6           0.59              0.59                 
N containing models:     3              4                 3                 
                      foliage_height_diversity hour_of_day day_of_season
Importance:           0.59                     0.59         0.4        
N containing models:     2                        2           4        
                     urban_boundary_distance
Importance:          <0.01                  
N containing models:     4

如何让 model.avg 正确读取公式中的预测变量？我在这里仅包含六个模型作为示例，但我想比较完整的 128 个模型集（并且我有其他具有更多预测变量的响应变量），因此单独输入它们是不可行的。

提前致谢。

编辑：可重现的示例

我花了一些时间来缩小问题范围。第一个示例 m.ave 显示了 for 循环中的问题。第二个示例 m.ave2 显示它使用键入的索引而不是使用变量。显然这只是预测变量的一小部分。

require(nlme)
require(MuMIn)

dat<-data.frame(median_Ta=rnorm(100), day_of_season=runif(100), hour_of_day=runif(100), pct_grey_cover=rnorm(100),
        foliage_height_diversity=rnorm(100), urban_boundary_distance=runif(100), tree_shannon_diversity=rnorm(100), 
        median_patch_size=rnorm(100))

f1<-"median_Ta ~ day_of_season + hour_of_day + pct_grey_cover + foliage_height_diversity + 
        urban_boundary_distance + tree_shannon_diversity + median_patch_size"

f1<-gsub("\s", "", f1) # remove whitespace
f1split <- strsplit(f1, split="~") # split predictors and response
response <- f1split[[1]][1]
predictors <- strsplit(f1split[[1]][2], split="+", fixed=TRUE)[[1]]

modelslist<-list()

combns <- combn(predictors, 6)
for (j in 1:7) {
    modelslist[[j]]<-gls(as.formula(paste(response,"~",paste(combns[,j], collapse="+"))), data=dat)
}

m.ave<-model.avg(modelslist[[2]], modelslist[[3]], modelslist[[4]],
        modelslist[[5]], modelslist[[6]], modelslist[[7]], modelslist[[8]])

summary(m.ave)

#compare....

modelslist2<-list()
modelslist2[[1]]<-gls(as.formula(paste(response,"~",paste(combns[,1], collapse="+"))), data=dat)
modelslist2[[2]]<-gls(as.formula(paste(response,"~",paste(combns[,2], collapse="+"))), data=dat)
modelslist2[[3]]<-gls(as.formula(paste(response,"~",paste(combns[,3], collapse="+"))), data=dat)
modelslist2[[4]]<-gls(as.formula(paste(response,"~",paste(combns[,4], collapse="+"))), data=dat)
modelslist2[[5]]<-gls(as.formula(paste(response,"~",paste(combns[,5], collapse="+"))), data=dat)
modelslist2[[6]]<-gls(as.formula(paste(response,"~",paste(combns[,6], collapse="+"))), data=dat)
modelslist2[[7]]<-gls(as.formula(paste(response,"~",paste(combns[,7], collapse="+"))), data=dat)

m.ave2<-model.avg(modelslist2[[1]], modelslist2[[2]], modelslist2[[3]], modelslist2[[4]],
        modelslist2[[5]], modelslist2[[6]], modelslist2[[7]])

summary(m.ave2)

Answer 1

这是 gls（在包 nlme 中）的 formula 方法中的错误。由于实际公式未存储在对象中的任何位置，因此它计算函数调用中的 "model" 参数。在 modellist 的元素的情况下，它们都是相同的，例如：

modelslist[[1]]$call$model
modelslist[[7]]$call$model

两个return

> formula(paste(response, "~", paste(combns[, j], collapse = "+")))

其中，当 eval 使用 j 的当前（最后）值时，所有 formula(modellist[[N]]) return 都是最后一个模型公式。

 all.equal(formula(modelslist[[1]]), formula(modelslist[[7]]))

returns

> TRUE

不用说，所有这些都让 model.avg 感到困惑，它使用公式来构建模型选择 table（这是后备方案，因为 gls 也缺少 terms ).

编辑：可能的解决方法

获得所需内容的更简单方法：

model.avg(dredge(..., m.lim = c(6,6)))

或者，如果您想进行预测：

modellist <- lapply(dredge(..., m.lim = c(6,6), evaluate = FALSE), eval)

但是，如果您想使用任意一组模型，请将每个 gls 模型对象中的 $call$model 元素替换为适当的公式，例如

combns <- combn(1:7, 6)
modellist <- vector("list", 7)
for (j in 1:7) {
    f <- reformulate(predictors[combns[, j]], response = response)
    fm <- gls(f, data = dat)
    fm$call$model <- f # assign the actual formula
    modellist[[j]] <- fm
}

R MuMIn model.avg() - 相对重要性不适用于粘贴的公式

R MuMIn model.avg() - relative importance not working with pasted formulae

r

mumin