楼主: 夏目贵志
6290 15

[Stata] 通过统计模拟学习Stata实用编程技巧(1): 异方差性 [推广有奖]

贵宾

学科带头人

96%

还不是VIP/贵宾

-

威望
1
论坛币
238580 个
通用积分
17117.6094
学术水平
851 点
热心指数
971 点
信用等级
711 点
经验
759232 点
帖子
4029
精华
1
在线时间
793 小时
注册时间
2012-7-15
最后登录
2017-9-16

相似文件 换一批

+2 论坛币
k人 参与回答

经管之家送您一份

应届毕业生专属福利!

求职就业群
赵安豆老师微信:zhaoandou666

经管之家联合CDA

送您一个全额奖学金名额~ !

感谢您参与论坛问题回答

经管之家送您两个论坛币!

+2 论坛币

大家好!今天我来挖一个大大的坑:通过统计模拟学习Stata实用编程技巧。这个系列呢,我也不知道要开几讲,其实也没计划,随性写点吧。不过至少还会有第二讲。

这个第一讲呢,目的是讲讲如何通过循环和宏变量收集数据并制表。那么用来作为例子的呢,是来研究一下异方差性对于估计量效率的影响。

具体地说,就是要比较同一个模型的同一个参数估计值在同方差时和异方差时的分布上的区别。但是,我们希望观察这个区别是否受样本量大小的影响。

考虑两个模型:

Y1=bX+u, uiid标准正态分布

Y2=bX+v, v是相互独立的正态分布,但是每个u_i的方差是x_i

我们要观察bOLS估计值在不同样本量下的分布。样本量采用20, 50,100, 200, 500。每个样本量做500次模拟。

这个模拟要求输出如下格式的表格(xxx表示数字):

statatrick1.png

我们如何用Stata编程实现呢?



首先,我们来用简单的语言描述一下这个程序需要干些什么

1.      按照上述要求生成数据

2.      跑上述两个模型

3.      储存每个模型b的估计值

4.      重复第1到第3步500次。

5.      重复第1到第4步5次,每次使用一个不同的样本量(20,50, 100, 200, 和500)。

好,那么为了实现上述功能,我们都有可能遇到什么问题呢?该用什么命令呢?

第一步:

生成数据的时候需要指定样本的大小。这个通过set obs 实现。生成变量使用generate命令。生成正态分布的随机变量使用rnormal()函数。

第二步:

跑回归模型使用regress命令。

第三步:

通过_b[x]来取得参数的估计值。使用local macro来保存这些值。注意:每一次模拟的值都需要保存,所以任何一次回归的结果都不能被之后的结果覆盖。也就是说,用于保存这些估计值的localmacro必须全部都有不同的名字。

第四步:

重复500次使用forvalues循环。

第五步:

重复5次,但是每次使用不同的样本量。这个需要使用foreach循环。



下面我们来看看写好的统计模拟的程序。

statatrick2.png

2行使用set seed来初始化伪随机数生成器,这样可以保证结果可以被他人重现。

3行和第4行把重复模拟的次数和样本量的选择保存在局部宏里,方便修改。

5行执行的是上述第5步。

6行注意不能使用clear all,不然会消除之前保存的所有结果。clear只消除数据。不消除其他保存的结果。

7行设定样本量。

8行注意:我们假定x的取值是给定的,非随机的。

9行执行上述第4步,开始模拟。

1014行生成数据。注意,b的真实值为3

1517行跑两个回归。注意,按照题目的要求,回归不包括常数项。

1618行保存参数的估计值。宏的名称包括两个变量,第一个表示样本量,第二个表示第几次模拟。这两个部分缺一不可。否则就无法保存每次回归的结果。原则是,有几次循环,保存结果的宏名字里就要有几个部分。这里有两个循环,一个foreach,一个forvalues,那么这里的宏名里就有两个变量。




下面我们来写保存和输出结果的程序。

statatrick3.png

这个部分是这一讲的关键所在,即,如何保存结果以便进一步处理或制表输出。

基本的思路是这样:

首先把所有要写进表里的数据保存在宏里。然后按照输出表格的格式建立一个数据集。接下来把保存的结果一个一个的写进表里。但是要注意的是我们需要保存所有的原始数据,而不是仅保留经过处理的数据。所以,最终的目标是保存一个数据集,让这个数据集可以通过table命令生成我们最终需要报告的表格。

我们这次的题目要求输出估计值的均值和标准差。而均值和标准差都是可以通过table命令计算的。这样我们就可以直接保存每次的估计值,而不用在保存之前就手动计算均值和标准差。

在上一步(程序的120行),我们已经将所有的估计值保存在宏里了。这里我们直接就可以建立新表。

21行清除现有数据。注意,千万不要clear all

2224行计算并设定我们需要的观测值数量。第22行通过数ssizes里的词数来确定我们有多少不同的样本量设定。对于每一个样本量设定,我们需要reps个观测值。这个乘积还需要乘以2,因为每次模拟我们有两个估计值需要保存(一个同方差的一个异方差的)。

2527行生成空变量。我们总共需要4个空变量。一个表示模型,一个表示样本量,一个表示模拟的次数,一个用来实际保存我们的估计值。

28行生成一个计数宏。因为我们每次写入这个数据的一行,所以我们需要知道我们下一步应该写入第几行。这个计数宏从1开始,到_N结束。每写入一行数据就给这个计数宏加1(见第36行)。

2939行实际写入数据。注意3235行使用的in条件。这就是我们需要一个计数宏的理由。

40行使用table命令输出结果。注意这一行的开头是用了noi命令,保证输出的结果会显示在屏幕上。




输出的结果如下所示。在Stata的结果页面选择Copy as table,然后粘贴到Excel里,简单的整理一下就可以粘贴到Word里输出啦!

statatrick4.png

最后呢,输出结果的长相是这样的:

statatrick5.png



最后留一个思考题:这个结果到底应该如何解读呢?异方差性对估计值效率的影响如何随着样本量增加变化呢?



二维码

扫码加我 拉你入群

请注明:姓名-公司-职位

以便审核进群资格,未注明则拒绝

关键词:Stata 统计模拟 异方差性 tata 异方差 技巧 统计

已有 3 人评分经验 论坛币 学术水平 热心指数 信用等级 收起 理由
我的素质低 + 100 + 50 + 5 精彩帖子
xddlovejiao1314 + 100 + 100 + 5 + 5 + 5 精彩帖子
niuniuyiwan + 60 + 5 + 5 + 5 精彩帖子

总评分: 经验 + 260  论坛币 + 150  学术水平 + 10  热心指数 + 15  信用等级 + 10   查看全部评分

本帖被以下文库推荐

沙发
夏目贵志 发表于 2015-8-8 11:09:34 |只看作者 |坛友微信交流群
完整的程序在这里
  1. clear all
  2. set seed 1172
  3. local reps "500"
  4. local ssizes "20 50 100 200 500"
  5. foreach ssize in `ssizes' {
  6.         clear
  7.         set obs `ssize'
  8.         gen x = runiform()*10
  9.         forvalues i=1/`reps' {
  10.                 cap drop u v y1 y2
  11.                 gen u = rnormal()
  12.                 gen v = rnormal()*x^2
  13.                 gen y1 = 3*x+u
  14.                 gen y2 = 3*x+v
  15.                 qui reg y1 x, noconst
  16.                 local b1_`ssize'_`i' = _b[x]
  17.                 qui reg y2 x, noconst
  18.                 local b2_`ssize'_`i' = _b[x]
  19.         }
  20. }
  21. clear
  22. local set = wordcount("`ssizes'")
  23. local rsize = `reps'*`set'*2
  24. set obs `rsize'
  25. foreach var in iter ssize model b {
  26.         gen `var' = .
  27. }
  28. local t = 1
  29. foreach ssize in `ssizes' {
  30.         forvalues i=1/`reps' {
  31.                 forvalues m = 1/2 {
  32.                         replace ssize = `ssize' in `t'
  33.                         replace iter = `i' in `t'
  34.                         replace b = `b`m'_`ssize'_`i'' in `t'
  35.                         replace model = `m' in `t'
  36.                         local t = `t'+1
  37.                 }
  38.         }
  39. }
  40. noi table ssize, contents(mean b sd b ) by(model) format(%5.3f)
复制代码


大家看了要是有什么错误或者有什么更好的做法请一定提出来!
已有 3 人评分经验 论坛币 学术水平 热心指数 信用等级 收起 理由
Sunknownay + 2 + 2 + 2 精彩帖子
xddlovejiao1314 + 50 + 20 + 1 + 1 + 1 精彩帖子
niuniuyiwan + 5 + 1 精彩帖子

总评分: 经验 + 50  论坛币 + 25  学术水平 + 3  热心指数 + 4  信用等级 + 3   查看全部评分

使用道具

藤椅
夏目贵志 发表于 2015-8-9 03:36:13 |只看作者 |坛友微信交流群
大家是觉得这样的内容太简单了吗?完全没有人气啊。。。
其实很多输出结果之类的问题都可以通过类似的方式解决。这样的的编程任务实际应用中是很多的。这里只是使用了一个很简单的例子。
要是没人气的话就腰斩了吧。。。
已有 1 人评分论坛币 热心指数 收起 理由
niuniuyiwan + 5 + 1 热心帮助其他会员

总评分: 论坛币 + 5  热心指数 + 1   查看全部评分

使用道具

感谢楼主写的那么用心,但是我完全没学过这些,所以看不懂
已有 2 人评分经验 论坛币 收起 理由
xddlovejiao1314 + 10 + 3 鼓励积极发帖讨论
niuniuyiwan + 5 精彩帖子

总评分: 经验 + 10  论坛币 + 8   查看全部评分

使用道具

报纸
liuding1111 发表于 2015-8-9 11:28:01 |只看作者 |坛友微信交流群
楼主辛苦了,楼主的劳动成果对于我们学习计量知识和Stata编程帮助很大,楼主要继续!
已有 2 人评分经验 论坛币 收起 理由
xddlovejiao1314 + 10 + 3 鼓励积极发帖讨论
niuniuyiwan + 5 精彩帖子

总评分: 经验 + 10  论坛币 + 8   查看全部评分

使用道具

地板
夏目贵志 发表于 2015-8-25 21:36:29 |只看作者 |坛友微信交流群
放在Stata版里效果不好。转来代码库吧。

使用道具

7
niuniuyiwan 在职认证  发表于 2015-8-25 21:55:01 |只看作者 |坛友微信交流群
这个系列非常好,也是花了不少心血的,感谢楼主,感谢分享。
已有 1 人评分经验 论坛币 收起 理由
xddlovejiao1314 + 10 + 3 鼓励积极发帖讨论

总评分: 经验 + 10  论坛币 + 3   查看全部评分

使用道具

8
BlackHawk123 在职认证  发表于 2015-8-28 11:02:17 |只看作者 |坛友微信交流群
感谢分享~~~
已有 1 人评分经验 论坛币 收起 理由
xddlovejiao1314 + 10 + 3 鼓励积极发帖讨论

总评分: 经验 + 10  论坛币 + 3   查看全部评分

使用道具

9
xddlovejiao1314 学生认证  发表于 2015-8-28 11:24:35 来自手机 |只看作者 |坛友微信交流群
夏目贵志 发表于 2015-8-8 11:07
大家好!今天我来挖一个大大的坑:通过统计模拟学习Stata实用编程技巧。这个系列呢,我也不知道要开几讲,其 ...
呵呵,好帖。谢谢分享。欢迎多来代码库分享好帖。
已有 1 人评分论坛币 收起 理由
niuniuyiwan + 4 我很赞同

总评分: 论坛币 + 4   查看全部评分

使用道具

10
tkt718 发表于 2015-9-5 00:26:24 来自手机 |只看作者 |坛友微信交流群
夏目贵志 发表于 2015-8-8 11:09
完整的程序在这里


感谢楼主分享!新手入门需要消化消化!
已有 2 人评分经验 论坛币 热心指数 收起 理由
xddlovejiao1314 + 10 + 3 鼓励积极发帖讨论
niuniuyiwan + 5 + 1 精彩帖子

总评分: 经验 + 10  论坛币 + 8  热心指数 + 1   查看全部评分

使用道具

您需要登录后才可以回帖 登录 | 我要注册

本版微信群
加好友,备注jltj
拉您入交流群

京ICP备16021002-2号 京B2-20170662号 京公网安备 11010802022788号 论坛法律顾问:王进律师 知识产权保护声明   免责及隐私声明

GMT+8, 2024-4-28 18:26