stata连享会

Stata: Logit 模型简介

2019-01-14  本文已影响96人  stata连享会

作者:杨柳 (西北大学) || 连玉君 (中山大学)

Stata 连享会: 知乎 | 简书 | 码云

Stata连享会 精彩推文1 || 精彩推文2

1. 引言

在实证研究中,我们会经常遇到被解释变量为 “是/否” 或者 “某事件发生/未发生”,此时,被解释变量只有 两种 取值,对应的数字编码可记为 01,称为二值变量0-1 变量。例如,我们想研究以下问题:

在这些情况下,采用线性回归模型研究被解释变量的影响因素是 不合适的 (可参见 Stata 连享会 推文 二元选择模型:probit 还是 logit?),需要使用 概率模型,其中最常用的是 二元 Logit 模型。接下来,小编就带大家一起来学习 二元 Logit 模型

2. 二元 Logit 模型

2.1 (0-1) 分布

对于某一个样本 i,由于被解释变量 y_{i} 的取值为 0 或者 1,我们可以将 y_{i} 视为随机变量 Y_{i} 的实现值:Y_{i} 取 1 的概率为 \pi_{i}Y_{i} 取 0 的概率为 1- \pi_{i}随机变量 Y_{i} 服从参数为 \pi_{i} 的 (0-1) 分布Y_{i} 的分布律为
\Pr[Y_{i}=y_{i}]=\pi_{i}^{y_{i}}(1-\pi_{i})^{1-y_{i}}, \quad y_{i}=0,1 \qquad (1) 显然,若 y_{i}=1,则 Y_{i} 的概率为 \pi_{i};若 y_{i}=0,则Y_{i} 的概率为 1-\pi_{i}

易于证明,Y_{i}期望方差分别为:
E(Y_{i})=\pi_{i} Var(Y_{i})=\pi_{i}(1-\pi_{i}) 可见,Y_{i} 的期望和方差都决定于 \pi_{i}任何影响概率的因素不但会影响观察值的均值,也会影响其方差。这就表明线性回归模型无法用于分析二值变量,因为线性回归模型 假设方差是固定不变的

2.2 Logit 变换

为了使上述模型更富有弹性,我们假设概率 \pi_{i} 受一系列变量的影响,设定为 \mathbf{x}_{i}^{'} 。一个非常直觉的想法是把二者之间的关系设定为线性函数: \pi_{i}=\mathbf{x}_{i}^{'}\boldsymbol{\beta} \qquad (1) 其中, \boldsymbol{\beta} 为系数向量。该模型通常称为线性概率模型,采用普通最小二乘法估计即可。其 主要缺陷 在于:由于等式左边的 \pi_{i} 表示概率,所以必须介于 0 和 1 之间,而右边的线性组合项则可能取任何值,所以在不对模型做严格约束的情况下,我们很难保证模型的预测值介于合理的范围内。

因此,我们必须对概率 \pi_{i} 进行变换以消除对其取值范围的约束,继而把变换后的数值设定为解释变量 x_{i} 的线性函数。处理过程包括两个步骤。

  • 第一步,我们依据概率 \pi_{i} 来定义胜算比 (odds)
    \Omega_{i}=\frac{\pi_{i}}{1-\pi_{i}} \qquad (2)y_{i}=1 的概率 \pi_{i}y_{i}=0 的概率 1-\pi_{i} 的比值。显然,胜算比可以取任意非负值,如此便可消除上限约束。
  • 第二步,取对数以计算 logit 或 log-odds
    logit(\Omega_{i})=\ln(\Omega_{i})=\ln(\frac{\pi_{i}}{1-\pi_{i}}) \qquad (3) 这样我们就可以去除下限约束。因为,随着概率 \pi_{i} 趋近于 0,logit 将趋近于 -\infty ; 而当概率 \pi_{i} 趋近于1 ,logit 将趋近于 +\infty 。因此,通过以上变换,logit(\Omega_{i}) 将概率� \pi_{i} 的取值范围从 (0, 1) 映射至整个实数轴。显然,如果概率为 0.5,胜算比为1,相应的 log-odds 为 0。log-odds 为负表示概率小于0.5,反之则表示概率大于0.5。

2.3 Logistic 模型

在完成了上述变换后,我们就可以定义 Logistic 回归模型了,此时 我们假设概率 \pi_{i} 的 logit 变换(而非概率 \pi_{i} 本身)服从线性模型,即
logit(\Omega_{i})=\ln(\frac{\pi_{i}}{1-\pi_{i}})=\mathbf{x}_{i}^{'}\boldsymbol{\beta} \qquad (4) 其中,\mathbf{x}_{i} 为解释变量构成的向量,\boldsymbol{\beta} 为系数向量。

由于 logit 变换是一一对应的,所以我们可以通过求取逆对数由 logit 反向得到概率值 (通常称为 antilogit )。由上式可解得:
\pi(\mathbf{x}_{i})=\frac{\exp(\mathbf{x}_{i}^{'}\boldsymbol{\beta})}{1+\exp(\mathbf{x}_{i}^{'}\boldsymbol{\beta})} \qquad (5) 进一步将被解释变量 {y}_{i} 表示为:
{y}_{i}=\pi(\mathbf{x}_{i})+\varepsilon_{i} \qquad (6) 其中,\varepsilon_{i} 为随机干扰项,有 两个 可能的取值。若 y_{i}=1,则 \varepsilon_{i}=1-\pi(\mathbf{x}_{i}),相应的概率为 \pi(\mathbf{x}_{i}); 若 y_{i}=0,则 \varepsilon_{i}=-\pi(\mathbf{x}_{i}),相应的概率为 1-\pi(\mathbf{x}_{i}) 。因此,\varepsilon_{i} 服从均值为 0,方差为 \pi(\mathbf{x}_{i})[1-\pi(\mathbf{x}_{i})] 的分布

  • 综合上面的介绍,可以看出当被解释变量是二值变量时:
    (1) 模型的条件均值必须限定于 0 和 1 之间;
    (2) 干扰项服从均值为 0,方差为 \pi(\mathbf{x}_{i})[1-\pi(\mathbf{x}_{i})] 的分布,而非正态分布,且其分布受所分析样本的具体情况的影响;
    (3) 分析线性模型的基本准则同样适用于分析 Logit 模型。

2.4 估计结果的解释与举例

在 Logistic 模型 (4) 式中,系数 \beta _{j} 的含义可解释如下:在其他变量不变的情况下,第 j 个变量 x_{j} 变动一个单位,\ln(\Omega_{i}) 的值将变动 \beta _{j} 个单位。问题在于,\ln(\Omega_{i}) 变动 \beta _{j} 个单位的经济含义很难解释。不过,我们可以对 (4) 式两边
同时进行指数运算,得到第 i 个观察值对应的胜算比:
\Omega_{i}=\exp(\mathbf{x}_{i}^{'}\boldsymbol{\beta}) \qquad (7) 假设第 j 个解释变量 x_{j} 增加 1 个单位,则胜算比变为:
\Omega_{i}(\mathbf{x}_{i}, {x}_{ij}+1)=\exp(\mathbf{x}_{i}^{'}\boldsymbol{\beta})\exp(\beta_{j}) \qquad (8) 结合 (7) 式和 (8) 式,可以得到奇比 (odds ratio)[^Note1]:
\frac{\Omega(\mathbf{x}, {x}_{j}+1)}{\Omega(\mathbf{x}, {x}_{j})}=\exp(\beta_{j}) \qquad (9) [^Note1]: 为了表述的方便,(9) 式中省略了下标 i

因此,\exp(\beta_{j}) 的含义 可以解释为:在其他变量不变的情况下,若x_{j} 增加 1 个单位,则 胜算比 (odds) 是原来的 \exp(\beta_{j})。若\exp(\beta_{j})>1,则表示在其他变量不变的情况下, x_{j} 增加 1 个单位时, 胜算比增大了,此时对应的\beta_{j}>0;若 \exp(\beta_{j})<1,则表示在其他变量不变的情况下, x_{j} 增加 1 个单位时, 胜算比减小了,此时对应的 \beta_{j}>0

在某些情况下,我们还需要分析 x_{j} 变动 1 个标准差 产生的影响: 在控制其他变量不变的情况下,若 x_{j} 增加 1 个标准差,则 胜算比 (odds) 是原来的 \exp(\beta_{j}\times s_{j}) 倍。

Stata 提供的 logit 命令可用于估计上面介绍的 二元logit 模型。这里,我们使用 Stata 附带的 auto.dta 数据 (1978年美国汽车数据) 来预测汽车产地 (进口= 1;国产= 0)。以 foreign (是否进口车) 作为被解释变量、以 mpg (每加仑汽油能够行驶的英里数)weight (汽车重量) 作为解释变量建立二元 logit 模型。Stata 中的命令和结果如下所示:

.  sysuse "auto.dta", clear
(1978 Automobile Data)

.  logit foreign mpg weight

Iteration 0:   log likelihood =  -45.03321  
Iteration 1:   log likelihood = -29.238536  
Iteration 2:   log likelihood = -27.244139  
Iteration 3:   log likelihood = -27.175277  
Iteration 4:   log likelihood = -27.175156  
Iteration 5:   log likelihood = -27.175156  

Logistic regression                             Number of obs     =         74
                                                LR chi2(2)        =      35.72
                                                Prob > chi2       =     0.0000
Log likelihood = -27.175156                     Pseudo R2         =     0.3966

------------------------------------------------------------------------------
     foreign |      Coef.   Std. Err.      z    P>|z|     [95% Conf. Interval]
-------------+----------------------------------------------------------------
         mpg |  -.1685869   .0919175    -1.83   0.067    -.3487418     .011568
      weight |  -.0039067   .0010116    -3.86   0.000    -.0058894    -.001924
       _cons |   13.70837   4.518709     3.03   0.002     4.851859    22.56487
------------------------------------------------------------------------------
  • 回归结果显示 mpgweight 的系数值显著为负,表明当控制其他变量不变时,mpgweight 变量的值分别增加 1 个单位时, 胜算比 (odds) 将分别变为原来的 \exp(-0.1685)=0.8448\exp(-0.0039)=0.9961 倍,表明 该车辆是进口车的概率将变小
  • 要获得这个变化的倍数值,即 奇比,只需在 logit 命令后附加 or 选项即可。Stata 中的命令和结果如下所示:
.  logit foreign mpg weight, or

Iteration 0:   log likelihood =  -45.03321  
Iteration 1:   log likelihood = -29.238536  
Iteration 2:   log likelihood = -27.244139  
Iteration 3:   log likelihood = -27.175277  
Iteration 4:   log likelihood = -27.175156  
Iteration 5:   log likelihood = -27.175156  

Logistic regression                             Number of obs     =         74
                                                LR chi2(2)        =      35.72
                                                Prob > chi2       =     0.0000
Log likelihood = -27.175156                     Pseudo R2         =     0.3966

------------------------------------------------------------------------------
     foreign | Odds Ratio   Std. Err.      z    P>|z|     [95% Conf. Interval]
-------------+----------------------------------------------------------------
         mpg |   .8448578   .0776572    -1.83   0.067     .7055753    1.011635
      weight |   .9961009   .0010077    -3.86   0.000     .9941279    .9980779
       _cons |   898396.7    4059594     3.03   0.002     127.9781    6.31e+09
------------------------------------------------------------------------------
Note: _cons estimates baseline odds.
  • Stata 还提供了另一条命令: logistic,可以得到完全相同的结果,命令如下(回归结果与上面相同,故省略):
.  logistic foreign mpg weight, or
  • 回归结果显示,mpgweight 系数对应的 奇比均小于 1,表明当车辆的 mpgweight 增加时,胜算比减小了,说明 该车辆是进口车的概率将减小

2.5 边际效应分析

从上面的例子中可以看到,虽然我们能从解释变量的系数值判断该车辆是进口车的概率值的变化情况,但我们仍然不能直观地看出概率值变化的具体数值。在实证分析中,我们往往需要计算出来解释变量的变化将对这辆车是进口车的概率的变化带来的边际影响值,此时,我们需要使用 margins 命令来计算边际效应。

  • 接下来,小编向大家介绍实证分析中常用的两种边际效应:
    (1) 平均边际效应。 即先分别计算在每个样本观测值上的边际效应,然后进行简单算术平均得到平均边际效应。
    (2) 样本均值处的边际效应。 即先分别计算各自变量的样本均值,然后计算在这一点处的边际效应。
  • (1) 平均边际效应

我们可以使用 margins 命令附加 dydx 选项来进行计算。Stata 中的命令和结果如下所示:

.  margins, dydx(mpg) //计算mpg变量对Pr(foreign)的平均边际效应

Average marginal effects                        Number of obs     =         74
Model VCE    : OIM

Expression   : Pr(foreign), predict()
dy/dx w.r.t. : mpg

------------------------------------------------------------------------------
             |            Delta-method
             |      dy/dx   Std. Err.      z    P>|z|     [95% Conf. Interval]
-------------+----------------------------------------------------------------
         mpg |  -.0197187   .0096987    -2.03   0.042    -.0387277   -.0007096
------------------------------------------------------------------------------

.  margins, dydx(weight) //计算weight变量对Pr(foreign)的平均边际效应

Average marginal effects                        Number of obs     =         74
Model VCE    : OIM

Expression   : Pr(foreign), predict()
dy/dx w.r.t. : weight

------------------------------------------------------------------------------
             |            Delta-method
             |      dy/dx   Std. Err.      z    P>|z|     [95% Conf. Interval]
-------------+----------------------------------------------------------------
      weight |  -.0004569   .0000571    -8.01   0.000    -.0005688   -.0003451
------------------------------------------------------------------------------
  • (2) 样本均值处的边际效应

我们可以使用 margins 命令附加 dydxatmeans 选项来进行计算。Stata 中的命令和结果如下所示:

.   margins, dydx(mpg) atmeans //计算mpg变量对Pr(foreign)在样本均值处的边际效应

Conditional marginal effects                    Number of obs     =         74
Model VCE    : OIM

Expression   : Pr(foreign), predict()
dy/dx w.r.t. : mpg
at           : mpg             =     21.2973 (mean)
               weight          =    3019.459 (mean)

------------------------------------------------------------------------------
             |            Delta-method
             |      dy/dx   Std. Err.      z    P>|z|     [95% Conf. Interval]
-------------+----------------------------------------------------------------
         mpg |  -.0223512   .0127037    -1.76   0.079      -.04725    .0025476
------------------------------------------------------------------------------

.   margins, dydx(weight) atmeans  //计算weight变量对Pr(foreign)在样本均值处的边际效应

Conditional marginal effects                    Number of obs     =         74
Model VCE    : OIM

Expression   : Pr(foreign), predict()
dy/dx w.r.t. : weight
at           : mpg             =     21.2973 (mean)
               weight          =    3019.459 (mean)

------------------------------------------------------------------------------
             |            Delta-method
             |      dy/dx   Std. Err.      z    P>|z|     [95% Conf. Interval]
-------------+----------------------------------------------------------------
      weight |  -.0005179   .0001389    -3.73   0.000    -.0007902   -.0002457
------------------------------------------------------------------------------

3. 二元 Logit 模型在实证分析中的应用举例

原文信息: Yizhen Gu, Elizabeth Deakin, Ying Long, The Effects of Driving Restrictions on Travel Behavior Evidence from Beijing, Journal of Urban Economics, 102 (2017), 106-122. DOI: /10.1016/j.jue.2017.03.001

3.1 研究背景与问题

汽车限号行驶政策已在国际上的很多城市实施过,许多研究质疑了这类政策的有效性,还有一些研究结果之间的结论存在差异。该篇文章将北京市作为案例城市,使用 微观的个人出行数据 来研究 汽车限号行驶政策 对居民的 出行行为 的影响,这些 出行行为 包括:

3.2 北京市汽车限号行驶政策

  • 周一限行 1 和 6
  • 周二限行 2 和 7
  • 周三限行 3 和 8
  • 周四限行 4 和 9
  • 周五限行 0 和 5

3.3 研究数据

该篇文献的研究数据来源于北京市政府组织调查的 2010 年的居民出行调查数据。调查采取 1% 的样本率,选择了北京市 1911 个交通小区的 1085 个区域,每个区域有 10 - 50 的家庭接受面对面的调查。居民出行调查的样本总量共计 46,900 家庭及 116,142 位家庭成员。居民出行调查内容包括:家庭成员一天的出行情况,家庭的信息和家庭成员的个人信息。在该样本中,有 71.1% 的家庭没有私人汽车, 26.7% 的家庭有一辆, 2.2% 的家庭至少有两辆车。

3.4 二元 Logit 模型建立与结果分析

  • (1) 基本模型

文献中建立的 二元 Logit 模型 如下所示:
log(\frac{P_{i,t}}{1-P_{i,t}})=\gamma _{1}Restrict_{i}+\gamma _{2}AdjacentRestrict_{i}+\gamma _{3}Restrict_{i}*z_{t} \\ +\gamma _{4}Restrict_{i}*x_{i}+\gamma _{5}x_{i}+\gamma _{6}z_{t} \qquad (10) 式中,P_{i,t} 为驾驶员 i 的第 t次出行使用小汽车方式的概率; Restrict_{i} 为虚拟变量,代表调查当日驾驶员的小汽车是否被限行;AdjacentRestrict_{i} 为虚拟变量,代表调查日是否与驾驶员的小汽车被限行的日子相邻; z_{t} 为代表一系列反映出行特征的协变量,例如出行目的、出行距离/时耗等;x_{i} 代表一系列反映驾驶员特征的协变量,例如性别、年龄、职业、收入、家庭距离地铁站的距离、家庭是否位于市中心的东城区或西城区、家庭是否位于长安街以北的地区等。

  • (2) 是否在工作日至少使用小汽车方式进行一次出行

基于二元 Logit 模型公式(10) ,将因变量中的 P_{i} 定义为 驾驶员 i 在工作日至少使用小汽车方式出行一次的概率,选择 两类驾驶员 样本进行回归分析,分别为以下两类,回归结果如下图所示:

图1 限行政策对驾驶员在工作日至少使用小汽车方式出行一次的边际效果.png

图中报告的结果是 平均边际效应值;括号内数字为 聚类标准误,以 驾驶员的家庭所在的交通小区 作为聚类的单元。从回归结果可以看到,限行政策使受限区域的驾驶员在工作日至少使用小汽车方式出行一次的概率显著降低了 8.1%-15.9%

为了进一步说明限行政策的边际效果在不同组别的驾驶员之间的区别,文献中按照性别机动/固定的工作时间高/低收入居住在长安街以北/以南地区是否有孩子对位于 受限区域且家庭中拥有一辆小汽车 的 5123 名驾驶员样本进行了 分组回归,结果如下图所示:

图2 限行政策分组回归的边际效果.png

图中报告的结果是各组的 平均边际效应值至少使用小汽车方式出行一次的样本均值调查当日被限行时至少使用小汽车方式出行一次的平均概率的下降情况;括号内数字为 聚类标准误,以 驾驶员的家庭所在的交通小区 作为聚类的单元。从回归结果可以看到,女性有机动工作时间高收入的驾驶员在 调查当日被限行时 至少使用小汽车方式出行一次的平均概率比 未限行时 的平均概率的 减小幅度 较大,分别为 27.3%、24% 与 18%

  • (3) 对限行政策的适应 / 调整行为

文献基于二元 Logit 模型公式 (10) ,将因变量中的 P_{i,t} 定义为 小汽车出行的概率,使用了 四组样本 进行回归分析,结果如下图所示:

图3 限行政策对选择小汽车出行方式的边际效果.png

从上图 Panel 1 中 I 与 II 列的回归结果可以看到,在 居住在受限区域且家庭拥有一辆车 的驾驶员进行的所有出行中,限行政策使得小汽车方式出行的平均概率显著的下降了近 10%;III 与 IV 列的回归结果显示,在 居住在非受限区域且家庭拥有一辆车的驾驶员进行的所有出行中,限行政策未能显著的影响小汽车方式出行的平均概率。

(4) 稳健性检验

文献中对 出行方式选择 的回归结果还进行了一系列的 稳健性检验,包括以下方法:

参考文献

关于我们

联系我们

往期精彩推文


欢迎加入Stata连享会(公众号: StataChina)
上一篇 下一篇

猜你喜欢

热点阅读