Написание собственной функции Probit в Stan - PullRequest
Новая
       15

Написание собственной функции Probit в Stan

1 голос
/ 16 мая 2019

Я пытаюсь кодировать пользовательскую функцию Probit в Stan, чтобы улучшить мое понимание языка 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);
    }
}

Error

ПАРСЕР ОЖИДАЕТСЯ: Ошибка в stanc (model_code = paste (program, collapse = "\ n"), имя_модели = имя_модели_папки: не удалось проанализировать модель 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 )

1 Ответ

2 голосов
/ 16 мая 2019

Полное сообщение об ошибке 

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.

Но написанные вами функции в любом случае не реализуют логарифмическую вероятность пробитной модели. Во-первых, стандартное отклонение ошибок не идентифицируется данными, поэтому вам не нужно 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));
...