在 stan 中使用数据增强进行 Probit 回归
Probit regression with data augmentation in stan
我正在尝试使用 stan 做一个带有数据扩充的概率模型。这是我们得到结果 y
0/1 的地方,它告诉我们潜在变量 ystar
的符号。这是我目前所知道的,但我不确定如何在 model
部分添加有关 y
的信息。有什么想法吗?
data {
int<lower=0> N; // number of obs
int<lower=0> K; // number of predictors
int<lower=0,upper=1> y[N]; // outcomes
matrix[N, K] x; // predictor variables
}
parameters {
vector[K] beta; // beta coefficients
vector[N] ystar; // latent variable
}
model {
vector[N] mu;
beta ~ normal(0, 100);
mu <- x*beta;
ystar ~ normal(mu, 1);
}
你可以
data {
int<lower=0> N; // number of obs
int<lower=0> K; // number of predictors
vector<lower=-1,upper=1> sign; // y = 0 -> -1, y = 1 -> 1
matrix[N, K] x; // predictor variables
}
parameters {
vector[K] beta; // beta coefficients
vector<lower=0>[N] abs_ystar; // latent variable
}
model {
beta ~ normal(0, 100);
// ignore the warning about a Jacobian from the parser
sign .* abs_ystar ~ normal(x * beta, 1);
}
就是说,没有理由在 Stan 中为二元概率模型进行数据扩充,除非某些结果缺失或其他原因。这样做更直接(并将参数 space 减少为 K 而不是 K + N)
data {
int<lower=0> N; // number of obs
int<lower=0> K; // number of predictors
int<lower=0,upper=1> y[N]; // outcomes
matrix[N, K] x; // predictor variables
}
parameters {
vector[K] beta; // beta coefficients
}
model {
vector[N] mu;
beta ~ normal(0, 100);
mu <- x*beta;
for (n in 1:N) mu[n] <- Phi(mu[n]);
y ~ bernoulli(mu);
}
如果你真的关心潜在效用,你可以通过 generated quantities
块中的拒绝采样来生成它,就像这样
generated quantities {
vector[N] ystar;
{
vector[N] mu;
mu <- x * beta;
for (n in 1:N) {
real draw;
draw <- not_a_number();
if (sign[n] == 1) while(!(draw > 0)) draw <- normal_rng(mu[n], 1);
else while(!(draw < 0)) draw <- normal_rng(mu[n], 1);
ystar[n] <- draw;
}
}
}
我正在尝试使用 stan 做一个带有数据扩充的概率模型。这是我们得到结果 y
0/1 的地方,它告诉我们潜在变量 ystar
的符号。这是我目前所知道的,但我不确定如何在 model
部分添加有关 y
的信息。有什么想法吗?
data {
int<lower=0> N; // number of obs
int<lower=0> K; // number of predictors
int<lower=0,upper=1> y[N]; // outcomes
matrix[N, K] x; // predictor variables
}
parameters {
vector[K] beta; // beta coefficients
vector[N] ystar; // latent variable
}
model {
vector[N] mu;
beta ~ normal(0, 100);
mu <- x*beta;
ystar ~ normal(mu, 1);
}
你可以
data {
int<lower=0> N; // number of obs
int<lower=0> K; // number of predictors
vector<lower=-1,upper=1> sign; // y = 0 -> -1, y = 1 -> 1
matrix[N, K] x; // predictor variables
}
parameters {
vector[K] beta; // beta coefficients
vector<lower=0>[N] abs_ystar; // latent variable
}
model {
beta ~ normal(0, 100);
// ignore the warning about a Jacobian from the parser
sign .* abs_ystar ~ normal(x * beta, 1);
}
就是说,没有理由在 Stan 中为二元概率模型进行数据扩充,除非某些结果缺失或其他原因。这样做更直接(并将参数 space 减少为 K 而不是 K + N)
data {
int<lower=0> N; // number of obs
int<lower=0> K; // number of predictors
int<lower=0,upper=1> y[N]; // outcomes
matrix[N, K] x; // predictor variables
}
parameters {
vector[K] beta; // beta coefficients
}
model {
vector[N] mu;
beta ~ normal(0, 100);
mu <- x*beta;
for (n in 1:N) mu[n] <- Phi(mu[n]);
y ~ bernoulli(mu);
}
如果你真的关心潜在效用,你可以通过 generated quantities
块中的拒绝采样来生成它,就像这样
generated quantities {
vector[N] ystar;
{
vector[N] mu;
mu <- x * beta;
for (n in 1:N) {
real draw;
draw <- not_a_number();
if (sign[n] == 1) while(!(draw > 0)) draw <- normal_rng(mu[n], 1);
else while(!(draw < 0)) draw <- normal_rng(mu[n], 1);
ystar[n] <- draw;
}
}
}