R中绘制预测值并比较随机效应plm模型中的交互项的最佳方法是什么？答案

【问题标题】：What is the best way in R to plot predicted values and compare interaction terms in a random effects plm model?R中绘制预测值并比较随机效应plm模型中的交互项的最佳方法是什么？
【发布时间】：2019-08-13 18:34:07
【问题描述】：

我有一个 plm 随机效应模型，其中包括连续变量和分类变量之间的交互作用（为了说明目的，假设分类变量是强烈的宗教信仰，而连续变量是对慈善事业的贡献）。

这是示例数据（IRL，数据要大得多）

ratio     contrib     relig    np_score     ID    year 
.4          3          1          11        1      1990  
0           7          0          8         2      1990
.9          7          1          6         1      1992
.7          6          1          10        1      1994
.1          2          0          4         2      1992
.3          9          0          8         2      1994

我想绘制预测结果，因为贡献增加取决于受访者是否认同强烈的宗教信仰。

m1 <- plm(ratio ~ contrib + relig + np_score + contrib*relig, 
            data = panel_dat, 
            index = c("ID", "year"),
            model = "random", random.method = "amemiya")

我尝试了interaction.plot如下：

interaction.plot(panel_dat$contrib, panel_dat$relig, predict(m1), col = 2:3, lty = 1)

但我得到了错误术语：

Error in tapply(response, list(x.factor, trace.factor), fun) : 
  arguments must have same length

我最终宁愿用 ggplot2 来做这件事。有什么建议？似乎应该有一个简单的答案...

【问题讨论】：

See here 提出一个人们可以帮助解决的 R 问题。这包括数据样本和所有必要的代码。
还可以查看 B. Bolker 的 glmmFAQ 以及可能的 DHARMa 包。后者目前缺乏科学论证（尽管看起来有些希望）。

标签： r plot interaction plm

【解决方案1】：

此示例代码生成下面的多线图（在 ggplot 中）。不确定轴是否正确...您可以根据需要进行调整。

# Create sample dataframe
df <- data.frame(
  id = c(1, 1, 1, 2, 2, 2),
  year = c(2000, 2001, 2002, 2000, 2001, 2002),  
  ratio = c(0.3, 0.6, 0.1, 0.4, 0.5, 0.6),
  contrib = c(310, 220, 230, 140, 0, 10),
  relig = c(1, 1, 1, 3, 3, 3)
)

# Build model and predictions
library(plm)
m1 <- plm(ratio ~ contrib + relig + contrib*relig, data = df, index = c('id','year'))
df$p_ratio <- predict(m1)

# Create plot
library(ggplot2)
ggplot(data=df, aes(x=contrib, y=p_ratio, group=relig))+geom_line(size=2, aes(color=relig))

【讨论】：