Статистика бутстрап-теста R-структуры из-за несферических возмущений

(1) Асимметрия и тяжелые хвосты. Как обычно в моделях линейной регрессии, асимптотическое обоснование вывода не зависит от нормальности, а также выполняется для любого другого распределения ошибок при нулевом математическом ожидании, гомоскедастичности и отсутствии корреляция (обычные предположения Гаусса-Маркова). Однако, если у вас есть подходящее асимметричное распределение для интересующих вас данных, вы можете повысить эффективность, основываясь на соответствующей модели. Например, пакет glogis предоставляет некоторые функции для тестирования структурных изменений и датирования на основе обобщенного логистического распределения, допускающего тяжелые хвосты и асимметрию. Windberger & Zeileis (2014, Eastern European Economics, 52, 66–88, doi:10.2753/EEE0012-8775520304) использовали его для отслеживания изменений асимметрии динамики инфляции с течением времени. (См. рабочий пример в разделе ?breakpoints.glogisfit.) Кроме того, если асимметрия сама по себе не представляет особого интереса, то преобразование журнала или sqrt также может быть достаточно хорошим, чтобы сделать данные более «нормальными».

(2) Гетероскедастичность и автокорреляция. Как обычно в моделях линейной регрессии, стандартные ошибки (или, в более широком смысле, ковариационная матрица) не согласуются при наличии гетероскедастичности и/или автокорреляции. Можно либо попытаться включить это явно в модель (например, модель AR), либо рассматривать его как неприятный термин и использовать ковариационные матрицы гетероскедастичности и автокорреляции (HAC) (например, Newey-West или квадратичный спектральный kernal Эндрюса HAC) . Функция Fstats() в strucchange позволяет подключать такие оценщики, например, из пакета sandwich. См. ?durab пример использования vcovHC().

(3) P-значения начальной загрузки и перестановки. «Скремблирование» временных рядов, описанное выше, больше похоже на применение перестановок (т. е. выборки без замены), а не на загрузку (т. е. выборку с заменой) . Первое возможно, если ошибки некоррелированы или взаимозаменяемы. Если вы регрессируете только на константу, вы можете использовать функцию maxstat_test() из пакета coin для выполнения теста supF. Статистика теста вычисляется несколько другим способом, однако можно показать, что он эквивалентен тесту supF в случае только констант (см. Zeileis & Hothorn, 2013, Statistical Papers, 54, 931–954, doi:10.1007/s00362-013-0503-4 ). Если вы хотите выполнить тест перестановки в более общей модели, вам придется выполнять перестановки «вручную» и просто сохранять статистику теста для каждой перестановки. В качестве альтернативы бутстрап можно применить, например, через пакет boot (где вам все равно нужно будет написать свою собственную небольшую функцию, которая вычисляет статистику теста из заданного бутстрап-образца). Есть также некоторые пакеты R (например, tseries), которые реализуют схемы начальной загрузки для зависимых серий.

Я добавляю второй ответ для анализа предоставленных смоделированных данных Thames. Что касается пунктов из моего первого общего методологического ответа: (1) В этом случае преобразование log() явно подходит для работы с крайней асимметрией наблюдений. (2) Поскольку данные гетероскедастичны, вывод должен быть основан на ковариациях HC или HAC. Ниже я использую оценку Newey-West HAC, хотя данные просто гетероскедастичны, но не автокоррелированы. Вывод с поправкой на HAC влияет на тест supF и доверительные интервалы для оценок контрольных точек. Сами точки останова и соответствующие специфические для сегмента перехваты оцениваются с помощью OLS, т. е. с учетом гетероскедастичности как мешающего члена. (3) Я не добавлял никакого бутстрепного или перестановочного вывода, поскольку в этом случае асимптотический вывод кажется достаточно убедительным.

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

library("strucchange")
set.seed(12)
Thames <- ts(c(rlnorm(120, 0, 1), rlnorm(120, 2, 2), rlnorm(120, 4, 1)),
  frequency = 12, start = c(1985, 1))

Затем мы вычисляем последовательность статистик Вальда/Ф с поправкой на HAC и оцениваем оптимальные контрольные точки (для m = 1, 2, 3,... разрывов) с помощью МНК. Чтобы проиллюстрировать, насколько лучше это работает для серии в журналах, а не в уровнях, показаны обе версии.

fs_lev <- Fstats(Thames ~ 1, vcov = NeweyWest)
fs_log <- Fstats(log(Thames) ~ 1, vcov = NeweyWest)
bp_lev <- breakpoints(Thames ~ 1) 
bp_log <- breakpoints(log(Thames) ~ 1) 

Визуализация ниже показывает временной ряд с подобранными точками пересечения в первой строке, последовательность статистики Вальда/Ф с 5%-ным критическим значением теста supF во второй строке и остаточную сумму квадратов и BIC для выбора количество точек останова в последней строке. Код для воспроизведения изображения находится в конце этого ответа.

Оба теста supF явно значимы, но на уровнях (sctest(fs_lev)) статистика теста составляет «всего» 82,79, а в логах (sctest(fs_log)) — 282,46. Кроме того, два пика, относящиеся к двум контрольным точкам, гораздо лучше видны при анализе данных в журналах.

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

confint(bp_lev, breaks = 2, vcov = NeweyWest)
## 
##          Confidence intervals for breakpoints
##          of optimal 3-segment partition: 
## 
## Call:
## confint.breakpointsfull(object = bp_lev, breaks = 2, vcov. = NeweyWest)
## 
## Breakpoints at observation number:
##   2.5 % breakpoints 97.5 %
## 1    NA         125     NA
## 2   202         242    263

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

confint(bp_log, breaks = 2, vcov = NeweyWest)
## 
##          Confidence intervals for breakpoints
##          of optimal 3-segment partition: 
## 
## Call:
## confint.breakpointsfull(object = bp_log, breaks = 2, vcov. = NeweyWest)
## 
## Breakpoints at observation number:
##   2.5 % breakpoints 97.5 %
## 1   107         119    121
## 2   238         240    250
## 
## Corresponding to breakdates:
##   2.5 %    breakpoints 97.5 %  
## 1 1993(11) 1994(11)    1995(1) 
## 2 2004(10) 2004(12)    2005(10)

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

par(mfrow = c(3, 2))
plot(Thames, main = "Thames")
lines(fitted(bp_lev, breaks = 2), col = 4, lwd = 2)
plot(log(Thames), main = "log(Thames)")
lines(fitted(bp_log, breaks = 2), col = 4, lwd = 2)
lines(confint(bp_log, breaks = 2, vcov = NeweyWest))
plot(fs_lev, main = "supF test")
plot(fs_log, main = "supF test")
plot(bp_lev)
plot(bp_log)