如何在 rmarkdown 中使用 stan
How to use stan in rmarkdown
我想在 rnotebook
中使用 rstan
获得模型的估计系数
我有以下 stan
块:
```{stan output.var="rats"}
data {
int<lower=0> N;
int<lower=0> T;
real x[T];
real y[N,T];
real xbar;
}
parameters {
real alpha[N];
real beta[N];
real mu_alpha;
real mu_beta; // beta.c in original bugs model
real<lower=0> sigmasq_y;
real<lower=0> sigmasq_alpha;
real<lower=0> sigmasq_beta;
}
transformed parameters {
real<lower=0> sigma_y; // sigma in original bugs model
real<lower=0> sigma_alpha;
real<lower=0> sigma_beta;
sigma_y = sqrt(sigmasq_y);
sigma_alpha = sqrt(sigmasq_alpha);
sigma_beta = sqrt(sigmasq_beta);
}
model {
mu_alpha ~ normal(0, 100);
mu_beta ~ normal(0, 100);
sigmasq_y ~ inv_gamma(0.001, 0.001);
sigmasq_alpha ~ inv_gamma(0.001, 0.001);
sigmasq_beta ~ inv_gamma(0.001, 0.001);
alpha ~ normal(mu_alpha, sigma_alpha); // vectorized
beta ~ normal(mu_beta, sigma_beta); // vectorized
for (n in 1:N)
for (t in 1:T)
y[n,t] ~ normal(alpha[n] + beta[n] * (x[t] - xbar), sigma_y);
}
generated quantities {
real alpha0;
alpha0 = mu_alpha - xbar * mu_beta;
}
```
我还有以下数据
```{r}
df <- read_delim("https://raw.githubusercontent.com/wiki/stan-dev/rstan/rats.txt",delim = " ")
y <- as.matrix(df)
x <- c(8,15,22,29,36)
xbar <- mean(x)
N <- nrow(y)
T <- ncol(y)
```
github 上的 documentation 显示 rats_fit <- stan(file = 'https://raw.githubusercontent.com/stan-dev/example-models/master/bugs_examples/vol1/rats/rats.stan')
,但由于我使用的是块,所以我没有可引用的文件。
我试过 stan(rats)
、summary(rats)
、print(rats)
,但其中 none 似乎有效。
第一个 RMarkdown 块在幕后调用 rats <- rstan::stan_model(model_code=the_text)
,因此为了从该后验分布中采样,您最终需要执行 rats_fit <- sampling(rats, data = list())
,其其余参数与 stan
。但是你必须先调用 library(rstan)
。
谢谢!在您的帮助下,我得出了以下结论
library(tidyverse)
df <- read_delim("https://raw.githubusercontent.com/wiki/stan-dev/rstan/rats.txt",delim = " ")
y <- as.matrix(df)
x <- c(8,15,22,29,36)
xbar <- mean(x)
N <- nrow(y)
T <- ncol(y)
library(rstan)
options(mc.cores = parallel::detectCores())
rats_fit <- rstan::sampling(rats,
data = list(y,x,xbar,N,T))
rstan::summary(rats_fit)
我想在 rnotebook
中使用rstan
获得模型的估计系数
我有以下 stan
块:
```{stan output.var="rats"}
data {
int<lower=0> N;
int<lower=0> T;
real x[T];
real y[N,T];
real xbar;
}
parameters {
real alpha[N];
real beta[N];
real mu_alpha;
real mu_beta; // beta.c in original bugs model
real<lower=0> sigmasq_y;
real<lower=0> sigmasq_alpha;
real<lower=0> sigmasq_beta;
}
transformed parameters {
real<lower=0> sigma_y; // sigma in original bugs model
real<lower=0> sigma_alpha;
real<lower=0> sigma_beta;
sigma_y = sqrt(sigmasq_y);
sigma_alpha = sqrt(sigmasq_alpha);
sigma_beta = sqrt(sigmasq_beta);
}
model {
mu_alpha ~ normal(0, 100);
mu_beta ~ normal(0, 100);
sigmasq_y ~ inv_gamma(0.001, 0.001);
sigmasq_alpha ~ inv_gamma(0.001, 0.001);
sigmasq_beta ~ inv_gamma(0.001, 0.001);
alpha ~ normal(mu_alpha, sigma_alpha); // vectorized
beta ~ normal(mu_beta, sigma_beta); // vectorized
for (n in 1:N)
for (t in 1:T)
y[n,t] ~ normal(alpha[n] + beta[n] * (x[t] - xbar), sigma_y);
}
generated quantities {
real alpha0;
alpha0 = mu_alpha - xbar * mu_beta;
}
```
我还有以下数据
```{r}
df <- read_delim("https://raw.githubusercontent.com/wiki/stan-dev/rstan/rats.txt",delim = " ")
y <- as.matrix(df)
x <- c(8,15,22,29,36)
xbar <- mean(x)
N <- nrow(y)
T <- ncol(y)
```
github 上的 documentation 显示 rats_fit <- stan(file = 'https://raw.githubusercontent.com/stan-dev/example-models/master/bugs_examples/vol1/rats/rats.stan')
,但由于我使用的是块,所以我没有可引用的文件。
我试过 stan(rats)
、summary(rats)
、print(rats)
,但其中 none 似乎有效。
第一个 RMarkdown 块在幕后调用 rats <- rstan::stan_model(model_code=the_text)
,因此为了从该后验分布中采样,您最终需要执行 rats_fit <- sampling(rats, data = list())
,其其余参数与 stan
。但是你必须先调用 library(rstan)
。
谢谢!在您的帮助下,我得出了以下结论
library(tidyverse)
df <- read_delim("https://raw.githubusercontent.com/wiki/stan-dev/rstan/rats.txt",delim = " ")
y <- as.matrix(df)
x <- c(8,15,22,29,36)
xbar <- mean(x)
N <- nrow(y)
T <- ncol(y)
library(rstan)
options(mc.cores = parallel::detectCores())
rats_fit <- rstan::sampling(rats,
data = list(y,x,xbar,N,T))
rstan::summary(rats_fit)