使用 Stata 和 R 绘制相同的随机数答案

【问题标题】：Draw same random numbers with Stata and R使用 Stata 和 R 绘制相同的随机数
【发布时间】：2020-02-03 04:55:12
【问题描述】：

我用Stata 和R 绘制相同的随机数。本质上，我想获得与R 中的sample 和Stata 中的rdiscrete 相同的随机数系列。不过，我已尝试为每种语言提供一个完整但小而可重复的示例。

我认为sample 函数与rdiscrete 函数做同样的事情，但我不确定。假设这些函数做同样的事情，我只需要它们返回相同的随机数。

我正在使用Stata 12。

这是我的R 代码：

set.seed(1234)

wave_of_cy  = 2
wave_obs = 20

fake_dat <- read.table(text = '
     nobs  p1   p2
      0   .20  .10
      1   .10  .15
      2   .10  .15
      3   .05  .10
      4   .05  .10
      5   .20  .05
      6   .10  .05
      7   .05  .05
      8   .05  .05
      9   .10  .20
', header = TRUE, stringsAsFactors = FALSE)

p_hrand  = fake_dat[, (wave_of_cy+1)]
pp_hrand = p_hrand / sum(p_hrand)

my_rdata = sample(nrow(fake_dat), wave_obs, prob=pp_hrand, replace = TRUE)
my_rdata

hrand    = fake_dat[my_rdata, 1]
hrand

这是我的Stata 代码：

clear
set seed 1234
global wave_of_cy  = 2
set obs 20
local wave_obs = _N

clear
input nobs p1 p2
0 .20 .10
1 .10 .15
2 .10 .15
3 .05 .10
4 .05 .10
5 .20 .05
6 .10 .05
7 .05 .05
8 .05 .05
9 .10 .20
end
list
save fake_dat

clear

use "fake_dat.dta", replace
putmata fake_data = (nobs p1 p2), replace

mata:
     p_hrand  = fake_data[., $wave_of_cy+1]
     pp_hrand = p_hrand :/ sum(p_hrand)
     my_rdata = rdiscrete(`wave_obs', 1, pp_hrand)
     my_rdata
     hrand    = fake_data[my_rdata, 1]
     hrand
end

【问题讨论】：

不幸的是，即使跨平台使用相同的###，您也无法复制随机种子。 Python 和 R、Matlab 和 Python、SAS 和 Stata 等都曾问过这个问题。它们都使用不同的算法。如果需要重复使用相同的数据，请尝试将数据保存到磁盘（csv、txt 等）或通过命令行 i/o 传递。
查看 www.random.org 了解跨平台解决方案。
实际上，@Parfait 的建议是人们在 99% 的时间里都会这样做。如果由于某种原因您确实不能或不会通过 IO 传递数据，则需要让 R & Stata 调用相同的外部代码，或者让 R 调用 Stata 代码，反之亦然。

标签： r stata random-seed

【解决方案1】：

如上所述，跨软件/语言的随机生成不容易复制，因为即使种子数相同，每种软件/语言也会运行不同的算法。为了重现相同的随机生成，您需要连接两个平台：

使用双语言 API（例如，rpy2 在 Python 中运行 R，reticulate 在 R 中运行 Python，或 twister 在 Matlab 中运行 Python 的 random.random()）
李>
运行像 C/C++ 这样的低级语言，以便在两个软件之间的应用层调用，例如 SAS and Stata;

这种方法在这里是可行的，因为 R 是用 C、Fortran 和 R 编写的，而 Stata（作为一种软件而不是语言）是用 C 编写的，因此两者都可以调用相同的随机数算法；
在任一平台上运行命令行并通过 i/o 文本处理导出/导入结果数据。

下面演示了最后一个选项。

R （在batch mode 中调用Stata，假定在最后一行end 之后没有空行）

setwd("C:\\Path\\To\\Working\\Directory")
# RUN DO SCRIPT WHICH OUTPUTS LOG OF SAME NAME
system("C:\\Path\\To\\StataMP-64.exe /e do myStataScript.do")

# READ IN LOG FILE TO CHARACTER VECTOR
stata_log <- readLines("myStataScript.log")

# EXTRACT NEEDED hrand OUTPUT LINES (N=20)
stata_data <- stata_log[(length(stata_log)-26):(length(stata_log)-7)]

# MATRIX BUILD OF EXTRACT AND RETURN SECOND ROW (TO MIRROR STATA'S RESULTS)
sapply(strsplit(stata_data, "\\|"), as.integer)[2,]
# [1] 9 9 1 9 0 9 4 1 0 2 2 2 0 6 2 7 1 5 3 1

Stata （调用Rscript自动执行）

首先在 R 示例脚本中添加所需的行：

setwd("C:\\Path\\To\\Working\\Directory")

... original code ...

# SAVE hrand DATA TO DISK
write.csv(data.frame(hrand), "RandomSeedDataSample.csv", row.names = FALSE)

然后运行Stata脚本：

* RUN R SCRIPT
shell "C:\Path\To\R\bin\Rscript.exe" "C:\Path\myRScript.R"

* IMPORT CSV FILE
import delimited using "C:\Path\To\Working\Directory\RandomSeedDataSample.csv", clear

* MATRIX BUILD (TO MIRROR R'S RESULTS)
putmata hrand = (hrand), replace

mata
    hrand
end

:         hrand
        1
     +-----+
   1 |  9  |
   2 |  3  |
   3 |  3  |
   4 |  3  |
   5 |  5  |
   6 |  3  |
   7 |  9  |
   8 |  2  |
   9 |  3  |
  10 |  4  |
  11 |  3  |
  12 |  4  |
  13 |  2  |
  14 |  8  |
  15 |  2  |
  16 |  6  |
  17 |  2  |
  18 |  2  |
  19 |  9  |
  20 |  2  |
     +-----+

【讨论】：