Расчет логарифмической вероятности набора наблюдений, отобранных из смеси двух нормальных распределений, с использованием R

Question

Я написал функцию для вычисления логарифмической вероятности набора наблюдений, отобранных из смеси двух нормальных распределений. Эта функция не дает мне правильного ответа.

Я не буду знать, из какого из двух распределений получен данный образец, поэтому функция должна суммировать по возможностям.

Эта функция принимает вектор пяти параметров модели в качестве первого аргумента (μ1, σ1, μ2, σ2 и p), где μi и σi - среднее значение и стандартное отклонение i-го распределения, а p - вероятность того, что выборка получена из первого распределение. Для второго аргумента функция принимает вектор наблюдений.

Я написал следующую функцию:

mixloglik <- function(p, v) {
    sum(log(dnorm(v, p[1], p[2])*p[5] + dnorm(v,p[3],p[4]))*p[5])
}

Я могу создать тестовые данные, для которых я знаю, что решение должно быть ~ -854.6359:

set.seed(42)
v<- c(rnorm(100), rnorm(200, 8, 2))
p <- c(0, 1, 6, 2, 0.5)

Когда я проверяю эту функцию на данных испытаний, я не получаю правильное решение

> mixloglik(p, v)
[1] -356.7194

Я знаю, что решение должно быть ~ -854.6359. Где я ошибаюсь в своей функции?

Answer 1

Правильное выражение для логарифмического правдоподобия выглядит следующим образом:

mixloglik <- function(p, v) {
  sum(log(p[5]*dnorm(v, p[1], p[2]) + (1 - p[5])*dnorm(v, p[3], p[4])))
}

Теперь попробуйте это:

set.seed(42)
v<- c(rnorm(100), rnorm(200, 8, 2))
p <- c(0, 1, 6, 2, 0.5)
mixloglik(p, v)
#[1] -854.6359

В подобных случаях лучший способ устранить ошибку - это перезапустить, переписав выражение на бумаге и перекодировать его.