#STATA

STATA本身有很多estimator是通过MLE方法估计的,例如logit, probit等。在这些模型之外,STATA同时提供了ML syntax来拓展可以估计maximum likelihood模型。下面我们举例子说明ML syntax的用法。

LF Estimator

当整体样本的Log likelihood可以通过每个样本点的log likelihood累加得到时,STATA认为这类模型符合线性约束(Linear Form Restriction),可以使用STATA ML syntax中的lf方法来估计这类模型。我们将使用lf方法来估计四种常见的模型: binary logit, binary probit, OLS, and mixed logit model.

Logit Model

Log Likelihood formula for Logit Model:
prob(y=1) = exp(xb)/(1+exp(xb))

Code:
‘’’

program logit_lf
    args lnfj theta 
    tempvar p 
    local y $ML_y1
    quietly gen `p' = exp(`theta')/(1+exp(`theta'))
    quietly replace `lnfj' = ln(`p') if `y' == 1
    quietly replace `lnfj' = ln(1-`p') if `y' == 0
end

sysuse auto, clear
ml model lf logit_lf (eq1: foreign = weight length)
ml maximize

‘’’

几个要注意的点:

+ args 是STATA提供的一种parse的工具,是把positional argument变成指定的local。例如,ml提供的第一个argument通过args变成了local lnfj,而在剩下的程序里,都将以`lnfj'的形式出现。
+ program 中的theta是所有x和系数的linear combination。这是使用lf方法最方便的地方,不需要考虑单独的系数大小,而是把系数的linear combination放在一起考虑。
+ lnfj 是每个observation 的log likelihood.
+ $ML_y1 是第一个公式中的 dependent variable,以global 的形式存在。
+ 所有program 中使用的新的变量要用tempvar
+ ml model lf logit_lf (eq1: foreign = weight length) 这一行是定义模型,即使用lf方法的logit_lf模型,模型中的y是foreign,而解释变量包括了weight和length。
+ ml maximize 这一行是对上述模型进行求解,计算出对应最大loglikelihood的参数。

Probit Model

类似地,我们可以通过lf方法估计probit模型:
‘’’

program probit_lf
    args theta lnfj
    local y $ML_y1
    quietly replace `lnfj' = ln(normal(`theta')) if `y' == 1
    quietly replace `lnfj' = ln(-normal(`theta')) if `y' == 0
end

‘’’

Linear Regression Model

线性回归模型(homoskedastic standard errors)也可以通过lf方法实现:
‘’’
program ols_lf
args lnfj theta std
local y $ML_y1
quietly replace lnfj' = ln(normalden(y’,theta',std’))
end

ml model lf ols_lf (eq1: mileage = weight length) (eq2:)
* ml model lf ols_lf (eq1: mileage = weight length) /eq2 
ml maximize 

‘’’

要注意的点:

  • Linear regression model和前面的logit probit模型不同的地方在于:线性回归模型的log likelihood不能通过一个theta(linear combination of variables)来表达,而是由两个系数同时决定的,一个是所有变量的线性组合,一个是残差项的standard deviation。因此,我们需要在写函数,和使用ml model时候都要做出相应的调整。首先,ml model 要加入第二个公式,即(eq2:)。这里,由于模型的假设包括了homoskedasticity,残差项的std不和任何解释变量相关,因此(eq2:)中不需要包括任何解释变量,只需要一个constant变量即可。
  • 当放松homoskedasticity的假设时,我们只需要稍微修改ml model中的eq2部分即可。

Mixed Binary Logit Model

Logit model虽然直观易估计,但是对于individual preference有比较强的限制。因此,我们可以估计一个mixed binary logit model,来解决这个问题。但是mixed logit model的log likelihood表达式没有closed form solution,因此只能通过simulation来解决,而STATA的ML syntax可以较好的提供simulation的方法。

// Simulate data for estimation 
/*
x~N(0,1)
beta1 = 1
beta2 ~ N(1,1)
epsilon~ extreme value distribution
y = (beta_1+ beta_2*x +epsilon>0)
*/

clear
set obs 1000
set seed 10101
gen x = rnormal()
gen beta_1 = 1
gen beta_2 = rnormal(1,1)
gen u = uniform()
gen epsilon = log(u)-log(1-u)
gen y = (beta_1 + beta_2*x+epsilon) > 0
forvalues i = 1/1000{
    quietly gen draw_`i' = uniform()
}

program mixed_logit
   args lnfj beta_1 beta_2 std 
   tempvar p sim_f sim_avg_f
   quietly gen `sim_avg_f' = 0
   quietly{
   forvalues i = 1/1000{
       gen `p' = exp(`beta'1+`beta_2'*x +  /// `std'*`draw_i'*x)/(1+exp(`beta'1+`beta_2'*x + `std'*`draw_i'*x))
       gen `sim_f' = `p' if $ML_y1 == 1
       replace `sim_f' = 1- `p' if $ML_y1 == 0
       replace `sim_avg_f` = `sim_avg_f` + `sim_f`/1000
   }
   replace `lnfj' = ln(sim_avg_f)
   }
end

ml model lf mixed_logit (beta_1:y = )(beta_2:x)(std:)
ml maximize

这里需要注意的是:

  • Log likelihood是一千次Simulation的平均。
  • 由于Log likelihood不能够通过所有变量的线性组合来表示,我们需要把系数拆分成三个部分,常数项的系数,变量X的系数的均值,变量X的std.这样的拆分不影响我们使用lf方法,因为样本的log likelihood仍然是由所有的样本的Log likelihood加总得到。

D0 Estimator

上述lf estimator只适用于样本的log likelihood仍然是由所有的样本的Log likelihood加总得到的情况,当上述条件不成立时,我们需要使用STATA提供的d0 estimator. D0 estimator的特点是我们需要在模型中提供整体样本的log likelihood,而不是每个单独样本点的log likelihood。下面我们用Conditional logit model来说明d0 estimator的用法。

Conditional Multinomial Logit Model

Model:
Individual: indexed by i
Choice situation: indexed by j
U_ij = XB+ epsilon, epsilon~extreme value distribution
Log Likelihood if choice j is chosen: exp(x_kb)/\sum(exp(x_jb))

Code:

cls
program drop _all
webuse choice,clear
gen japan = car==2
gen europe = car ==3 

program clogit_sch 
version 14
/* 
model: lnfi = exp(xb)/(\sum exp(xb))
*/
args todo b lnf 
tempvar xb e_xb sum_e_xb lnfj
// local group_id $ML_id 
local y $ML_y1
mleval `xb' = `b', eq(1)
sort id
quietly{
    gen `e_xb' = exp(`xb')
    bysort id: egen `sum_e_xb' = total(`e_xb')
    gen double `lnfj'   = ln(`e_xb'/`sum_e_xb')
    mlsum `lnf' = `lnfj' if `y' == 1
}
end


global ML_id id 
ml model d0 clogit_d0 (eq1: choice = dealer japan europe) 
ml check 
ml max

需要注意的是:

  • mleval: 计算系数b'带来的linear combination的xb’。
  • mlsum: 计算所有相关的loglikelihood的总和。
  • mlcheck:用来对ml model进行debug。

和STATA提供改的asclogit的结果进行对比

webuse choice,clear   
set more off
gen japan = car==2
gen europe = car ==3 
cd C:/Users/yan/Dropbox/college_entrance_exam/program/2018.12
program drop _all
global ML_id id 
ml model d0 clogit_sch (eq1: choice = dealer japan europe,noconstant) 
ml check 
ml max  
asclogit choice dealer, case(id) alternatives(car)

这里需要指出,asclogit会在后台生成choice specific fixed effect,而我们的estimator需要自己生成这些变量,加入到模型中,但是好处是我们的模型更加flexible。

另外几个很有用的ml command是 ml query, ml report, ml init,和 ml graph。ml query是显示目前处理的模型,包括了每条equation,变量,和系数的初始值。ml report会显示目前的系数vector, gradient vector, negative Hessian, and maximization direction。ml init用来设置模型参数的初始值,例如 ml init 1 2 -2, copy。Ml graph可以用来检测潜在的convergence issue。Newton-Rhapson 是默认的maximization routine,可以通过technique()选项来尝试其他方法,例如bhhh, dfp, and bfgs等。

[TOC]

版本

2018.06.24 初稿

2018.07.19 第二稿

2018.08.03 终稿

前言

描述性表格(Summary Statistics Table),作为每篇论文中最先亮相的表格,它的设计和内容在很大程度上决定了读者对文章的第一印象。Latex虽然可以制作出专业的表格,但是tex环境中输入各类型的描述性数字枯燥且容易出错。那么如何快速高效的制作出专业的描述性表格呢?

Stata中可以用来做描述性表格的命令众多,例如,tabulate, table, tabstat等。而本文将介绍的主要工具为Tabout命令。Tabout和上述提到的命令相比,能够制作更多类型的表格,且能够快捷的输出表格的tex code,从而方便的将表格插入到完整的tex文档中。

本文将介绍如何通过使用Tabout命令,设计出一个简洁而又内容丰富的描述性表格。

Tabout

Tabout的核心思想是将描述性表格分成了两类:frequency table和summary table. Frequency table, 顾名思义就是变量取值的频率分布表,而summary table则是输出例如均值,中位数等样本统计数字的表格。两种表格都可以分为oneway, 和twoway两种类型。在Tabout的逻辑下,可以通过一套语法,细调部分选项,来生成这几种表格。下面我们通过举例来说明tabout的完整逻辑,以及语法细节。

Frequency Table

Twoway Frequency Table

我们使用STATA自带的cancer.dta数据来制作frequency table。数据中的两个变量,drug, died,分别记录了病人使用的药的种类,以及是否死亡。下图统计了使用不同药物种类的人数,以及使用各个药物后的死亡人数的分布。我们将分别对输出这个表格的Stata Code进行详细解释。

Imgur

1
2
3
4
5
6
7
#delimit ;  //change delimiter to ;
tabout died drug using freq_two.tex, replace style(tex)
c(freq col cum) // column content
clab(Freq Col_Pct Cum_Pct) f(0 0 0) // column head and content formating
title(Table: Twoway Frequency Table) fn(Source:auto.dta) // title and footnote
topf(table_top.tex) botf(table_end.tex) ; // tex codes
#delimit cr
  • 第一行 tabout died drug using freq_one.tex, replace style(tex) 是基本设置,即选择所要输出的tex file,以及设定输出的格式为tex。 值得讨论的是 died drug这两个变量的顺序。 Tabout将变量分成列变量,和行变量,并且规定最后一个变量为行变量。行变量的意思很直白,就是每个变量的值在表格中以行的形式出现。 这个设定的好处是我们可以设置多个列变量,同时和行变量进行cross tabulation。

  • 第二行c()是frequency table的核心选项,即表格中要输出的内容:freq, col, cum分别代表了频数,列比,和累计百分比。

  • 第三行是输出类容的格式设定。clab的含义是输出内容的标题,即频数列应该被标记为Freq, 而列比应该被标记为Col_Pctf选项设置了每一类型的数值应当保留的小数点位数。 值得指出的是,tabout在行标题上的逻辑也较为独特,需要理解后才能很好掌握。 Tabout设置了3级行标题,第一级标题是行变量的名称,第二级行标题是行变量的具体取值,第三季标题是tabulate的内容,例如在这里是频数,列比例,和累计比例。这三级标题的设置可以通过h1, h2和h3来设置。 clab同样可以设置第三级标题,区别是clab可以改变所有panel里面的第三级标题,而h3只能改变第一个panel。

  • 第四行为表格标题,以及脚注的设定。

  • 第五行是整个使用tabout的workflow的核心。tabout输出的虽然是tex code,但是并不能直接用来编译。不能编译的结果就是我们必须手动的添加topf 以及botftabout的输出结果中。而当我们需要经常的调整表格内容时,每一次都需要手动添加这些内容显然非常的繁琐。而topf 以及botf使这一手动过程可以自动化。 tofbotf的内容是非常简单的tex代码,只是用来可以直接编译生成的latex表格代码,具体内容可见下面的code block。

tof文件的内容为:

1
2
3
4
\documentclass[leqno,11pt]{article}
\usepackage{booktabs}
\usepackage{tabularx}
\begin{document}

botf文件的内容是:

1
\end{document}

执行上述Stata命令之后,我们会得到一个tex文档。编译tex文档即可得到如下图所示的表格。

Oneway Frequency Table

虽然入门有难度,但是Tabout的好处是在理解了它的逻辑之后,稍微改动一些选项就可以生成其他类型的表格。例如,当我们需要生成oneway frequency table时,只需要在选项中加入oneway即可。

1
2
3
4
5
6
7
#delimit ;  //change delimiter to ;
tabout died using freq_one.tex, replace style(tex) oneway
c(freq col cum) // column content
clab(Freq Col_Pct Cum_Pct) f(0 2 2) // column head and content formating
title(Table:Oneway Frequency Table) fn(Source:auto.dta) // title,footnote
topf(table_top.tex) botf(table_end.tex) ; // tex codes
#delimit cr

Imgur

Summary Table

类似的,我们只需要微调选项既可以使用tabout生成summary table。 下面我们举例说明。

Twoway Summary Table

1
2
3
4
5
6
sysuse auto, clear
tabout rep78 foreign using table10.tex, replace ///
style(tex) font(italic) c(mean weight) f(0c) sum ///
twidth(9) h1(Car type (mean weight in lbs.)) h3(nil) ///
title(Table 10: Simple twoway summary table of means) ///
fn(auto.dta)

第三行是微调的核心,加入了sum选项,即表明这里需要生成的是summary table, 而不是frequency table。另外c()选项中的mean weight说明了我们需要的数字是车辆的平均重量。上述代码生成的表格如下图。

Imgur

Oneway Summary Table

类似的,我们也可以生成oneway summary table。具体可见下面的代码和图,不再赘述。

1
2
3
4
5
6
7
8
sysuse auto, clear
tabout foreign rep78 using table12.tex, replace ///
style(tex) font(bold) twidth(13) sum npos(tufte) ///
c(mean mpg mean weight mean length median price ///
median headroom) f(1c 1c 1c 2cm 1c) h2(Mean Median) ///
h2c(3 2) clab(MPG Weight_(lbs) Length_(in) Price ///
Headroom_(in)) title(Table 12: Oneway summary table ///
\\ with multiple summary measures) fn(Source: auto.dta)

Imgur

[TOC]

版本

2018.06.13 初稿
2018.07.15 第二稿

前言

STATA是回归分析的最常用的工具。当我们进行了大量的回归分析之后通常需要解决三个问题。怎样才能直观的展示需要关注的系数?如何才能方便的对比不同回归中的系数?怎样才能生成论文中可以直接使用的高质量的回归表格?本教程将试图对这三个问题给出自己的理解。

本教程将使用STATA中的Estout Package来回答以上三个问题。Estout功能非常强大,其中包括了estout, esttab, estadd等命令。Estout的具体细节介绍可以参考官方介绍文档。 但是由于它的选项非常丰富,学习起来有一定的难度,因此本文旨在对Estout的使用方法进行总结,来达到降低使用estout制作表格的难度的目的。

Estout介绍

我们将首先介绍如何使用Estout输出最基本的回归表格,接着介绍如何把描述性数字以及回归中控制的变量类型加入到表格中,最后是如何生成pdf格式的表格。

基本回归表格

使用Estout制作基本的表格非常简单,只需要在regress命令后使用esttab即可。 下面我们以STATA自带的auto数据来展示esttab的效果。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
. sysuse auto
(1978 Automobile Data)

. eststo: quietly regress price weight mpg
(est1 stored)

. eststo: quietly regress price weight mpg foreign
(est2 stored)

. esttab

--------------------------------------------
(1) (2)
price price
--------------------------------------------
weight 1.747** 3.465***
(2.72) (5.49)

mpg -49.51 21.85
(-0.57) (0.29)

foreign 3673.1***
(5.37)

_cons 1946.1 -5853.7
(0.54) (-1.73)
--------------------------------------------
N 74 74
--------------------------------------------
t statistics in parentheses
* p<0.05, ** p<0.01, *** p<0.001

. eststo clear

我们用了一个命令就得到了一个简单的回归表格,而这个回归表格已经能够满足清晰的对比不同回归中的相同变量的系数的目的。Esttab之所以如此简单,强大,因为它其实是更为复杂的estout命令的一个wrapper,也就是说当我们使用esttab命令时,其实已经使用了STATA设计好的初始设置。例如,回归系数,标准误差的小数点默认为三位,以及表格中会加入观测值数量等。当然,这个表格还不能满足我们的所有需求。当我们希望更为丰富的表格内容时该如何做呢?

加入描述性数据,以及控制变量类型

可以加入到回归表格中的描述性数据分为两类:回归本身返回的数据,例如R-Square, 观测值数量等;需要另外计算加入到回归表格中的数据,例如因变量的均值,控制变量类型等。

加入回归返回的数据

将回归本身返回的数据加入到表格的方法相对简单,只需要使用esttab中的选项stats即可。例如,当我们需要加入R Square和回归的观测值时,可以加入stats(r2 N, labels(“R Square” “Num of Obs”))。其中,r2,N分别代表了要加入的数据,而labels中的内容为数据的标签。具体效果可见下图。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
    
.esttab, stats(r2 N, labels("R Square" "Num of Obs"))

--------------------------------------------
(1) (2)
price price
--------------------------------------------
weight 1.747** 3.465***
(2.72) (5.49)

mpg -49.51 21.85
(-0.57) (0.29)

foreign 3673.1***
(5.37)

_cons 1946.1 -5853.7
(0.54) (-1.73)
--------------------------------------------
R Square 0.293 0.500
Num of Obs 74 74
--------------------------------------------
t statistics in parentheses
* p<0.05, ** p<0.01, *** p<0.001

加入其他信息

如果所需要加入回归表格中的信息不是回归中返回值时,我们需要用到estout package中的estadd命令。例如,当需要加入因变量的均值时,我们可以使用以下命令

1
2
sum price 
estadd r(mean)

我们也常常需要在回归表格中标注所控制的变量的类型,这一步骤也同样的可以使用estadd来实现。例如,当我们希望在上述部分回归中控制车辆是否是外国品牌时,我们可以在回归表格中加入一行来显示每个回归是否控制了该变量。

1
estadd local  Foreign_FE = " ", replace

完整代码,以及输出的表格效果可见下图。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
    sysuse auto,clear

reg price trunk headroom length weight
estimates store est1
sum price
estadd r(mean)
estadd local Foreign_FE = " ", replace

reg price trunk headroom length weight i.foreign
estimates store est2
sum price
estadd scalar ymean = r(mean)
estadd local Foreign_FE = "X", replace

. esttab, stats(Foreign_FE ymean r2 N, ///
> labels("Foreign FE" "Dep Mean" "R Square" "Num of Obs"))

--------------------------------------------
(1) (2)
price price
--------------------------------------------
trunk 114.1 63.31
(1.04) (0.68)

headroom -711.6 -606.2
(-1.60) (-1.62)

length -101.7* -89.03*
(-2.41) (-2.52)

weight 4.753*** 5.780***
(4.24) (6.04)

1.foreign 3526.8***
(5.51)

_cons 11488.5* 5340.2
(2.53) (1.35)
--------------------------------------------
Foreign FE X
Dep Mean 6165.3 6165.3
R Square 0.372 0.566
Num of Obs 74 74
--------------------------------------------
t statistics in parentheses
* p<0.05, ** p<0.01, *** p<0.001

制作tex文档

上述表格都是在STATA窗口中的显示效果,并不是我们的最终目的。我们希望能够将表格制作和文章写作打通起来,实现完全自动化,从进行回归分析,到将表格加入到文章中,不需要进行任何手动的复制粘贴。将这一过程自动化的目的,是在不断修改回归的过程中,减少人为出错的概率。为了达到这个目的,我们需要使用estout 中的using .tex, 以及prehead, postfoot这几个选项。using x.tex将输出结果更改为tex格式,而prehead, postfoot分别包括了使得tex文档能够直接编译成pdf的tex代码。具体代码可见下图。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
#delimit ; 

esttab est1 est2 using estout_eg.tex, replace style(tex) booktabs
keep(trunk headroom length weight 1.foreign)
stats(Foreign_FE ymean r2 N, labels("Foreign FE" "Dep Mean" "R-Squared" "N")
fmt(0 3 3 0)) // content options
b(3) se(3) star(* 0.1 ** 0.05 *** 0.01) // content formating options
mlabels(,none) numbers // column heading options
label varlabels(1.foreign "Foreign") // row head formating options
prehead( // tex code needed to compile document
\documentclass[leqno,11pt]{article}
\usepackage{booktabs}
\usepackage{tabularx}
\begin{document}
\begin{table}[h]
\def\sym#1{\ifmmode^{#1}\else\(^{#1}\)\fi}
\caption{Using Estout to Output Regression Coefficients}
\begin{center}
\begin{tabular}{l c c}
\toprule
)
postfoot( // more tex code needed to compile document
\end{tabular}
\end{center}
\end{table}
\end{document}
)
;
#delimit cr

STATA在执行完上述代码后,会生成一个estout_eg.tex的文档。通过Latex编译这个文档,我们可以得到如下图所示的PDF文档

代码下载链接
完整版本PDF下载链接

版权声明

Your browser is out-of-date!

Update your browser to view this website correctly. Update my browser now

×