如何从xtlogit模型重现平均边际效应

问题描述

我对从随机效应logit模型（使用xtlogit在Stata中运行）中再现平均边际效应感兴趣。我了解如何使用Delta方法从logit模型重现平均边际效应。例如，在下面的代码中，我成功地再现了age中报告的margins的平均边际效应。

*** Stata code
* download data
webuse union,clear

* calculate delta and copy variable of interest - age
sum age
gen xdelta = r(sd)/1000
clonevar age_ = age

* run logit model and calculate average marginal effect using margins
logit union age_
margins,dydx(age_)

* calculate average marginal effect by hand - mean of xme equals result from margins above
predict p1
replace age_ = age_+xdelta
predict p2
gen xme = (p2 - p1) / xdelta
sum xme

* calculate average marginal effect at fixed value of age using margins
margins,at(age=(16))

* calculate average marginal effect by hand - mean of p3 equals result from margins above
replace age_ = 16
predict p3
sum p3

我正在苦苦挣扎的是再现xtlogit模型的平均边际效应。

*** Stata code
* download data and designate panel variables
webuse union,clear
xtset idcode year

* run xtlogit model
xtlogit union age

* calculate average marginal effects - can't figure out how to reproduce these estimates :(
margins,dydx(*)
margins,at(age=(16))

对于弄清楚如何重现xtlogit边际效应的任何帮助，我们深表感谢。谢谢！

-----已编辑，以更清楚地表明我对复制margins报告的估算值感兴趣

解决方法

有两种方法可以做到这一点，但本质上问题可以归结为事实

$$ \ Pr（y_ {it} = 1 \ vert x_ {it}）= \ int \ Lambda（u_i + x_ {it}'\ beta）\ cdot \ varphi（0，\ sigma_u ^ 2）du_i $$

其中$ \ varphi（）$是正常密度。在您的代码中，您实际上将随机效果$ u_i $设置为零（predict(pu0)就是这样做的）。这会将RE设置为平均值，这可能与您的想法不符。当然，{u} {i}不会被xtlogit,re观察甚至估计，因此，如果要复制predict(pr)的作用，则需要将随机效应进行积分，以使用估计方差得出无条件概率

在Stata中执行此操作的一种方法是使用用户编写的integrate命令执行一维数值积分，如下所示：

webuse union,clear
xtset idcode year
xtlogit union age,nolog
margins,at(age=(16)) predict(pr)
margins,dydx(*) at(age=16) predict(pr)
capture ssc install integrate
/* phat at age 16 */
integrate,f(invlogit(x - 3.0079682 + .01929225*16)*normalden(x,2.477537654945496)) l(-10) u(10) vectorise quadpts(1000)
/* ME at age 16 */
integrate,f(invlogit(x - 3.0079682 + .01929225*16)*(1-invlogit(x - 3.0079682 + .01929225*16))*(.01929225)*normalden(x,2.477537654945496)) l(-10) u(10) vectorise quadpts(1000)

如果使用实际系数（例如_b [_cons]和e（sigma_u））而不是将值粘贴到其中，则可能会获得更高的精度。

也许可以使用Mata，Python甚至是模拟方法来执行此操作，但效率更高，但我将由您自己决定。您还可以使用xtmelogit获得RE的后验模态估计。