Stata:蒙特卡洛模拟新命令-simulate2

发布时间:2022-04-25 阅读 227

Stata连享会   主页 || 视频 || 推文 || 知乎 || Bilibili 站

温馨提示: 定期 清理浏览器缓存,可以获得最佳浏览体验。

New! lianxh 命令发布了:
随时搜索推文、Stata 资源。安装:
. ssc install lianxh
详情参见帮助文件 (有惊喜):
. help lianxh
连享会新命令:cnssc, ihelp, rdbalance, gitee, installpkg

课程详情 https://gitee.com/lianxh/Course

⛳ Stata 系列推文:

PDF下载 - 推文合集

作者:李坤 (华南师范大学)
邮箱kunli1001@163.com


目录


1. 引言

蒙特卡洛模拟 (MC) 是指当根据总体的分布函数很难求出想要的数字特征时,从总体中生成大量的样本,使用样本的数字特征估计总体的数字特征。Stata 有两种常用的方法进行 MC 模拟,第一种是使用 postfile 命令,另一种是使用simulate 命令。

在本文中,我们为大家介绍一个新的高效的 MC 模拟命令 simulate2,以及另外一个可以进行并行实验的命令 psimulate2,这两个命令可以大大简化 MC 模拟的程序并提高运行的速度。

2. 命令介绍

simulate2 命令安装:

ssc install simulate2, replace

需要注意,simulate2 命令需要 Stata16 及以上才能正常使用;psimulate2 命令需要 Stata15 及以上才能正常使用,目前仅适用于 Windows 系统。

simulate2 命令语法:

simulate2 [exp_list] , reps(#) [options1 options2 ] : command

其中,[exp_list] 是指将每一次 command 的执行结果按照 [exp_list] 格式提取出来,若不是 [exp_list] 格式则会报错。command 代表每一次蒙特卡洛模拟执行的程序,可以是 Stata 自带或者外部命令,也可以是用户自己编写的程序。rep(#) 代表做 # 次蒙特卡洛模拟。

[options1] 主要选项如下:

  • nodots:在 MC 模拟过程中不在屏幕上打点;
  • dots(#):每隔 # 次模拟,在屏幕上打印一个点;
  • noisily:将每次 MC 的结果都显示在屏幕上;
  • trace:追踪命令运行的过程;
  • nolegend:不显示 MC 的信息;
  • verbose:显示 MC 信息。

[options2] 主要选项如下:

  • aving(filename, …):将结果保存到一个文件夹中;
  • seed(options):设定随机种子数;
  • seedsave(options):保存用过的种子。

psimulate2 命令语法:

psimulate2 [exp_list] , reps(#) parallel(#2, options3) [options2] : command

其中,parallel(#2,()) 中的 #2 是设置并行运行 Stata 实例的数量,建议使用的实例数量不要超过可用的 CPU 内核。options3 主要选项如下:

  • exe(string):设置 Stata.exe 的路径,若 Stata.exe 位于非常规文件夹中或具有非常规文件名,可能会查找失败;
  • temppath(string):设置保存临时文件的替代路径;
  • processors(integer):设置每个 Stata 实例允许的最大处理数,仅适用于 Stata MP 版本;
  • simulate:使用 simulate 而不是 simulate2,如果使用 Stata15 运行 psimulate2,则默认使用 simulate
  • seed(#):设置随机种子数;
  • seedstream(integer):种子流,可以为第一个实例设置初始种子流编号。例如 parallel(3, seedstream(4)) 表示实例 1 将使用种子流编号为 4、实例 2 将使用种子流编号为 5,依次类推。该选项只能在 psimulate2 中使用。

3. 具体案例

3.1 simulate2 命令

我们先通过 1000 次样本容量为 100 的 MC 模拟实验来估计 OLS 的参数和标准误,并记录时间。在生成的 dta 文件中,每一行的数据代表一次蒙特卡洛模拟,每一列的数据代表一个基于样本计算出的统计值,最后可以通过计算上述数据的均值来近似估计回归方程的系数和标准误。

在上述过程中,我们将种子保存到名为 seed_frame 的 frame 中。在 Stata16 之前的版本中,Stata 的内存只能存储一份数据集,这在实际中非常的不方便。frame 则可以允许在 Stata 内存中同时存储、操作多份数据,并且可以将这些数据连接起来,极大的提升了数据处理的效率。

* 定义程序 testsiuml,并以 r() 形式储存结果
program define testsimul, rclass 
    version 17       // 使用的Stata的版本
    syntax anything
    clear
    set obs `anything' 
    gen x = rnormal(1,4)   
    gen y = 2 + 3*x + rnormal()  // 假设真实的模型参数为 a=2 b=3
    reg y x                      // 用 reg y x 命令来拟合模拟
    matrix b = e(b)              // 提取参数估计值
    matrix se = e(V)             // 提取标准误的估计值
    ereturn clear
    return scalar b = b[1,1]
    return scalar V = se[1,1]
    return local time "`c(current_time)'"
end

* 选取 100 个样本进行 1000 次 MC 模拟        
simulate2 time = r(time) b = r(b) V = r(V), reps(1000) seedsave(seed_frame, frame): testsimul 100

现在我们可以利用保存的种子重新进行前 500 次实验:

. simulate2 time = r(time) b = r(b) V = r(V), reps(500) seed(seed_frame seed, frame): testsimul 100

对第二个 500 次重复实验,将起始种子编号设置为 501:

. simulate2 time = r(time) b = r(b) V = r(V), reps(500) ///
>     seed(seed_frame seed, frame start(501)): testsimul 100

与上述的步骤相反,我们也可以先进行 500 次抽样保存种子,再进行 500 次抽样追加保存种子。最后对所有 1000 次抽样进行实验,并将比较结果保存在 frame 中。

. simulate2 time = r(time) b = r(b) V = r(V), reps(500) seed(123) ///
>     seedsave(seed_frame, frame) saving(first500, frame): testsimul 100

. * 对于第二次运行,没有设置新种子,而是使用了选项 append 
. simulate2 time = r(time) b = r(b) V = r(V), reps(500) ///
>     seedsave(seed_frame, frame append) saving(second500, frame): testsimul 100  

. simulate2 time = r(time) b = r(b) V = r(V), reps(1000) ///
>     seed(seed_frame seed, frame): testsimul 100

分别比较用两种方法计算得到的前 500 次重复模拟结果与后 5​​00 次重复模拟的结果:

. sum if _n <= 500

    Variable |        Obs        Mean    Std. dev.       Min        Max
-------------+---------------------------------------------------------
           b |        500    3.001443     .024907   2.921489   3.071222
           V |        500    .0006364    .0001251    .000357   .0010889
        time |          0

. frame first500: sum

    Variable |        Obs        Mean    Std. dev.       Min        Max
-------------+---------------------------------------------------------
           b |        500    3.001443     .024907   2.921489   3.071222
           V |        500    .0006364    .0001251    .000357   .0010889
        time |          0

. sum if _n > 500

    Variable |        Obs        Mean    Std. dev.       Min        Max
-------------+---------------------------------------------------------
           b |        500    2.999219    .0241378   2.930022   3.078007
           V |        500    .0006446    .0001321   .0003689   .0012382
        time |          0

. frame second500: sum

    Variable |        Obs        Mean    Std. dev.       Min        Max
-------------+---------------------------------------------------------
           b |        500    2.999219    .0241378   2.930022   3.078007
           V |        500    .0006446    .0001321   .0003689   .0012382
        time |          0

可以看出两个命令的结果是相同的,表明最终结果与进行实验的方式是无关的。同时观察数据可以看出,估计出的值都接 3,与我们设定的回归方程的系数非常接近。

3.2 psimulate2 命令

psimulate2 命令是 simulate2 的并行版本。它的工作原理是将复制的数据分成相等的块,每个块在单独的 Stata 实例上运行。为此 psimulate2 创建了一个 do 文件和一个批处理文件。然后使用批处理文件来启动一个运行相应 do 文件的新 Stata 实例。正在运行的实例充当父实例,另一个是子实例。

psimulate2 在 Stata 上输出与 simulatesimulate2 的输出不同。它显示已经运行完成的百分比、经过的时间和预期的剩余时间以及完成时间。因此,我们可以将上面的实验使用 psimulate2 进行并行化模拟。例如我们要同时运行两个实例,即实例 1 运行第一个 500 次重复,实例 2 运行第二个 500 次重复,采用下述的代码:

. psimulate2, reps(500) seed(123) parallel(2, temppath("D:\Stata17\do")): testsimul 200	 
. sum  // 查看结果       

    Variable |        Obs        Mean    Std. dev.       Min        Max
-------------+---------------------------------------------------------
           V |        500    .0003177    .0000478   .0002151   .0005266
           b |        500    2.999787    .0187097   2.946175    3.04934
        time |          0

parallel(2, temppath("D:\Stata17\do")) 含义为设置两个并行实例,temppath() 指定保存临时文件的路径。仅当 Windows 不允许从临时文件夹运行批处理文件时,才需要此选项。运行第一个和第二个 500 次重复具有相同的随机数状态,但种子流不同,因此随机抽取也不同。可以看出,该命令最后的结果与使用 simulate2 命令计算出的结果都是接近的,但是运行速度更快,语法结构也更加简洁。

当我们需要将一个或多个 psimulate2 嵌套在循环中或按顺序调用,它们将使用相同的初始种子,可能造成最终的结果出现偏误。因此为了避免这种情况,但还能按顺序运行 psimulate2,可以使用以下代码来让 psimulate2 返回 r() 中最后一个实例的最后一个种子状态。

. psimulate2, reps(100) seed(_current) p(2) seedsave(seed, frame): testsimul 200 
. set rng `r(rng_current)'
. set rngstate `r(rngstate)'
. psimulate2, reps(100) seed(_current) p(2) seedsave(seed, frame): testsimul 200
. sum

    Variable |        Obs        Mean    Std. dev.       Min        Max
-------------+---------------------------------------------------------
           V |        100    .0003271    .0000516   .0002419   .0004835
           b |        100    3.001315    .0181299   2.961705   3.052167
        time |          0

psimulate2, seed() 中使用 current,可以将使用 current 父实例的种子作为所有子实例的初始种子,这样每个子实例仍将具有不同的种子流,以确保随机数抽取不同。

3.3 simulate 与 simulate2 的区别

simulate2 遵循 simulate 的语法并简化了执行蒙特卡洛模拟的编程任务。相比于 simulate,它的新颖之处在于允许其调用的程序将宏 (字符串) 返回到 r()e() 来扩展命令,然后将它们保存到使用的数据集中。使用该命令时,还可以产生一个当下时间的数据集,记录该结果产生的时间。

此外,它还具有高级选项,可将随机数状态 (种子、种子流和类型) 分配给特定的模拟重复,种子可以保存到 Stata 数据集 (dta) 或 frame 中,以便无缝复制模拟结果,模拟结果也可以直接保存到 frame 中。

4. 结语

通过上面的介绍,向大家展示了 simulate2 以及 psimulate2 的具体使用方法。可以发现当需要进行 MC 模拟时,simulate2psimulate2 命令拥有很大的优势,不仅语法简洁,而且使用起来也很方便。我们希望通过这篇文章,能够让大家对这两个命令有一个更加深入的认识。

5. 相关推文

Note:产生如下推文列表的 Stata 命令为:
lianxh MC 模拟, m
安装最新版 lianxh 命令:
ssc install lianxh, replace

相关课程

免费公开课

最新课程-直播课

专题 嘉宾 直播/回看视频
最新专题 文本分析、机器学习、效率专题、生存分析等
研究设计 连玉君 我的特斯拉-实证研究设计-幻灯片-
面板模型 连玉君 动态面板模型-幻灯片-
面板模型 连玉君 直击面板数据模型 [免费公开课,2小时]
  • Note: 部分课程的资料,PPT 等可以前往 连享会-直播课 主页查看,下载。

课程主页

课程主页

关于我们

  • Stata连享会 由中山大学连玉君老师团队创办,定期分享实证分析经验。
  • 连享会-主页知乎专栏,700+ 推文,实证分析不再抓狂。直播间 有很多视频课程,可以随时观看。
  • 公众号关键词搜索/回复 功能已经上线。大家可以在公众号左下角点击键盘图标,输入简要关键词,以便快速呈现历史推文,获取工具软件和数据下载。常见关键词:课程, 直播, 视频, 客服, 模型设定, 研究设计, stata, plus, 绘图, 编程, 面板, 论文重现, 可视化, RDD, DID, PSM, 合成控制法

连享会小程序:扫一扫,看推文,看视频……

扫码加入连享会微信群,提问交流更方便

✏ 连享会-常见问题解答:
https://gitee.com/lianxh/Course/wikis

New! lianxhsongbl 命令发布了:
随时搜索连享会推文、Stata 资源,安装命令如下:
. ssc install lianxh
使用详情参见帮助文件 (有惊喜):
. help lianxh