0

我正在尝试测试两个边际效应之间的差异。我可以让 R 来计算效果,但我找不到任何资源来解释如何测试它们的差异。

我查看了边距文档和其他边际效果包,但无法找到可以测试差异的东西。

data("mtcars")

mod<-lm(mpg~as.factor(am)*disp,data=mtcars)

(marg<-margins(model = mod,at = list(am = c("0","1"))))
 at(am)     disp    am1
      0 -0.02758 0.4518
      1 -0.05904 0.4518
summary(marg)
 factor     am     AME     SE       z      p   lower   upper
    am1 1.0000  0.4518 1.3915  0.3247 0.7454 -2.2755  3.1791
    am1 2.0000  0.4518 1.3915  0.3247 0.7454 -2.2755  3.1791
   disp 1.0000 -0.0276 0.0062 -4.4354 0.0000 -0.0398 -0.0154
   disp 2.0000 -0.0590 0.0096 -6.1353 0.0000 -0.0779 -0.0402

我想进行一个测试来确定每一行 marg 中的边际效应是否显着不同;即,边际效应图中的斜率不同。这似乎是正确的,因为置信区间不重叠——表明位移的影响对于 am=0 和 am=1 是不同的。

我们在下面的评论中讨论我们可以使用 emmeans 测试对比,但这是对 am=0 和 am=1 的平均响应的测试。

emm<-emmeans(mod,~ as.factor(am)*disp)
emm
 am disp emmean    SE df lower.CL upper.CL
  0  231   18.8 0.763 28     17.2     20.4
  1  231   19.2 1.164 28     16.9     21.6
cont<-contrast(emm,list(`(0-1)`=c(1,-1)))
cont
 contrast estimate   SE df t.ratio p.value
 (0-1)      -0.452 1.39 28 -0.325  0.7479 

这里 p 值很大,表明当 am=0 时的平均响应与 am=1 时的平均响应没有显着差异。

这样做是否合理(例如测试两种方法的差异)?

smarg<-summary(marg)
(z=as.numeric((smarg$AME[3]-smarg$AME[4])/sqrt(smarg$SE[3]^2+smarg$SE[4]^2)))
[1] 2.745
2*pnorm(-abs(z))
[1] 0.006044

该 p 值似乎与非重叠置信区间的分析一致。

4

3 回答 3

3

如果我理解您的问题,可以使用以下方法回答emtrends

library(emmeans)

emt = emtrends(mod, "am", var = "disp")

emt  # display the estimated slopes

## am disp.trend      SE df lower.CL upper.CL
##  0    -0.0276 0.00622 28  -0.0403  -0.0148
##  1    -0.0590 0.00962 28  -0.0787  -0.0393
##
## Confidence level used: 0.95 

pairs(emt) # test the difference of slopes

## contrast estimate     SE df t.ratio p.value
## 0 - 1      0.0315 0.0115 28 2.745   0.0104
于 2019-06-08T02:07:42.903 回答
1

对于“斜率在统计上是否不同,表明位移的影响对于 am=0 与 am=1 不同?”的问题。问题,您可以直接从 lm() 拟合的摘要中获取与比较相关的 p 值。


> summary(mod)

Call:
lm(formula = mpg ~ as.factor(am) * disp, data = mtcars)

Residuals:
    Min      1Q  Median      3Q     Max 
-4.6056 -2.1022 -0.8681  2.2894  5.2315 

Coefficients:
                     Estimate Std. Error t value Pr(>|t|)    
(Intercept)         25.157064   1.925053  13.068 1.94e-13 ***
as.factor(am)1       7.709073   2.502677   3.080  0.00460 ** 
disp                -0.027584   0.006219  -4.435  0.00013 ***
as.factor(am)1:disp -0.031455   0.011457  -2.745  0.01044 *  
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 2.907 on 28 degrees of freedom
Multiple R-squared:  0.7899,    Adjusted R-squared:  0.7674 
F-statistic: 35.09 on 3 and 28 DF,  p-value: 1.27e-09

请注意,该as.factor(am)1:disp术语的 p 值为 0.01044,与 Russ Lenth 答案中的输出相匹配pairs(emt)

(作为答案发布,因为声誉不足以作为评论发布)

于 2020-10-24T18:28:04.973 回答
0

不确定,但您可能正在查看边际效应的对比或成对比较?您可以使用emmeans包执行此操作:

library(margins)
library(emmeans)
library(magrittr)

data("mtcars")
mod <- lm(mpg ~ as.factor(am) * disp, data = mtcars)

marg <- margins(model = mod, at = list(am = c("0", "1")))
marg
#> Average marginal effects at specified values
#> lm(formula = mpg ~ as.factor(am) * disp, data = mtcars)
#>  at(am)     disp    am1
#>       0 -0.02758 0.4518
#>       1 -0.05904 0.4518

emmeans(mod, c("am", "disp")) %>% 
  contrast(method = "pairwise")
#>  contrast                    estimate   SE df t.ratio p.value
#>  0,230.721875 - 1,230.721875   -0.452 1.39 28 -0.325  0.7479

emmeans(mod, c("am", "disp")) %>% 
  contrast()
#>  contrast            estimate    SE df t.ratio p.value
#>  0,230.721875 effect   -0.226 0.696 28 -0.325  0.7479 
#>  1,230.721875 effect    0.226 0.696 28  0.325  0.7479 
#> 
#> P value adjustment: fdr method for 2 tests

或者干脆使用summary()

library(margins)

data("mtcars")
mod <- lm(mpg ~ as.factor(am) * disp, data = mtcars)

marg <- margins(model = mod, at = list(am = c("0", "1")))
marg
#> Average marginal effects at specified values
#> lm(formula = mpg ~ as.factor(am) * disp, data = mtcars)
#>  at(am)     disp    am1
#>       0 -0.02758 0.4518
#>       1 -0.05904 0.4518

summary(marg)
#>  factor     am     AME     SE       z      p   lower   upper
#>     am1 1.0000  0.4518 1.3915  0.3247 0.7454 -2.2755  3.1791
#>     am1 2.0000  0.4518 1.3915  0.3247 0.7454 -2.2755  3.1791
#>    disp 1.0000 -0.0276 0.0062 -4.4354 0.0000 -0.0398 -0.0154
#>    disp 2.0000 -0.0590 0.0096 -6.1353 0.0000 -0.0779 -0.0402

reprex 包(v0.3.0)于 2019 年 6 月 7 日创建

于 2019-06-07T16:15:57.823 回答