AIC, BIC значения ARIMA с ограниченными коэффициентами в R

Различные способы задания одной и той же модели AR (или MA), оцениваемой функцией arima() в пакете forecast в R, дают различные значения BIC (байесовского информационного критерия).

Почему это происходит?

Рассмотрим две модели:

(1) AR (1)
(2) AR(2) с коэффициентом на AR2, ограниченным нулем

На бумаге эти две модели одинаковы. Однако их оценка может отличаться (?). Не знаю, почему они дают одинаковые оценки коэффициентов, равные логарифмические значения правдоподобия и равные значения AIC - но разные значения BIC.

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

Это оправданно, или это ошибка ???

Интересно, в чем разница и как вычисляется BIC в случае (2). Я хотел бы иметь возможность ... воспроизвести результаты, так что мне нужно понять, как здесь все работает.

Ниже я привожу воспроизводимый пример. Выполнив его в R, посмотрите на напечатанные значения BIC, а также AICc-они отличаются между моделями. 
library(forecast)
T=1000; seed=1; set.seed(seed); x=rnorm(T) 
model1=arima(x,order=c(1,0,0)                 ,method="CSS-ML",transform.pars=FALSE)
model2=arima(x,order=c(2,0,0),fixed=c(NA,0,NA),method="CSS-ML",transform.pars=FALSE)
print(model1)
print(model2)
То же самое относится к моделям AR(p) и MA(q), которые я не обсуждаю явно, чтобы упростить их. Было бы здорово, если бы кто-нибудь объяснил, почему это происходит. Спасибо!
1 2
time-series

1 ответ:

Вычисление AICc и BIC выполняется в функции forecast:::print.Arima, в то время как AIC возвращается функцией arima(). Если вы посмотрите на код для forecast:::print.Arima, то увидите следующее: 

npar <- length(x$coef) + 1
nstar <- length(x$residuals) - x$arma[6] - x$arma[7] * x$arma[5]
bic <- x$aic + npar * (log(nstar) - 2)
aicc <- x$aic + 2 * npar * (nstar/(nstar - npar - 1) - 1)
Обратите внимание, что npar не учитывает не оцененные коэффициенты (т. е. коэффициенты, ограниченные заданными значениями). Он предполагает, что все коэффициенты в x$coef были оценены. Это можно было бы исправить, используя 
npar <- length(x$coef[x$mask]) + 1

Я исправил версию пакета github , поэтому CRAN версия будет обновлена в следующем выпуске.

2