Stata连享会 主页 || 视频 || 推文 || 知乎 || Bilibili 站
温馨提示: 定期 清理浏览器缓存,可以获得最佳浏览体验。
New!
lianxh
命令发布了:
随时搜索推文、Stata 资源。安装:
. ssc install lianxh
详情参见帮助文件 (有惊喜):
. help lianxh
连享会新命令:cnssc
,ihelp
,rdbalance
,gitee
,installpkg
⛳ Stata 系列推文:
作者:王硕 (山东大学)
邮箱:shweung@163.com
编者按:本文整理自 Types Of Transformations For Better Normal Distribution,特此致谢!
目录
以概率分布为核心的研究大都聚焦于正态分布。例如一般线性模型假定
在本推文中,我们将给出 2 个判断指标以及 5 种类型的数据转换方案,以更好地辨别和适应正态分布。
偏度衡量随机变量概率分布的不对称性,是相对于平均值不对称程度的度量。通过对偏度系数的测量,我们能够判定数据分布的不对称程度以及方向,其计算公式为
峰度是研究数据分布陡峭或平滑的统计量,通过对峰度系数的测量,我们能够判定数据相对于正态分布而言是更陡峭 / 平缓,其计算公式为
在查看数据偏度和峰度统计量时,可以借助 Stata 中的 tabstat
命令,并加入 count mean p50 skewness kurtosis
等选项呈现出观测值、均值、中位数、峰度和偏度等指标。以 Stata 软件自带的数据文件 nlsw88.dta 中的 wage 工资变量为例进行演示,该数据包含了 1988 年 2246 个美国妇女的资料,其命令的语法结构如下:
. sysuse nlsw88, clear
. tabstat wage, stat (count mean p50 skewness kurtosis)
Variable | N Mean p50 Skewness Kurtosis
----------+---------------------------------------------
wage | 2246 7.766949 6.27227 3.096199 15.85446
--------------------------------------------------------
从呈现结果可以看出,偏度为 3.096,均值 7.767 高于中位数 6.272,因此很明显工资变量是一个右偏的变量。峰度为 15.854,远大于标准正态分布的峰值 3,表明工作变量较为集中地分布在均值附近。
图形可以便于我们理解连续变量分布,从而直观地评估偏度和峰度,以建立对数据的初步了解。Stata 中的直方图可以较好地展示数据的偏度和峰度特征,其命令的语法结构如下:
. sysuse nlsw88,clear
. sum wage, d
. local p50 = r(p50)
. local mean = r(mean)
. twoway (hist wage, vertical color("136 187 203") ///
> xline(`mean', lcolor("220 146 20")) xline(`p50', lcolor(black))) ///
> (kdensity wage, lpattern(dash)), ///
> ylabel(0(0.035)0.14) graphregion(color(white)) ///
> xtitle(Wage,place(right)) ytitle(Density, place(top)) ///
> legend(ring(0) col(1) position(2.5))
. gr save wage0.gph, replace
直方图横轴表示数据类型,纵轴表示分布情况。从上图可以看出,均值 (橙色虚线) 位于中位数 (黑色虚线) 右侧,且右侧分布较多明显的极大值,这同样可以得出工资变量是一个右偏的变量的结论。从峰度来看,数据集中于 3-6 范围内,且高出同均值和方差的正态分布曲线,呈现出 “尖峰” 样式。
首先,采用由开发者提供的命令链接进行下载安装。
. ssc install jb6, replace
. sysuse nlsw88, clear
. jb6 wage
Jarque-Bera normality test: 1.9e+04 Chi(2) 0
Jarque-Bera test for Ho: normality: (wage)
从上述结果可以看出,该命令的原假设 H0 为数据服从正态分布,而结果中的
对数变换即将原始数据
. sysuse nlsw88,clear
. gen Logwage=log(wage)
. twoway (hist wage, vertical color("136 187 203")) ///
> (kdensity wage, lpattern(dash)), ///
> ylabel(0(0.035)0.14) graphregion(color(white)) ///
> xtitle(Wage,place(right)) ytitle(Density,place(top)) ///
> legend(ring(0) col(1) position(1))
. gr save wage.gph, replace
. twoway (hist Logwage,vertical color("96 133 149")) ///
> (kdensity Logwage, lpattern(dash)), ///
> graphregion(color(white)) ///
> xtitle(Logwage,place(right)) ytitle(Density,place(top)) ///
> legend(ring(0) col(1) position(1))
. gr save Logwage.gph, replace
. gr combine wage.gph Logwage.gph, ysize(3)
上图是原始数据 (左) 和对数转换数据 (右) 的比较。可以看出转换后的数据的偏度明显降低且峰度与正态分布接近。
平方根变换即将原始数据
. sysuse nlsw88,clear
. gen Sqrtwage = sqrt(wage)
. twoway (hist wage, vertical color("136 187 203")) ///
> (kdensity wage, lpattern(dash)), ///
> ylabel(0(0.035)0.14) graphregion(color(white)) ///
> xtitle(Wage,place(right)) ytitle(Density,place(top)) ///
> legend(ring(0) col(1) position(1))
. gr save wage.gph, replace
. twoway (hist Sqrtwage, vertical color("96 133 149")) ///
> (kdensity Sqrtwage, lpattern(dash)), ///
> graphregion(color(white)) ///
> xtitle(Sqrtwage,place(right)) ytitle(Density,place(top)) ///
> legend(ring(0) col(1) position(1))
. gr save Sqrtwage.gph, replace
. gr combine wage.gph Sqrtwage.gph, ysize(3)
上图是原始数据 (左) 和平方根转换数据 (右) 的比较。从上图可以看出,无论是从偏度还是峰度角度,平方根转换效果均低于对数转换效果。
在倒数转换中,
. sysuse nlsw88,clear
. gen Reciwage=1/wage
. twoway (hist wage, vertical color("136 187 203")) ///
> (kdensity wage, lpattern(dash)), ///
> ylabel(0(0.035)0.14) graphregion(color(white)) ///
> xtitle(Wage,place(right)) ytitle(Density,place(top)) ///
> legend(ring(0) col(1) position(1))
. gr save wage.gph, replace
. twoway (hist Reciwage, vertical color("96 133 149")) ///
> (kdensity Reciwage, lpattern(dash)), ///
> graphregion(color(white)) ///
> xtitle(Reciwage,place(right)) ytitle(Density,place(top)) ///
> legend(ring(0) col(1) position(1))
. gr save Reciwage.gph, replace
. gr combine wage.gph Reciwage.gph, ysize(3)
上图是原始数据 (左) 和倒数转换数据 (右) 的比较。从上图可以看出,无论是从偏度还是峰度角度,倒数转换未发生明显变动,其正态效果皆弱于对数转换和倒数转换。
Box-Cox 变换是一种广义幂变换方法,可以明显地改善数据的正态性和方差齐性,其计算公式为:
其中,
通过求解 boxcox
命令获得参数 boxcox
命令的语法结构如下:
boxcox depvar [indepvars] [if] [in] [weight] [, options]
depvar
和 indepvars
:被解释变量和解释变量。options
选择项如下:
model(lhsonly)
:指定只对因变量进行 Box-Cox 变换,而不对自变量进行变换。这种情况下,Box-Cox 转换的指数 model(rhsonly)
:指定只对解释变量进行 Box-Cox 变换,而不对被解释变量进行变换。这种情况下,Box-Cox 转换的指数 model(lambda)
:指定 Box-Cox 转换的指数 boxcox
命令的默认选项,它同时对因变量和自变量进行 Box-Cox 变换。model(theta)
:指定 Box-Cox 转换的指数 notrans
:指定不对数据进行任何变换,仅计算出基于因变量和自变量的联合最大似然估计 (MLE) 得到的 Box-Cox 变换的指数 下面仅对 wage 工资变量为例进行转换说明:
. sysuse nlsw88, clear
. boxcox wage, model(lhsonly)
Number of obs = 2,246
LR chi2(0) = 0.00
Log likelihood = -6108.0084 Prob > chi2 = .
----------------------------------------------------------------
wage | Coefficient Std. err. z P>|z| [95% conf. interval]
-----+-----------------------------------------------------------
/theta| -.2241236 .0288636 -7.76 0.000 -.2806952 -.167552
-----------------------------------------------------------------
Estimates of scale-variant parameters
----------------------------
| Coefficient
-------------+--------------
Notrans |
_cons | 1.502707
-------------+--------------
/sigma | .3727356
----------------------------
-----------------------------------------------------
Test Restricted LR statistic
H0: log likelihood chi2 Prob > chi2
-----------------------------------------------------
theta = -1 -6449.0767 682.14 0.000
theta = 0 -6138.868 61.72 0.000
theta = 1 -7117.2949 2018.57 0.000
-----------------------------------------------------
第一个表是 Box-Cox 变换参数
第二个表包含尺度变量参数的估计值,会显示各个自变量的参数估计值和模型标准差。本例只是对因变量进行了变化,该表实际意义不大;
第三个表包含了 Box-Cox 变换中三个常用函数的限制性似然比检验输出结果,即倒数变换 (
接下来,对 Box-Cox 变换参数
. sysuse nlsw88, clear
. gen Boxwage=(wage^(-0.2241)-1)/(-0.2241)
. twoway (hist wage, vertical color("136 187 203")) ///
> (kdensity wage, lpattern(dash)), ///
> ylabel(0(0.035)0.14) graphregion(color(white)) ///
> xtitle(Wage,place(right)) ytitle(Density,place(top)) ///
> legend(ring(0) col(1) position(1))
. gr save wage.gph, replace
. twoway (hist Boxwage, vertical color("96 133 149")) ///
> (kdensity Boxwage, lpattern(dash)), ///
> graphregion(color(white)) ///
> xtitle(Boxwage,place(right)) ytitle(Density,place(top)) ///
> legend(ring(0) col(1) position(10))
. gr save Boxwage.gph, replace
. gr combine wage.gph Boxwage.gph, ysize(3)
上图是原始数据 (左) 和 Box-Cox 转换数据 (右) 的比较。从上图可以看出,Box-Cox 转换的效果较好且与对数转化结果相似,实际上,对数转化模型是 Box-Cox 转换在模型选择时的一种情况 (
Yeo-Johnson 转换也是一种广义幂变换方法,通过构建一组单调函数对随机变量进行数据变换,其特点在于其可被应用于包含 0 值和负值的样本中,因此其也被认为是 Box-Cox 变换在实数域的推广。当随机变量非负时,其表达式为:
当随机变量为负时,其表达式为:
由于 Stata 中暂时没有专门处理 Yeo-Johnson 转换的命令,我们借助 Python 中的 scipy.stats 模块中的 yeojohnson 函数进行处理,相应命令为:
. sysuse nlsw88, clear
. python
---- python (type end to exit) ---
>>> from sfi import Data
>>> from scipy.stats import yeojohnson
>>> Wage=Data.get(va = 'wage')
>>> YJwage,lam = yeojohnson(Wage)
>>> print(YJwage,lam)
>>> Data.addVarFloat("YJwage")
>>> Data.store("YJwage",None,val=YJwage)
>>> end
-----------------------------------
. twoway (hist wage, vertical color("136 187 203")) ///
> (kdensity wage, lpattern(dash)), ///
> ylabel(0(0.035)0.14) graphregion(color(white)) ///
> xtitle(Wage,place(right)) ytitle(Density,place(top)) ///
> legend(ring(0) col(1) position(1))
. gr save wage.gph, replace
. twoway (hist YJwage, vertical color("96 133 149")) ///
> (kdensity YJwage, lpattern(dash)), ///
> graphregion(color(white)) ///
> xtitle(YJwage,place(right)) ytitle(Density,place(top)) ///
> legend(ring(0) col(1) position(10))
. gr save YJwage.gph, replace
. gr combine wage.gph YJwage.gph, ysize(3)
上图是原始数据 (左) 和 Yeo-Johnso 转换数据 (右) 的比较。从上图可以看出,Yeo-Johnso 转换的效果与 Box-Cox 转换非常相似,也是一种较好的转化方式。
通过前述讨论,我们有以下几方面的发现:
Note:产生如下推文列表的 Stata 命令为:
lianxh 对数 python交互, m
安装最新版lianxh
命令:
ssc install lianxh, replace
免费公开课
最新课程-直播课
专题 | 嘉宾 | 直播/回看视频 |
---|---|---|
⭐ 最新专题 | 文本分析、机器学习、效率专题、生存分析等 | |
研究设计 | 连玉君 | 我的特斯拉-实证研究设计,-幻灯片- |
面板模型 | 连玉君 | 动态面板模型,-幻灯片- |
面板模型 | 连玉君 | 直击面板数据模型 [免费公开课,2小时] |
⛳ 课程主页
⛳ 课程主页
关于我们
课程, 直播, 视频, 客服, 模型设定, 研究设计, stata, plus, 绘图, 编程, 面板, 论文重现, 可视化, RDD, DID, PSM, 合成控制法
等
连享会小程序:扫一扫,看推文,看视频……
扫码加入连享会微信群,提问交流更方便
✏ 连享会-常见问题解答:
✨ https://gitee.com/lianxh/Course/wikis
New!
lianxh
和songbl
命令发布了:
随时搜索连享会推文、Stata 资源,安装命令如下:
. ssc install lianxh
使用详情参见帮助文件 (有惊喜):
. help lianxh