Краткое повторение ⏭

Предыдущие заметки серии:

В предыдущих заметках говорилось о важности своевременного определения наличия рецессии вообще и в текущий момент в частности. Предполагается, что текущий политический контекст будет способствовать затягиванию официального объявления рецессии в США. Не являясь политологом и не имея желания погружаться в детали противостояния республиканцев и демократов, я фокусирую свое внимание на экономитрической модели раннего или, правильнее сказать, своевременного обнаружения рецессии, публикуемой на сайте ФРС США. Такая модель призвана сигнализировать о наступлении рецессии за долго до официального объявления профильным ведомством The National Bureau of Economic Research (NBER).

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

Теория 📚

В данном разделе будет более-менее подробно разобрана модель динамического изменения факторов на основе цепи Маркова или в английском варианте dynamic-factor markov-switching model. В предыдущей заметке уже говорилось о том, что складывается ощущение работы экономической системы как-будто в двух режимах. Поддавшись силе интуиции можно предположить, что существует некая скрытая переменная, которая характеризует переключение между двумя или более режимами. Таким образом, модель может быть представлена в виде комбинации составляющей скрытой переменной и стандартной видимой регрессионной составляющей в т.ч. обобщенной регрессии, например, логистической или авторегрессионной модели. Предполагаю, что видимая составляющая читателю хорошо знакома и не требует подробных пояснений.

Скрытая составляющая может быть представлена в виде Марковского процесса (Марковской цепи событий) т.е. принимается предположение, что каждое последующее состояние системы зависит только от текущего состояния, количество возможных состояний (режимов) – конечно, а вероятности перехода между состояниями задаются матрицей. Для двух режимов такую матрицу можно проиллюстрировать следующим образом:

Другими словами, задав вероятности [p_{11}] и [p_{22}] можно получить матрицу перехода для всей системы: $$ P = \begin{bmatrix}p_{11} & 1-p_{11}\\ 1-p_{22} & p_{22}\end{bmatrix} $$ Вероятность состояния будущего периода в момент времени [t] можно записать многократным матричным умножением (возведением в степень), если не подмешивать какие-либо дополнительные факторы: $$ S_{t} = S_{0}P^{t} $$ Переходя к смешанной модели можно записать следующее выражение для случая двух режимов [s_{t}], одного предиктора [x] и простой линейной регрессии: $$ y_{t} = a_{1}x + b_{1} + \epsilon_{t},\quad when\quad s_{t} = 1\\ y_{t} = a_{2}x + b_{2} + \epsilon_{t},\quad when\quad s_{t} = 2 $$ , где:

  • [y_{t}] – целевая переменная, которая предсказывается
  • [a] и [b] – коэффициенты модели регрессии для состояния [s_{t}] от предикторов
  • [\epsilon_{t}] – случайная ошибка модели

Другими словами, имеются две модели линейной регрессии, где одна модель сменяет другую и обратно. Ничего особо сложного, если не ударятся в обобщения и большую абстракцию математического описания. Пусть модель описывает данные достаточно хорошо т.е. ошибка модели будет иметь нормальное распределение вокруг нуля без признаков гомоскедастичности т.е. дисперсия ошибки модели не будет меняться в зависимости от значений предикторов. В силу того, что выше в качестве примера была выбрана модель линейной регрессии функция распределения ошибки может быть записана: $$ f(y_{t}, x_{t}, a, b, s_{t}) = \frac{1}{\sigma_{s} \sqrt{2\pi}}exp\Big[-\frac{y_{t}-a_{s}x_{t}-b_{s}}{2\sigma_{s}}\Big] $$ Для оценки матрицы [P] используется логистическая функция связи с параметром [\eta]: $$ p_{11}=\frac{1}{1+exp(-\eta_{1})},\quad p_{22}=\frac{1}{1+exp(-\eta_{2})},\quad p_{12}=1-p_{11},\quad p_{21}=1-p_{22} $$ Параметры [\eta] вычисляется методом максимального правдоподобия наравне с прочими коэффициентами модели [a_{s}, b_{s}, \sigma_{s}].

Далее, пользуясь формулой полной вероятности, можно записать комбинированную условную плотность вероятности для момента [t], где перемешивается модель линейной регрессии с моделью Марковских цепей. Условной такая вероятность – является потому как она зависит от предыдущего состояния системы.

$$ \hat f(y_{t})=S_{1,t-1}p_{11}f_{1t}+S_{1,t-1}p_{12}f_{2t}+S_{2,t-1}p_{21}f_{1t}+S_{2,t-1}p_{22}f_{2t} $$ В силу того, что ранее были определены только два режима и сумма состояний в каждый момент времени должна быть равна единице т.е. включается хотя бы один и только один режим из двух возможных: $$ \sum_{j=1}^{2}S_{jt}=1 $$ Теперь можно записать апостериорную вероятность состояния системы [S_{t}=1] в момент времени [t] корректированную с учетом предыдущего состояния [S_{t-1}=1]: $$ S_{1,t}= \frac{S_{1,t-1}p_{11}f_{1t}}{\hat f(y_{t})} $$ Аналогичное выражение можно записать для вероятности получить состояние [S_{t}=2] в момент времени [t], корректированную с учетом предыдущего состояния [S_{t-1}=1]: $$ S_{2,t}= \frac{S_{1,t-1}p_{12}f_{2t}}{\hat f(y_{t})} $$

Аналогично можно записать вероятности переходов для всей матрицы [P], которая для двух режимов имеет размерность – 2 x 2. Таким образом, получилось выразить каждое текущее состояние через предыдущее: $$ S_{jt} \to S_{i,t-1} $$

С практической точки зрения расчет вероятностей производится сразу для двух режимов, а в качестве прогноза используется уравнение линейной регрессии того режима, который имеет более высокую вероятность. Вероятность получения того или иного режима называется filtered probability или фильтрованная вероятность.

Как уже упоминалось, поиск коэффициентов модели осуществляется методом максимального правдоподобия, который ищет параметры [a_{s}, b_{s}, \sigma_{s}, \eta_{s}]. Параметр [\eta_{s}] пересчитывается в матрицу [P], что в совокупности c [a_{s}, b_{s}] характеризует полученную модель. Аналогичные выкладки можно более компактно записать с использованием матричных выражений, но, на мой взгляд, в этом случае теряется наглядность того, что из себя представляют коэффициенты модели. Тем не менее, реальные расчеты в R разумнее производить именно с использованием матричных выражений, что позволяет легко адаптировать модель под большее число режимов и предикторов.

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

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

Расчеты 🧮

Отличной новостью – является то, что нет необходимости самому писать алгоритм данной модели, а можно воспользоваться готовым вариантом, упакованным в пакет R MSwM. Буду использовать этот пакет наравне с теми, которые я использовал в предыдущей заметке:

library(MSwM) # пакет который умеет считать марковскую модель с переключением
library(data.table) # быстрый пакет для работы с табличными данными
library(collapse) # быстрый пакет для работы с табличными данными
library(stringi) # пакет для работы с текстовой информацией 
library(DT) # пакет создания интерактивных таблиц
library(echarts4r) # интерактивные графики
library(firatheme) # палитра, которую я использую для непрерывных индикаторов
library(emphatic) # пакет для раскрашивания матриц и табличек для наглядности 

# Сохраняем палитру в отдельную переменную
my_pal <- palette.colors(palette = "Tableau") %>% unname() 

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

# Загрузка подготовленных данных
indctrs0 <- fst::read_fst("data/indctrs.fst", as.data.table = TRUE)

preds <- c("PAYEMS", "INDPRO", "W875RX1", "CMRMTSPL")

rec_dt <- indctrs0  %>%  
  # Получение относительных изменений
  .[, stri_c(preds, "_ch") := lapply(.SD, \(x) shift((x - shift(x))/x)), .SDcols = preds]  %>% 
  # Получение сдвига на один период
  .[, stri_c(c(preds, stri_c(preds, "_ch")), "lg") := lapply(.SD, shift), .SDcols = c(preds, stri_c(preds, "_ch"))] |> 
  fsubset(date > as.Date("1970-01-01") & date < as.Date("2022-01-01")) |>
  na_omit()

datatable(rec_dt[1:30], style = 'bootstrap4', extensions = 'Responsive', 
          options = list(pageLength = 6), caption = "Индикаторы рецессии")

Спецификация модели регрессии с марковским переключением должна включать базовую модель. В данном случае будет использована модель логистической регрессии, которая описывалась в предыдущей заметке. Также необходимо указать количество режимов k = 2 и выбрать предикторы, которые имеют переключения sw, что делается с помощью TRUE/FALSE вектора в том же порядке в каком определены регрессоры (предикторы) базовой модели, включая свободный член.

Отдельно стоит затронуть параметр control, который позволяет сделать настройки для оптимизатора. Пакет MSwM используют собственный оптимизатор, расчитывающий матрицы Гессе (матрица вторых частных производных), поэтому параметры сложно назвать привычными и описания этих параметров в справке не сильно добавляет понимания. Я посмотрел исходный код по диагонали, но не смог найти интуитивного и простого пояснения. Могу сказать, что модель очень сильно зависит от выбранных начальных значений и для того чтобы получить хороший результаты нужно существенно увеличить перебор этих самых начальных значений, что собственно я и сделал. Такой подход существенно увеличил время вычислений, но позволил получить спецификацию модели, которая на выходе дает стабильно-воспроизводимый результат, качество которого можно оценивать через метрику AIC. В данном случае AIC ~ 140 и дальнейшее увеличение maxiter параметров не улучшает этот результат. Еще раз отмечу, что важна не столько характеристика качества модели AIC, но сколько воспроизводимость получения приблизительно одного и того же значения при многократных запусках расчетов.

# Формула логистической регрессии
frml <- formula(stri_c("USRECD ~", stri_c(stri_subset(names(rec_dt), regex = "ch"), collapse = "+")))

# Оценка коэффициентов логистической регрессии
rec_glm <- glm(data = rec_dt, frml, family="binomial")

# Оценка коэффициентов логистической регрессии с марковским переключением
rec_msm <- msmFit(rec_glm, k = 2, sw = rep(TRUE, 9), family = "binomial", 
                  control = list(tol = 1e-16, maxiterOuter = 200, # настройка оптимизатора
                                 maxiter = 1000, maxiterInner = 200))
# Таблица с коэффициентами модели
rec_msm@Coef |> 
  qDT(row.names.col = "regime") |>
  melt(id.vars = "regime") |> 
  dcast(variable ~ regime) |>
  frename(\(x)stri_c("regime_", x), cols = 2:3) |>
  datatable(style = 'bootstrap4', extensions = 'Responsive', 
            options = list(pageLength = 6), 
            caption = "Коэффициенты модели марковского переключения") |> 
  formatRound(columns = 2:3, digits = 2)

После нескольких минут расчета модель выдала результат в виде параметров логистической модели для двух режимов. Любопытно заметить, что знаки перед коэффициентами PAYEMS_ch и CMRMTSPL_ch при смене режима инвертируются т.е. предикторы в разных режимах предсказывают разные вероятности рецессии и это ровно то чего добивается модель: получить две существенно разные модели для двух режимов.

Марковская матрица вероятностей переходов между режимами получилась следующая:

TM <- round(rec_msm@transMat, 3)

hl_mat(TM, scale_color_fira(continuous = TRUE)) |> 
  as_html(style = "color:White;font-size:15px")
        [,1]  [,2]
[1,]   0.907 0.015
[2,]   0.093 0.985

Или в виде поясняющей диаграммы переходов с конкретными значениями вероятности:

Прогноз 🔮

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

В свете вышесказанного возникает вопрос: модель динамического изменения факторов на основе цепи Маркова создавалась прежде всего для прогнозирования или нет?

Для ответа на этот вопрос приведу цитату из работы3 Marcelle Chauvet

A dynamic factor model with regime switching is proposed as an empirical characterization of business cycles. The approach integrates the idea of comovements among macroeconomic variables and asymmetries of business cycle expansions and contractions. The first is captured with an unobservable dynamic factor and the second by allowing the factor to switch regimes. The model is estimated by maximizing its likelihood function and the empirical results indicate that the combination of these two features leads to a successful representation of the data relative to extant literature. This holds for within and out-of-sample and for both revised and real time data.

Или в переводе:

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

Под внешними данными как раз подразумевается способность осуществлять прогнозированием на те периоды времени когда сведений о наступлении рецессии еще нет. Значит логично будет использовать эту модель для прогнозирования наступления рецессии. Проблема заключается в том, что метода predict для пакета MSwM на языке R не существует и я не один, кто этим интересуется интересовался. Можно подумать, что это какое-то недоразумение и разработчики пакета R просто не настолько многочисленны чтобы делать все возможные опции. Хорошо, в этом случае можно посмотреть на аналогичную реализацию модели для самого популярного языка на планете – Python для пакета statmodels, но вот незадача: в этом пакете также нет метода predict, но есть соответствующие ветки обсуждения тут и тут.

Не хочется ударятся в конспирологию, но почему такая замечательная и важная модель, которую использует ФРС для ориентировки капитала всего мира относительно бизнес цикла, не имеет функции прогнозирования в самых популярных и очевидных языков для эконометрики? 🤯

На этой замечательной мысли можно было бы остановиться, но я не просто так распылялся: расписывал теорию и методику расчета коэффициентов модели. Кроме того, R – это прежде всего открытый софт, а следовательно можно найти исходный код модели и немного модифицировать его для получения самопального метода predict. Собственно именно в исходном коде я подсматривал как осуществляется оптимизация максимального правдоподобия. Теперь осталось найти метод для модели логистической регрессии, который называется .MSM.glm.msmFilter и заимствовать донорский код для своих манипуляций:

# Функция прогнозирования для модели msmFit (параметр object) и out-of-sample данных (параметр newdata)
pred_msm <- function(object, newdata){
  # object <- rec_msm # строчка для отладки модели
  
  model=object["model"] # непосредственно модель 
  k=object["k"] # количество режимов
  family=model$family # функция связи (логистическая функция)
  
  Coef=object["Coef"] # коэффициенты модели в виде матрицы (2 режима x 9 регрессеров)
  P=object["transMat"] # 2 х 2 матрица вероятностей переходов между режимами 
  
  # Некоторые преобразования, необходимые для расчетов
  newdata_prep <- get_vars(newdata, vars = "_ch", regex = TRUE) |>
    tfm(`(Intercept)`= 1) %>%  # свободный член в регрессии
    .[, .SD, .SDcols = names(Coef)] # сортировка предикторов в соответствии с очередностью в модели
  
  # Предварительные расчеты 
  # nr=length(model$model[,1]) # исходный код
  nr=fnrow(newdata_prep) # исправленный код
  
  terms=model.matrix(model)
  
  # Расчет условного среднего для двух моделей логистической регрессии
  #CondMean=family$linkinv(as.matrix(terms)%*%t(as.matrix(Coef))) # исходный код
  CondMean=family$linkinv(as.matrix(newdata_prep)%*%t(as.matrix(Coef))) # исправленный код
  # error= as.matrix(model$model[,1,drop=FALSE])%*%matrix(rep(1,k),nrow=1)-CondMean # исходный код
  
  # Расчет условной вероятности рецесии для двух моделей логистической регрессии
  # Likel=object@Likelihood(as.matrix(model$model[,1,drop=FALSE])%*%matrix(rep(1,k),nrow=1), mu=CondMean) # исходный код
  rec_prob <- dbinom(matrix(rep(1,k),nrow=1), size = 1, prob = CondMean) # исправленный код
  
  ###Расчет фильтрованной вероятности####
  fProb=matrix(data=0,nrow=nr,ncol=k)
  margLik=matrix(data=0,nrow=nr,ncol=1)
  
  # Расчет совокупной вероятности (нормализирующего коэффициента) для первого значения
  # margLik[1,1]=sum ((P%*%matrix(object["iniProb"],ncol=1)) * t(Likel[1,,drop=FALSE])) # исходный код
  margLik[1,1]=sum ((P%*%matrix(c(1,1),ncol=1))* t(rec_prob[1,,drop=FALSE])) # исправленный код
  
  # Расчет фильтрованной вероятности для первого значения
  # fProb[1,]= ((P%*%matrix(object["iniProb"],ncol=1))*t(Likel[1,,drop=FALSE]))/margLik[1,1] # исходный код
  fProb[1,]= ((P%*%matrix(c(1,1),ncol=1))*t(rec_prob[1,,drop=FALSE]))/margLik[1,1] # исправленный код
  for (i in 2:nr){
    # Расчет смешанной вероятности
    # Расчет совокупной вероятности (нормализирующего коэффициента) для периода с учетом предыдущего
    margLik[i,1]=sum ((P%*%t(fProb[i-1,,drop=FALSE])) * t(rec_prob[i,,drop=FALSE]))
    # Расчет фильтрованной вероятности для периода с учетом предыдущего
    fProb[i,]= ((P%*%t(fProb[i-1,,drop=FALSE])*t(rec_prob[i,,drop=FALSE]))/margLik[i,1])
  }
  # Вывод результата
  return(data.table(fProb = fProb, rec_prob = rec_prob))
}

res <- pred_msm(rec_msm, rec_dt)

datatable(res[1:30], style = 'bootstrap4', extensions = 'Responsive', 
            options = list(pageLength = 6), caption = "Прогнозные значения") |> 
  formatRound(columns = 1:4, digits = 4)

Теперь можно построить график прогноза с использованием фильтрованной вероятности из модели для данных, которые модель “знает”

cbind(res, rec_dt) |>
  tfm(RECPROUSM156N = RECPROUSM156N/100) |> # приведения прогноза Marcelle Chauvet к процентам
  cbind(fltr = rec_msm@Fit@filtProb) |> # фильтрованная вероятность из модели
  tfm(mdl = fifelse(fltr.V2 < fltr.V1, fProb.V1, fProb.V2)) |> # фильтрация прогноза
  e_charts(date) |> 
  e_area(USRECD, symbol= 'none', lineStyle = list(width = 0), 
         name = "Период рецессии", color = "rgba(225, 87, 89, .5)") |> 
  e_line(RECPROUSM156N, symbol= 'none', color = my_pal[1], 
         name = "Вероятность рецессии Marcelle Chauvet") |>
  e_line(mdl, symbol= 'none', color = my_pal[2], lineStyle = list(type = 'dashed', width = 1),
         name = "Вероятность рецессии авторская") |>
  e_datazoom(x_index = 0, type = "slider") |>
  e_title(text = "Оценка вероятности рецесии в США") |> 
  e_y_axis(name = "Вероятность рецесии", formatter = e_axis_formatter("percent")) |> 
  e_x_axis(name = "Месяц") |> 
  e_draft("InvestCookies.ru", size = "80px", opacity = 0.2) |> 
  e_legend(padding = 30) |> 
  e_tooltip(trigger = "axis", axisPointer = list(type = "line"), 
            backgroundColor = "rgba(255,255,255,0.7)", 
            formatter = e_tooltip_pointer_formatter("percent")) 

В целом авторская модель ничуть не хуже справляется с прогнозированием по сравнению с официально публикуемой на сайте ФРС США, но есть одна неприятная тонкость. Дело в том, что фильтрованная вероятность рассчитывается в алгоритме, отталкиваясь от конкретного значения отсутствия/наличия рецеcсии на каждый момент времени. Это происходит в строчке кода: Likel= object@Likelihood(as.matrix(model$model[,1,drop=FALSE]) %*% matrix(rep(1,k),nrow=1), mu=CondMean) метода .MSM.glm.msmFilter. Следовательно невозможно вычислить значение фильтрованной вероятности без конкретного значения предсказываемой переменной. Звучит как полный бред, но это так: красивый график фильтрации можно получить, показав модели фактические значения. Как вариант, можно показывать модели значения только за предыдущий месяц, но это никак не решает проблемы когда задержка в публикации наличия рецессии составлет пол года или больше.

Я бы не хотел останавливаться на мысли о том, что модель абсолютно бестолковая, но собираюсь построить еще один график для прогнозных значений, для которых не нужно показывать фактические значения. В силу того, что модель рассчитывает вероятность от предыдущего состояния важно понимать от какой точки производится расчет и в данном случае начальной точкой будет далекий 1970 год, а завершением будет текущий месяц 2022 года. Другими словами, модель рассчитывает вероятность рецессии на каждом шаге исходя из значений предикторов и также предыдущего состояния т.е. прогнозные значения представляют собой цепь связанных событий ⛓

new_dt <- indctrs0  %>%  
  # Получение относительных изменений
  .[, stri_c(preds, "_ch") := lapply(.SD, \(x) shift((x - shift(x))/x)), .SDcols = preds]  %>% 
  # Получение сдвига на один период
  .[, stri_c(c(preds, stri_c(preds, "_ch")), "lg") := lapply(.SD, shift), .SDcols = c(preds, stri_c(preds, "_ch"))] |> 
  fsubset(date > as.Date("1970-01-01")) |>
  na_omit()

pred <- pred_msm(rec_msm, new_dt)

cbind(pred, new_dt)[, n := .I]  |> 
  tfm(mdl = fifelse(fProb.V2 > fProb.V1, rec_prob.V1, rec_prob.V2)) |> # фильтрация прогноза
  e_charts(date) |> 
  e_area(USRECD, symbol= 'none', lineStyle = list(width = 0), 
         name = "Период рецессии", color = "rgba(225, 87, 89, .5)") |> 
  e_line(mdl, symbol= 'none', color = my_pal[2], lineStyle = list(type = 'dashed', width = 1),
         name = "Вероятность рецессии авторская") |>
  e_loess(mdl ~ n, lineStyle = list(type = "dotted", width = 1), model_args = list(span = .005),
          emphasis = list(focus = 'self'), showSymbol = FALSE, color = my_pal[1], 
          name = "Сглаженная вероятность рецессии авторская") |>
  e_datazoom(x_index = 0, type = "slider") |>
  e_title(text = "Оценка вероятности рецесии в США") |> 
  e_y_axis(name = "Вероятность рецесии", formatter = e_axis_formatter("percent")) |> 
  e_x_axis(name = "Месяц") |> 
  e_draft("InvestCookies.ru", size = "80px", opacity = 0.2) |> 
  e_legend(padding = 30) |> 
  e_tooltip(trigger = "axis", axisPointer = list(type = "line"), 
            backgroundColor = "rgba(255,255,255,0.7)", 
            formatter = e_tooltip_pointer_formatter("percent")) |> 
  e_mark_area(data = list(list(xAxis = as.Date("2022-01-01"), name = "Прогноз"), 
                          list(xAxis = as.Date("2022-07-01"))), 
              itemStyle = list(color = my_pal[1], opacity = .2))

Модель уже не такая уверенная в своих предсказаниях, тем не менее, картинка местами неплохо указывает на рецессию. В качестве характеристик модели можно использовать стандартные метрики. Функцию для вычисления метрик уже приводилась в предыдущей заметке и следовательно не требует дополнительных пояснений:

# Функция которая считает полезные метрики модели
get_metrics <- function(dt){
  mtrx <- table(dt)
  n <- fnrow(dt)
  accuracy <- round(fsum(diag(mtrx))/n, 2)
  precision = round(diag(mtrx) / colSums(mtrx) , 2)
  recall = round(diag(mtrx) / rowSums(mtrx), 2)
  f1 = round(2 * precision * recall / (precision + recall) , 2)
  
  qDF(mtrx) |> cbind(data.frame(metrics = c("->", "->"), precision, recall, f1)) |> 
    hl(scale_color_fira(continuous = TRUE), cols = 1:2, calc_scale = "each") |>
    as_html(style = "color:White")
}

cbind(pred, rec_dt)[date < as.Date("2022-01-01")]  |> 
  tfm(mdl = fifelse(fProb.V2 > fProb.V1, rec_prob.V1, rec_prob.V2)) |> 
  tfm(pred  = fifelse(mdl > .5, 1, 0)) |> 
  fselect(USRECD, pred) |>
  get_metrics()
      0  1 metrics precision recall   f1
0   539  0      ->       0.9   1.00 0.95
1    57 33      ->       1.0   0.37 0.54

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

tfm(rec_dt, fit = rec_glm$fitted.values) |>
  tfm(resid = fit - USRECD,
      pred = fifelse(fit > .26, 1, 0)) |> # предельное значения для определения рецессии = 26%
  fselect(actual = USRECD, prediction = pred)|> 
  get_metrics()
      0  1 metrics precision recall   f1
0   502 37      ->      0.93   0.93 0.93
1    37 47      ->      0.56   0.56 0.56

Отклик (Recall) модели существенно ухудшился, но именно эта характеристика – является ключевой с моей точки зрения потому как важнее знать о проблемах в экономике чем об их отсутствии применительно к инвестициям в целом и к инвестиционному портфелю в частности 💰

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

Выводы 🍪

По заведенной традиции сделаю краткие выводы по заметке в целом:

  1. Я достаточно подробно разобрал модель динамического изменения факторов на основе цепи Маркова, которая используется на сайте ФРС для прогнозирования наступления рецессии
  2. Модель показывает слабый показатель отклика (чувствительности) на период рецессии: метрика оказалась даже хуже простой логистической регрессии, которую я приводил в качестве некой базовой модели
  3. Модель прогнозирует отсутствие рецессии в США для первых шести месяцев 2022 года
  4. И все-таки эта модель какая-то бестолковая для предсказания рецессии и я полагаю, что можно найти какую-то более вменяемую альтернативу

Простой способ узнать о новых публикациях – подписаться на Telegram-канал:


  1. Hamilton Regime Switching Model using R code ↩︎

  2. The Markov Switching Dynamic Regression Model ↩︎

  3. AN ECONOMETRIC CHARACTERIZATION OF BUSINESS CYCLE DYNAMICS WITH FACTOR STRUCTURE AND REGIME SWITCHING ↩︎