如何使用plm包比较R中的2个模型？

Question

我正在使用R中的plm包运行一个固定效应模型，我想知道如何比较哪个模型更适合。

例如，这是我构建的两个模型的代码：

library(plm)

eurofix <- plm(rlogmod ~ db+gdp+logvix+gb+i+logtdo+fx+ld+euro+core, 
               data=euro, 
               model="within")

eurofix2 <- plm(rlogmod ~ db+gdp+logvix+gb+i+logtdo+ld+euro+core, 
                data=euro,
                model="within")

我知道在普通的 lm 调用中，可以通过运行 anova 测试来比较两个模型，但是在这种情况下似乎不起作用。我总是得到以下错误：

Error in UseMethod("anova") : 
  no applicable method for 'anova' applied to an object of class "c('plm', 'panelmodel')"

有人知道如何使用plm包吗？Wald检验是否适用？

- NuclearPenguins

---

我最终决定只使用Wald检验。结果发现，这在plm包中运行良好。 - NuclearPenguins

谢谢，我正想问同样的问题，直到看到这个。你能分享一下实现Wald检验比较两个模型的代码吗？我遇到了些麻烦。 - daanoo

是的，在线性回归模型中的方差分析（ANOVA）等同于Wald检验。请参考此链接：https://stats.stackexchange.com/questions/131401/how-to-get-anova-table-with-robust-standard-errors。 - Juan Antonio Roldán Díaz

数据集包中的欧元数据（默认可用）没有支持这个问题的结构。 - IRTFM

在@42上，有一个关于交叉验证的问题，几乎相同，并包括数据和所有工作示例，可以在这里找到https://stats.stackexchange.com/questions/351746/comparing-groups-in-repeated-measures-fe-models-with-a-nested-error-component - Eric Fail

Answer 1

以下代码回答了Cross Validated中类似的问题。那里的问题也涉及到plm例程中的测试（联合）假设。将代码应用于您的问题应该很简单。

library(plm)  # Use plm
library(car)  # Use F-test in command linearHypothesis
library(tidyverse)
data(egsingle, package = 'mlmRev')
dta <- egsingle %>% mutate(Female = recode(female, .default = 0L, `Female` = 1L))
plm1 <- plm(math ~ Female * (year), data = dta, index = c('childid', 'year', 'schoolid'), model = 'within')

# Output from `summary(plm1)` --- I deleted a few lines to save space.
# Coefficients:
#                 Estimate Std. Error t-value Pr(>|t|)    
# year-1.5          0.8842     0.1008    8.77   <2e-16 ***
# year-0.5          1.8821     0.1007   18.70   <2e-16 ***
# year0.5           2.5626     0.1011   25.36   <2e-16 ***
# year1.5           3.1680     0.1016   31.18   <2e-16 ***
# year2.5           3.9841     0.1022   38.98   <2e-16 ***
# Female:year-1.5  -0.0918     0.1248   -0.74     0.46    
# Female:year-0.5  -0.0773     0.1246   -0.62     0.53    
# Female:year0.5   -0.0517     0.1255   -0.41     0.68    
# Female:year1.5   -0.1265     0.1265   -1.00     0.32    
# Female:year2.5   -0.1465     0.1275   -1.15     0.25    
# ---

xnames <- names(coef(plm1)) # a vector of all independent variables' names in 'plm1'
# Use 'grepl' to construct a vector of logic value that is TRUE if the variable
# name starts with 'Female:' at the beginning. This is generic, to pick up
# every variable that starts with 'year' at the beginning, just write
# 'grepl('^year+', xnames)'.
picked <- grepl('^Female:+', xnames)
linearHypothesis(plm1, xnames[picked])

# Hypothesis:
# Female:year - 1.5 = 0
# Female:year - 0.5 = 0
# Female:year0.5 = 0
# Female:year1.5 = 0
# Female:year2.5 = 0
# 
# Model 1: restricted model
# Model 2: math ~ Female * (year)
# 
#   Res.Df Df Chisq Pr(>Chisq)
# 1   5504                    
# 2   5499  5  6.15       0.29

- semibruin

---

你认为有可能以某种方式重写-c(1:5)块，使代码更加通用吗？（同样在这里） - Eric Fail

@EricFail 好建议。我用正则表达式匹配替换了计数(-c(1:5))。 - semibruin

Answer 2

我也曾经为此苦恼，但最终在一位博士朋友的帮助下找到了以下解决方案。使用你的示例，参考下面的样本解决方案。

使用AIC标准比较面板模型，具体如下：

library(plm)

eurofix <- plm(rlogmod ~ db+gdp+logvix+gb+i+logtdo+fx+ld+euro+core, 
               data=euro, 
               model="within")

eurofix2 <- plm(rlogmod ~ db+gdp+logvix+gb+i+logtdo+ld+euro+core, 
                data=euro,
                model="within")

# AIC = log(RSS/N) + 2K/N  for linear models
# AIC = log(RSS/n) + 2K/n  for panel models

Sum1 <- summary(eurofix)
RSS1 <- sum(Sum1$residuals^2)
K1 <- max(eurofix$assign)
N1 <- length(eurofix$residuals)
n1 <- N1 - K1 - eurofix$df.residual

AIC_eurofix = log(RSS1/n1) + (2*K1)/n1

Sum2 <- summary(eurofix2)
RSS2 <- sum(Sum2$residuals^2)
K2 <- max(eurofix2$assign)
N2 <- length(eurofix2$residuals)
n2 <- N2 - K2 - eurofix2$df.residual

AIC_eurofix2 = log(RSS2/n2) + (2*K2)/n2

较小的AIC值是首选模型！

- Glen Wade

Answer 3

你有使用过 plm 函数中的 anova() 吗？根据你提出问题的方式，我无法确定。如果没有使用过，可以试一下。

如果你的问题更多是关于选择一种方法来帮助你在模型之间进行裁决而不是技术性的问题，那么答案实际上取决于你如何定义“合适”。如果两个模型唯一的区别是第一个模型中包含了 fx，那么可以使用几个统计测试来评估你的模型最小化平方误差的程度（例如 R^2）或由于残差的非随机分布而导致的缺陷（例如 VIF）。

如果您想知道是否包括fx会产生一个相对于过度拟合而言更适合您的数据的模型，可以考虑使用BIC。我通常更喜欢BIC，因为它比其他模型拟合统计量（如AIC）更积极地惩罚额外的参数。具有最低BIC的模型往往是最好的拟合模型（尽管您还应该使用Wald测试/F测试进行确认，特别是当嵌套模型是其理想用例时）。您应该能够使用plm获取模型对象上的BIC值，如下所示：

anova(model1, model2)

如果那不起作用，我发现lme4包的函数非常有用：

BIC(model1, model2)

如果我误解了问题，请告诉我 - 并让我们知道你的发现！

- J.Q