Skip to main content

数据科学中的 R 语言

第 60 章 模型的边际效应

本章介绍模型的边际效应，主要围绕marginaleffects宏包，本章的内容也是来源该宏包的说明文档。

60.1 边际效应

边际效应，测量的是某一个预测因子变化一个单位与伴随的响应变量的变化之间的关联。用数学语言表示，就是回归方程对x求偏导。

假定我们建立的回归方程是一个二次函数， \[ y = -x^2 \]

那么对x的偏导数 \[ \frac{\partial y}{\partial x} = -2x \]

可以看到，此时的边际效应就是曲线的斜率

  1. 当\(x<0\)，斜率为正，x增加y也增加，即，边际效应为正
  2. 当\(x=0\)，斜率为0，在这个位置上边际效应为0
  3. 当\(x>0\)，斜率为负，x增加y也减少，在这个位置上边际效应为负

60.2 marginaleffects function

最简单的线性模型，每个因子的边际效应就是预测因子的系数，与因子的取值无关。但是复杂点模型，因子边际效应不仅仅与因子的取值有关，而且还与其它因子的值也有关。

我们下面用企鹅数据来说明。

我们先构建一个二元变量fat_penguin（是否为胖子）, 1表示是，0表示不是。并建立logitisc回归模型

dat <- penguins %>%
  drop_na() %>% 
  mutate(
    fat_penguin = if_else(body_mass_g > median(body_mass_g), 1, 0)
  )

mod <- glm(
  fat_penguin ~ bill_length_mm + flipper_length_mm + species,
  data = dat, 
  family = binomial(link = "logit")
)
mod
## 
## Call:  glm(formula = fat_penguin ~ bill_length_mm + flipper_length_mm + 
##     species, family = binomial(link = "logit"), data = dat)
## 
## Coefficients:
##       (Intercept)     bill_length_mm  flipper_length_mm   speciesChinstrap  
##          -42.0591             0.3485             0.1408            -5.0386  
##     speciesGentoo  
##            0.8165  
## 
## Degrees of Freedom: 332 Total (i.e. Null);  328 Residual
## Null Deviance:       461.3 
## Residual Deviance: 167.3     AIC: 177.3
mfx <- marginaleffects(mod, type = "response") 
head(mfx)
## 
##            Term Contrast Estimate Std. Error    z Pr(>|z|)    S    2.5 % 97.5 %
##  bill_length_mm    dY/dX   0.0164    0.00800 2.05  0.04000  4.6 0.000751 0.0321
##  bill_length_mm    dY/dX   0.0336    0.01201 2.80  0.00511  7.6 0.010091 0.0572
##  bill_length_mm    dY/dX   0.0806    0.02099 3.84  < 0.001 13.0 0.039427 0.1217
##  bill_length_mm    dY/dX   0.0339    0.00634 5.35  < 0.001 23.4 0.021466 0.0463
##  bill_length_mm    dY/dX   0.0482    0.01334 3.61  < 0.001 11.7 0.022064 0.0743
##  bill_length_mm    dY/dX   0.0154    0.00744 2.07  0.03813  4.7 0.000845 0.0300
## 
## Columns: rowid, term, contrast, estimate, std.error, statistic, p.value, s.value, conf.low, conf.high, predicted, predicted_hi, predicted_lo, fat_penguin, bill_length_mm, flipper_length_mm, species

marginaleffects() 函数对数据框dat的每一行观测给出了边际效应估计，最后输出一个数据框 。注意到，我们的模型有3个预测因子（两个连续变量，一个离散变量），一个预测因子对应一个与原数据框等长的数据框，因此最终返回的结果是原来数据框长度的3倍。

mfx %>% 
  count(term, contrast)
## 
##               Term           Contrast
##  bill_length_mm    dY/dX             
##  flipper_length_mm dY/dX             
##  species           Chinstrap - Adelie
##  species           Gentoo - Adelie   
## 
## Columns: term, contrast, n

边际效应对连续变量非常适合。但离散变量的边际效应，不太好理解，因此采用对照方法，具体为，以某一层级为基线，那么从基线切换到其它层级，伴随着响应变量的变化，就为离散变量的边际效应。

60.3 平均边际效应 

summary(mfx)
## 
##               Term                       Contrast Estimate Std. Error      z
##  bill_length_mm    mean(dY/dX)                      0.0268    0.00591  4.538
##  flipper_length_mm mean(dY/dX)                      0.0108    0.00241  4.506
##  species           mean(Chinstrap) - mean(Adelie)  -0.4076    0.05786 -7.044
##  species           mean(Gentoo) - mean(Adelie)      0.0658    0.10956  0.601
##  Pr(>|z|)    2.5 %  97.5 %
##    <0.001  0.01525  0.0384
##    <0.001  0.00613  0.0156
##    <0.001 -0.52100 -0.2942
##     0.548 -0.14892  0.2805
## 
## Columns: term, contrast, estimate, std.error, statistic, p.value, conf.low, conf.high