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)
##   rowid     type           term contrast       dydx   std.error statistic
## 1     1 response bill_length_mm    dY/dX 0.01643670 0.007567926  2.171890
## 2     2 response bill_length_mm    dY/dX 0.03363322 0.011641242  2.889144
## 3     3 response bill_length_mm    dY/dX 0.08056324 0.021049901  3.827250
## 4     4 response bill_length_mm    dY/dX 0.03388757 0.006397870  5.296696
## 5     5 response bill_length_mm    dY/dX 0.04820226 0.013242635  3.639930
## 6     6 response bill_length_mm    dY/dX 0.01543305 0.007016936  2.199400
##        p.value    conf.low  conf.high fat_penguin bill_length_mm
## 1 2.986394e-02 0.001603842 0.03126957           0           39.1
## 2 3.862921e-03 0.010816810 0.05644964           0           39.5
## 3 1.295827e-04 0.039306194 0.12182029           0           40.3
## 4 1.179166e-07 0.021347979 0.04642717           0           36.7
## 5 2.727123e-04 0.022247175 0.07415735           0           39.3
## 6 2.784952e-02 0.001680104 0.02918599           0           38.9
##   flipper_length_mm species
## 1               181  Adelie
## 2               186  Adelie
## 3               195  Adelie
## 4               193  Adelie
## 5               190  Adelie
## 6               181  Adelie

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

mfx %>% 
  count(term, contrast)
##                term           contrast   n
## 1    bill_length_mm              dY/dX 333
## 2 flipper_length_mm              dY/dX 333
## 3           species Chinstrap - Adelie 333
## 4           species    Gentoo - Adelie 333

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

60.3 平均边际效应 

summary(mfx)
## Average marginal effects 
##                Term           Contrast   Effect Std. Error z value   Pr(>|z|)
## 1    bill_length_mm              dY/dX  0.02684   0.005921  4.5328 5.8212e-06
## 2 flipper_length_mm              dY/dX  0.01084   0.002416  4.4871 7.2211e-06
## 3           species Chinstrap - Adelie -0.40759   0.058147 -7.0097 2.3887e-12
## 4           species    Gentoo - Adelie  0.06581   0.109543  0.6008    0.54798
##       2.5 %   97.5 %
## 1  0.015234  0.03844
## 2  0.006106  0.01558
## 3 -0.521561 -0.29363
## 4 -0.148889  0.28051
## 
## Model type:  glm 
## Prediction type:  response