在 Stan 中编写自定义 Probit 函数
Writing a custom Probit function in Stan
我正在尝试在 Stan 中编写自定义 Probit 函数,以提高我对 Stan 语言和可能性的理解。到目前为止,我已经编写了普通 pdf 的对数,但收到一条错误消息,当我尝试编写可能性时,我发现它无法理解。我做错了什么?
斯坦模型
functions {
real normal_lpdf(real mu, real sigma) {
return -log(2 * pi()) / 2 - log(sigma)
- square(mu) / (2 * sigma^2);
}
real myprobit_lpdf(int y | real mu, real sigma) {
return normal_lpdf(mu, sigma)^y * (1 - normal_lpdf(mu, sigma))^(1-y);
}
}
data {
int N;
int y[N];
}
parameters {
real mu;
real<lower = 0> sigma;
}
model {
for (n in 1:N) {
target += myprobit_lpdf(y[n] | mu, sigma);
}
}
错误
预期解析器:
stanc(model_code = paste(program, collapse = "\n"), model_name = model_cppname, 错误:
由于上述错误,无法解析 Stan 模型 'Probit_lpdf'。
R代码模拟数据
## DESCRIPTION
# testing a Probit model
## DATA
N <- 2000
sigma <- 1
mu <- 0.3
u <- rnorm(N, 0, 2)
y.star <- rnorm(N, mu, sigma)
y <- ifelse(y.star > 0,1, 0)
data = list(
N = N,
y = y
)
## MODEL
out.stan <- stan("Probit_lpdf.stan",data = data, chains = 2, iter = 1000 )
完整的错误信息是
SYNTAX ERROR, MESSAGE(S) FROM PARSER:
Probabilty functions with suffixes _lpdf, _lpmf, _lcdf, and _lccdf,
require a vertical bar (|) between the first two arguments.
error in 'model2a7252aef8cf_probit' at line 7, column 27
-------------------------------------------------
5: }
6: real myprobit_lpdf(real y, real mu, real sigma) {
7: return normal_lpdf(mu, sigma)^y * (1 - normal_lpdf(mu, sigma))^(1-y);
^
8: }
-------------------------------------------------
这告诉您 normal_lpdf
函数除了三个输入和一个分隔第一个和第二个的竖线。
给您的函数取一个与 Stan 语言中已有的函数相同的名称也不是一个好主意,例如 normal_lpdf
.
但是您编写的函数并没有实现概率模型的 log-likelihood。首先,误差的标准差没有被数据识别出来,所以你不需要sigma
。那么,正确的表达式应该是
real Phi_mu = Phi(mu);
real log_Phi_mu = log(Phi_mu);
real log1m_Phi_mu = log1m(Phi_mu);
for (n in 1:N)
target += y[n] == 1 ? log_Phi_mu : log1m_Phi_mu;
尽管这只是一种缓慢的方式
target += bernoulli_lpmf(y | Phi(mu));
我正在尝试在 Stan 中编写自定义 Probit 函数,以提高我对 Stan 语言和可能性的理解。到目前为止,我已经编写了普通 pdf 的对数,但收到一条错误消息,当我尝试编写可能性时,我发现它无法理解。我做错了什么?
斯坦模型
functions {
real normal_lpdf(real mu, real sigma) {
return -log(2 * pi()) / 2 - log(sigma)
- square(mu) / (2 * sigma^2);
}
real myprobit_lpdf(int y | real mu, real sigma) {
return normal_lpdf(mu, sigma)^y * (1 - normal_lpdf(mu, sigma))^(1-y);
}
}
data {
int N;
int y[N];
}
parameters {
real mu;
real<lower = 0> sigma;
}
model {
for (n in 1:N) {
target += myprobit_lpdf(y[n] | mu, sigma);
}
}
错误
预期解析器: stanc(model_code = paste(program, collapse = "\n"), model_name = model_cppname, 错误: 由于上述错误,无法解析 Stan 模型 'Probit_lpdf'。
R代码模拟数据
## DESCRIPTION
# testing a Probit model
## DATA
N <- 2000
sigma <- 1
mu <- 0.3
u <- rnorm(N, 0, 2)
y.star <- rnorm(N, mu, sigma)
y <- ifelse(y.star > 0,1, 0)
data = list(
N = N,
y = y
)
## MODEL
out.stan <- stan("Probit_lpdf.stan",data = data, chains = 2, iter = 1000 )
完整的错误信息是
SYNTAX ERROR, MESSAGE(S) FROM PARSER:
Probabilty functions with suffixes _lpdf, _lpmf, _lcdf, and _lccdf,
require a vertical bar (|) between the first two arguments.
error in 'model2a7252aef8cf_probit' at line 7, column 27
-------------------------------------------------
5: }
6: real myprobit_lpdf(real y, real mu, real sigma) {
7: return normal_lpdf(mu, sigma)^y * (1 - normal_lpdf(mu, sigma))^(1-y);
^
8: }
-------------------------------------------------
这告诉您 normal_lpdf
函数除了三个输入和一个分隔第一个和第二个的竖线。
给您的函数取一个与 Stan 语言中已有的函数相同的名称也不是一个好主意,例如 normal_lpdf
.
但是您编写的函数并没有实现概率模型的 log-likelihood。首先,误差的标准差没有被数据识别出来,所以你不需要sigma
。那么,正确的表达式应该是
real Phi_mu = Phi(mu);
real log_Phi_mu = log(Phi_mu);
real log1m_Phi_mu = log1m(Phi_mu);
for (n in 1:N)
target += y[n] == 1 ? log_Phi_mu : log1m_Phi_mu;
尽管这只是一种缓慢的方式
target += bernoulli_lpmf(y | Phi(mu));