【发布时间】:2021-11-13 09:44:20
【问题描述】:
我想提取 Stan 拟合的预测值(在生成量块中)并将它们与实际观察结果进行比较,但我找不到简单的解决方案。下面是一个简单的逻辑回归模型是如何做到的:
library(rstan)
library(tidyverse)
library(boot)
rstan_options(auto_write = TRUE)
options(mc.cores = parallel::detectCores())
T <- 40
set.seed(123)
x <- sort(runif(T, 0, 10))
alpha <- 1
beta <- 0.2
logit_p <- alpha + beta * x
p <- inv.logit(logit_p)
y <- rbinom(T, 1, p)
model_code <- "
data {
int<lower=0> N;
vector[N] x;
int<lower=0,upper=1> y[N];
}
parameters {
real alpha;
real beta;
}
model {
y ~ bernoulli_logit(alpha + beta * x);
}
generated quantities {
vector[N] z;
for (n in 1:N)
z[n] = bernoulli_logit_rng(alpha + beta * x[n]);
}"
model_data <- list(
N = T,
x = x,
y = y
)
stan_run <- stan(
data = model_data,
model_code = model_code
)
posterior <- rstan::extract(stan_run)
df <- as.data.frame(posterior$z)
df <-df %>% summarise(across(everything(.),
~ ifelse(length(.[which(. == 1)]) > length(.[which(. == 0)]), 1, 0)))
我不知道我的方法是否正确。有谁知道任何直接的方法吗?
【问题讨论】:
-
我在分配
model_data对象时收到错误“对象x未找到”。x和y是什么? -
@scrameri 我只是想展示我在拟合后提取后验的方法,因此 x 和 y 是什么无关紧要,您可以模拟您自己选择的任何 x 和 y 并将它们提供给模型。
-
好的,但即使在模拟 x 和 y 之后我也无法拟合模型(语法错误:分配中的尺寸不匹配,第 18 行)。您的问题似乎更多是关于如何从 data.frame 中提取和总结某些东西,所以也许最好提供一个示例
posterior对象,或者提供一个使用 rstan 和模拟 x 和 y 的可工作的、可重现的示例:stackoverflow.com/questions/5963269/…跨度> -
@scrameri 我假设任何知道答案的人甚至根本不需要运行代码(我给出了多余的代码),因为这个问题非常笼统,与模型或数据(可以是任何数据或模型)。但是现在你问我编辑了代码,现在它应该可以工作了。
-
是的,但是需要知道
posterior的结构是什么样的,这就是为什么现在有这个可重现的例子很好。不过,dputposterior对象的一小部分就足够了。