Модели смесей в Stan - векторизация

Сейчас я изучаю Стэна и хотел реализовать простую модель смеси.

В справочнике (stan-reference-2.14.0) уже есть решение:

data {
  int<lower=1> K; // number of mixture components
  int<lower=1> N; // number of data points
  real y[N]; // observations
}
parameters {
  simplex[K] theta; // mixing proportions
  real mu[K]; // locations of mixture components
  real<lower=0> sigma[K]; // scales of mixture components
}
model {
  real ps[K]; // temp for log component densities
  sigma ~ cauchy(0, 2.5);
  mu ~ normal(0, 10);
  for (n in 1:N) {
    for (k in 1:K) {
      ps[k] = log(theta[k])
      + normal_lpdf(y[n] | mu[k], sigma[k]);
    }
    target += log_sum_exp(ps);
  }
}

На следующей странице описывается, что векторизация внешнего цикла невозможна. Однако мне было интересно, сохраняется ли параллелизация внутреннего цикла.

И вот я попробовал следующую модель:

data {
  int<lower=1> K; // number of mixture components
  int<lower=1> N; // number of data points
  real y[N]; // observations
}

parameters {
  simplex[K] theta; // mixing proportions
  vector[K] mu; // locations of mixture components
  vector<lower=0>[K] sigma; // scales of mixture components
}

model {
  vector[K] ps;//[K]; // temp for log component densities
  vector[K] ppt;
  sigma ~ cauchy(0, 2.5);
  mu ~ normal(0, 10);
  for (n in 1:N) {
    ppt = log(theta);
    /*
    for (k in 1:K) {
      ps[k] = ppt[k] + //log(theta[k])
        normal_lpdf(y[n] | mu[k], sigma[k]);
    }
    */
    ps = ppt + normal_lpdf(y[n] | mu, sigma);
    target += log_sum_exp(ps);
  }
}

... и эта модель дает неправильные оценки (в отличие от исходной модели).

data("faithful")
erupdata <- list(
  K = 2,
  N = length(faithful$eruptions),
  y = faithful$eruptions
)

fiteruptions <- stan(file = 'mixturemodel.stan', data = erupdata, iter = 1000, chains = 1)

Мне интересно, что я неправильно понимаю в спецификации модели. Я хотел бы понять разницу, которую обеспечивает синтаксис (среди прочего, разницу между vector[K] и real[K]), и, возможно, получить более глубокое понимание Стэна.


rstan stan mixture-model mixture
person Drey    schedule 03.02.2017    source источник

Ответы (2)


arrow_upward
2
arrow_downward

Вторая программа определяет другую плотность. normal_lpdf возвращает одно скалярное значение, которое представляет собой сумму PDF-файлов журнала по контейнерам данных / параметров.

В руководстве есть глава о матрицах / векторах и массивах.

Вы хотите поднять определение ppt выше для повышения эффективности.

person Bob Carpenter    schedule 04.02.2017

arrow_upward
0
arrow_downward

Я не являюсь пользователем Stan и не объясняю разницу между vector[K] и real[K]. Но согласно некоторым ресурсам (здесь и здесь), предыдущие распределения места смешения должны быть определены как ordered[K] вместо vector[K], чтобы избежать проблемы неидентифицируемости / заменяемости априорных значений.

Для более глубокого понимания вышеупомянутой проблемы, Идентификация моделей байесовской смеси настоятельно рекомендуется. Рабочая модель скопирована из этого блога:

data {
 int<lower = 0> N;
 vector[N] y;
}

parameters {
  ordered[2] mu;
  real<lower=0> sigma[2];
  real<lower=0, upper=1> theta;
}

model {
 sigma ~ normal(0, 2);
 mu ~ normal(0, 2);
 theta ~ beta(5, 5);
 for (n in 1:N)
   target += log_mix(theta,
                     normal_lpdf(y[n] | mu[1], sigma[1]),
                     normal_lpdf(y[n] | mu[2], sigma[2]));
}
person Minh Vu    schedule 13.12.2018