Как сделать выборку из суммы двух распределений: биномиального и пуассоновского - PullRequest
2 голосов
/ 01 апреля 2020

Есть ли способ предсказать значение из суммы двух распределений? Я получаю синтаксическую ошибку в rstan, когда пытаюсь оценить y здесь: y ~ binomial (,) + poisson ()


library(rstan)

BH_model_block <- "
data{
  int y; 
  int a; 
}

parameters{
  real <lower = 0, upper = 1> c;
  real <lower = 0, upper = 1> b;
}

model{
  y ~ binomial(a,b)+ poisson(c);
}
"
BH_model <- stan_model(model_code = BH_model_block)
BH_fit <- sampling(BH_model,
                   data = list(y = 5,
                               a = 2), 
                   iter= 1000)

Выдает эту ошибку:

SYNTAX ERROR, MESSAGE(S) FROM PARSER:

  error in 'model2c6022623d56_457bd7ab767c318c1db686d1edf0b8f6' at line 13, column 20
  -------------------------------------------------
    11: 
    12: model{
    13:   y ~ binomial(a,b)+ poisson(c);
                           ^
    14: }
  -------------------------------------------------

PARSER EXPECTED: ";"
Error in stanc(file = file, model_code = model_code, model_name = model_name,  : 
  failed to parse Stan model '457bd7ab767c318c1db686d1edf0b8f6' due to the above error.

Ответы [ 2 ]

2 голосов
/ 03 апреля 2020

Stan не поддерживает целочисленные параметры, поэтому технически вы не можете этого сделать. Для двух реальных переменных это будет выглядеть так:

parameters {
  real x;
  real y;
}
transformed parameters {
  real z = x + y;
}
model {
  x ~ normal(0, 1);
  y ~ gamma(0.1, 2);
}

Тогда вы получите распределение суммы для z. Если переменные являются дискретными, они не будут компилироваться.

Если вам не нужно z в модели, то вы можете сделать это в блоке сгенерированных количеств,

generated quantities {
  int x = binomial_rng(a, b);
  int y = poisson_rng(c);
  int z = x + y;
}

Недостатком этого является то, что ни одна из переменных не доступна в блоке модели. Если вам нужны дискретные параметры, они должны быть обособлены, как описано в главе руководства пользователя по скрытым дискретным параметрам (также в главе о смесях и НММ). С Пуассоном это не так просто, потому что поддержка не ограничена. Если ожидания двух дискретных распределений невелики, вы можете сделать это примерно с помощью al oop над вероятными значениями.

Из примера в исходном посте видно, что z - это данные. Это немного другая маргинализация по сравнению с x и y, но вы суммируете только по x и y, так что x + y = z, так что комбинаторика значительно уменьшается.

1 голос
/ 07 апреля 2020

Альтернативой является замена Бинома Пуассоном и использование аддитивности Пуассона:

BH_model_block <- "
data{
  int y; 
  int a; 
}

parameters{
  real <lower = 0, upper = 1> c;
  real <lower = 0, upper = 1> b;
}

model{
  y ~ poisson(a * b + c);
}
"

Это отличается тем, что если b не мало, Бином имеет меньшую дисперсию, чем Пуассон, но, может быть, в любом случае есть чрезмерная дисперсия?

...