Stan 中指数随机变量的模拟 (RStan Package/Interface)
Simulations of Exponential Random Variables in Stan (RStan Package/Interface)
我正在尝试使用 RStan 代码模拟指数随机变量。
模拟指数随机变量的 R 代码如下所示:
A <- rexp(1000, 4)
B <- rexp(1000, lambda2)
C <- rexp(1000, lambda3)
D <- rexp(1000, lambda4)
E <- rexp(1000, lambda5)
F <- rexp(1000, lambda6)
A + B = AB
C + D = CD
mean(AB)
mean(CD)
如你所见,在R中模拟指数随机变量对我来说非常简单。以A <- rexp(1000, 4)
为例:它生成一个包含1000个结果的随机样本,其中lambda = 4。我可以然后对这些模拟值做统计分析(求均值等)。
现在我想在使用 RStan 时做同样的事情。使用 RStudio 的 R 笔记本功能,可以 "insert" 不同类型的代码,其中之一是 Stan:
我有以下 Stan 和 R 代码:
{stan output.var="exponential"}
generated quantities{
real total;
real number;
number = 0.0;
total = 0.0;
while(total < 1000) {
total += exponential_rng(4);
number += 1.0;
}
number -= 1.0;
}
{r}
simul2 <- sampling(exponential, algorithm="Fixed_param")
{r}
print(simul2, pars=c("number"), digits = 5)
但是,当我执行此代码(具体来说,print
命令)时,我得到以下输出:
我不明白为什么我会得到如此荒谬的统计数据(平均 4000!?)。我的目标是能够获得与我在 R 中进行分析时相同的统计数据,正如我在上面所展示的那样;换句话说,我的目标是获得与如果我完成 A <- rexp(1000, 4)
时相同的值,在这种特定情况下,以及更普遍的 X <- rexp(1000, lambda)
。
显然,我的 Stan 代码不正确,所以我将不胜感激,请大家花时间解释一下正确的使用方法。
那是告诉你 number
的后验分布,这只是你的计数器。你的 Stan 代码做的事情与你的 R 代码不同,但是如果你传递 N
和 lambdas 作为数据,那么做你的 R 代码的模拟并不难:
data {
int<lower=1> N;
real<lower=0> lambda2;
real<lower=0> lambda3;
real<lower=0> lambda4;
}
generated quantities {
real mean_AB;
real mean_CD;
{
vector[N] A;
vector[N] B;
vector[N] C;
vector[N] D;
for (n in 1:N) {
A[n] = exponential_rng(4);
B[n] = exponential_rng(lambda2);
C[n] = exponential_rng(lambda3);
D[n] = exponential_rng(lambda4);
}
mean_AB = mean(A + B);
mean_CD = mean(C + D);
}
}
我正在尝试使用 RStan 代码模拟指数随机变量。
模拟指数随机变量的 R 代码如下所示:
A <- rexp(1000, 4)
B <- rexp(1000, lambda2)
C <- rexp(1000, lambda3)
D <- rexp(1000, lambda4)
E <- rexp(1000, lambda5)
F <- rexp(1000, lambda6)
A + B = AB
C + D = CD
mean(AB)
mean(CD)
如你所见,在R中模拟指数随机变量对我来说非常简单。以A <- rexp(1000, 4)
为例:它生成一个包含1000个结果的随机样本,其中lambda = 4。我可以然后对这些模拟值做统计分析(求均值等)。
现在我想在使用 RStan 时做同样的事情。使用 RStudio 的 R 笔记本功能,可以 "insert" 不同类型的代码,其中之一是 Stan:
我有以下 Stan 和 R 代码:
{stan output.var="exponential"}
generated quantities{
real total;
real number;
number = 0.0;
total = 0.0;
while(total < 1000) {
total += exponential_rng(4);
number += 1.0;
}
number -= 1.0;
}
{r}
simul2 <- sampling(exponential, algorithm="Fixed_param")
{r}
print(simul2, pars=c("number"), digits = 5)
但是,当我执行此代码(具体来说,print
命令)时,我得到以下输出:
我不明白为什么我会得到如此荒谬的统计数据(平均 4000!?)。我的目标是能够获得与我在 R 中进行分析时相同的统计数据,正如我在上面所展示的那样;换句话说,我的目标是获得与如果我完成 A <- rexp(1000, 4)
时相同的值,在这种特定情况下,以及更普遍的 X <- rexp(1000, lambda)
。
显然,我的 Stan 代码不正确,所以我将不胜感激,请大家花时间解释一下正确的使用方法。
那是告诉你 number
的后验分布,这只是你的计数器。你的 Stan 代码做的事情与你的 R 代码不同,但是如果你传递 N
和 lambdas 作为数据,那么做你的 R 代码的模拟并不难:
data {
int<lower=1> N;
real<lower=0> lambda2;
real<lower=0> lambda3;
real<lower=0> lambda4;
}
generated quantities {
real mean_AB;
real mean_CD;
{
vector[N] A;
vector[N] B;
vector[N] C;
vector[N] D;
for (n in 1:N) {
A[n] = exponential_rng(4);
B[n] = exponential_rng(lambda2);
C[n] = exponential_rng(lambda3);
D[n] = exponential_rng(lambda4);
}
mean_AB = mean(A + B);
mean_CD = mean(C + D);
}
}