Как bnlearn рассчитывает BIC непрерывных данных?

Я использую пакет bnlearn в R и хотел бы знать, как пакет вычисляет BIC-g (BIC в распределении по Гауссу).

Давайте создадим структуру, и я смогу найти показатель BIC следующим образом.

library(bnlearn)
X = iris[, 1:3]
names(X) = c("A", "B", "C")
Network = empty.graph(names(X))
bnlearn::score(Network, X, type="bic-g")

bnlearn предоставляет мне более подробную информацию о том, как можно рассчитать этот балл,

bnlearn::score(Network, X, type="bic-g", debug=TRUE)

И это приводит к

----------------------------------------------------------------
* processing node A.
  > loglikelihood is -184.041441.
  > penalty is 2.505318 x 2 = 5.010635.
----------------------------------------------------------------
* processing node B.
  > loglikelihood is -87.777815.
  > penalty is 2.505318 x 2 = 5.010635.
----------------------------------------------------------------
* processing node C.
  > loglikelihood is -297.588727.
  > penalty is 2.505318 x 2 = 5.010635.
[1] -584.4399

Я знаю, как рассчитать BIC для дискретных данных в байесовских сетях, ссылаясь на здесь. Но я не знаю, как это можно обобщить на совместный гауссовский (многомерный нормальный) случай.

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

bnlearn::score(Network, X, type="loglik-g", debug=TRUE)

Но я хотел бы знать, как я могу конкретно рассчитать вероятности и штрафы, учитывая данные.

Я нашел материал, объясняющий Laplace Approximation (см. стр. 57), но я не мог связать это.

Кто-нибудь, чтобы помочь мне?


r bnlearn
person moreblue    schedule 26.02.2019    source источник
comment
логлик рассчитывается с sum(sapply(X, function(i) dnorm(i, mean(i), sd(i), log=TRUE))), а штрафной срок с 0.5*log(nrow(X)) * 6, который равен 0,5 * log количество наблюдений * количество параметров   -  person user20650    schedule 28.02.2019

Ответы (1)


arrow_upward
1
arrow_downward

БИК рассчитывается как

BIC = -2* logLik + nparams* log(nobs)

но в bnlearn это масштабируется на -2 (см.?score), чтобы дать

BIC = logLik -0,5* nparams* log(nobs)

Итак, для вашего примера без ребер вероятность рассчитывается с использованием предельных средних и ошибок (или, в более общем случае, для каждого узла количество параметров задается путем суммирования 1 (перехват) + 1 (остаточная ошибка) + количество родителей ), например

library(bnlearn)
X = iris[, 1:3]
names(X) = c("A", "B", "C")
Network = empty.graph(names(X))

(ll = sum(sapply(X, function(i) dnorm(i, mean(i), sd(i), log=TRUE)))) 
#[1] -569.408
(penalty = 0.5* log(nrow(X))* 6)
#[1] 15.03191

ll - penalty
#[1] -584.4399

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

Network = set.arc(Network, "A", "B")

Нам нужны компоненты логарифмического правдоподобия из узла A и C

(llA = with(X, sum(dnorm(A, mean(A), sd(A), log=TRUE))))
#[1] -184.0414
(llC = with(X, sum(dnorm(C, mean(C), sd(C), log=TRUE))))
#[1] -297.5887

и мы получаем условные вероятности B из линейной регрессии

m = lm(B ~ A, X)
(llB = with(X, sum(dnorm(B, fitted(m), stats::sigma(m), log=TRUE))))
#[1] -86.73894

Предоставление

(ll = llA + llB + llC)
#[1] -568.3691
(penalty = 0.5* log(nrow(X))* 7)
#[1] 17.53722
ll - penalty
#[1] -585.9063 

#  bnlearn::score(Network, X, type="bic-g", debug=TRUE)
# ----------------------------------------------------------------
# * processing node A.
#    loglikelihood is -184.041441.
#    penalty is 2.505318 x 2 = 5.010635.
# ----------------------------------------------------------------
# * processing node B.
#    loglikelihood is -86.738936.
#    penalty is 2.505318 x 3 = 7.515953.
# ----------------------------------------------------------------
# * processing node C.
#    loglikelihood is -297.588727.
#    penalty is 2.505318 x 2 = 5.010635.
# [1] -585.9063
person user20650    schedule 01.03.2019