1#!/usr/local/bin/python3.8
2# coding: utf-8
3
4# DO NOT EDIT
5# Autogenerated from the notebook markov_regression.ipynb.
6# Edit the notebook and then sync the output with this file.
7#
8# flake8: noqa
9# DO NOT EDIT
10
11# ## Markov switching dynamic regression models
12
13# This notebook provides an example of the use of Markov switching models
14# in statsmodels to estimate dynamic regression models with changes in
15# regime. It follows the examples in the Stata Markov switching
16# documentation, which can be found at
17# http://www.stata.com/manuals14/tsmswitch.pdf.
18
19import numpy as np
20import pandas as pd
21import statsmodels.api as sm
22import matplotlib.pyplot as plt
23
24# NBER recessions
25from pandas_datareader.data import DataReader
26from datetime import datetime
27
28usrec = DataReader("USREC",
29                   "fred",
30                   start=datetime(1947, 1, 1),
31                   end=datetime(2013, 4, 1))
32
33# ### Federal funds rate with switching intercept
34#
35# The first example models the federal funds rate as noise around a
36# constant intercept, but where the intercept changes during different
37# regimes. The model is simply:
38#
39# $$r_t = \mu_{S_t} + \varepsilon_t \qquad \varepsilon_t \sim N(0,
40# \sigma^2)$$
41#
42# where $S_t \in \{0, 1\}$, and the regime transitions according to
43#
44# $$ P(S_t = s_t | S_{t-1} = s_{t-1}) =
45# \begin{bmatrix}
46# p_{00} & p_{10} \\
47# 1 - p_{00} & 1 - p_{10}
48# \end{bmatrix}
49# $$
50#
51# We will estimate the parameters of this model by maximum likelihood:
52# $p_{00}, p_{10}, \mu_0, \mu_1, \sigma^2$.
53#
54# The data used in this example can be found at https://www.stata-
55# press.com/data/r14/usmacro.
56
57# Get the federal funds rate data
58from statsmodels.tsa.regime_switching.tests.test_markov_regression import fedfunds
59
60dta_fedfunds = pd.Series(fedfunds,
61                         index=pd.date_range("1954-07-01",
62                                             "2010-10-01",
63                                             freq="QS"))
64
65# Plot the data
66dta_fedfunds.plot(title="Federal funds rate", figsize=(12, 3))
67
68# Fit the model
69# (a switching mean is the default of the MarkovRegession model)
70mod_fedfunds = sm.tsa.MarkovRegression(dta_fedfunds, k_regimes=2)
71res_fedfunds = mod_fedfunds.fit()
72
73res_fedfunds.summary()
74
75# From the summary output, the mean federal funds rate in the first regime
76# (the "low regime") is estimated to be $3.7$ whereas in the "high regime"
77# it is $9.6$. Below we plot the smoothed probabilities of being in the high
78# regime. The model suggests that the 1980's was a time-period in which a
79# high federal funds rate existed.
80
81res_fedfunds.smoothed_marginal_probabilities[1].plot(
82    title="Probability of being in the high regime", figsize=(12, 3))
83
84# From the estimated transition matrix we can calculate the expected
85# duration of a low regime versus a high regime.
86
87print(res_fedfunds.expected_durations)
88
89# A low regime is expected to persist for about fourteen years, whereas
90# the high regime is expected to persist for only about five years.
91
92# ### Federal funds rate with switching intercept and lagged dependent
93# variable
94#
95# The second example augments the previous model to include the lagged
96# value of the federal funds rate.
97#
98# $$r_t = \mu_{S_t} + r_{t-1} \beta_{S_t} + \varepsilon_t \qquad
99# \varepsilon_t \sim N(0, \sigma^2)$$
100#
101# where $S_t \in \{0, 1\}$, and the regime transitions according to
102#
103# $$ P(S_t = s_t | S_{t-1} = s_{t-1}) =
104# \begin{bmatrix}
105# p_{00} & p_{10} \\
106# 1 - p_{00} & 1 - p_{10}
107# \end{bmatrix}
108# $$
109#
110# We will estimate the parameters of this model by maximum likelihood:
111# $p_{00}, p_{10}, \mu_0, \mu_1, \beta_0, \beta_1, \sigma^2$.
112
113# Fit the model
114mod_fedfunds2 = sm.tsa.MarkovRegression(dta_fedfunds.iloc[1:],
115                                        k_regimes=2,
116                                        exog=dta_fedfunds.iloc[:-1])
117res_fedfunds2 = mod_fedfunds2.fit()
118
119res_fedfunds2.summary()
120
121# There are several things to notice from the summary output:
122#
123# 1. The information criteria have decreased substantially, indicating
124# that this model has a better fit than the previous model.
125# 2. The interpretation of the regimes, in terms of the intercept, have
126# switched. Now the first regime has the higher intercept and the second
127# regime has a lower intercept.
128#
129# Examining the smoothed probabilities of the high regime state, we now
130# see quite a bit more variability.
131
132res_fedfunds2.smoothed_marginal_probabilities[0].plot(
133    title="Probability of being in the high regime", figsize=(12, 3))
134
135# Finally, the expected durations of each regime have decreased quite a
136# bit.
137
138print(res_fedfunds2.expected_durations)
139
140# ### Taylor rule with 2 or 3 regimes
141#
142# We now include two additional exogenous variables - a measure of the
143# output gap and a measure of inflation - to estimate a switching Taylor-
144# type rule with both 2 and 3 regimes to see which fits the data better.
145#
146# Because the models can be often difficult to estimate, for the 3-regime
147# model we employ a search over starting parameters to improve results,
148# specifying 20 random search repetitions.
149
150# Get the additional data
151from statsmodels.tsa.regime_switching.tests.test_markov_regression import ogap, inf
152
153dta_ogap = pd.Series(ogap,
154                     index=pd.date_range("1954-07-01", "2010-10-01",
155                                         freq="QS"))
156dta_inf = pd.Series(inf,
157                    index=pd.date_range("1954-07-01", "2010-10-01", freq="QS"))
158
159exog = pd.concat((dta_fedfunds.shift(), dta_ogap, dta_inf), axis=1).iloc[4:]
160
161# Fit the 2-regime model
162mod_fedfunds3 = sm.tsa.MarkovRegression(dta_fedfunds.iloc[4:],
163                                        k_regimes=2,
164                                        exog=exog)
165res_fedfunds3 = mod_fedfunds3.fit()
166
167# Fit the 3-regime model
168np.random.seed(12345)
169mod_fedfunds4 = sm.tsa.MarkovRegression(dta_fedfunds.iloc[4:],
170                                        k_regimes=3,
171                                        exog=exog)
172res_fedfunds4 = mod_fedfunds4.fit(search_reps=20)
173
174res_fedfunds3.summary()
175
176res_fedfunds4.summary()
177
178# Due to lower information criteria, we might prefer the 3-state model,
179# with an interpretation of low-, medium-, and high-interest rate regimes.
180# The smoothed probabilities of each regime are plotted below.
181
182fig, axes = plt.subplots(3, figsize=(10, 7))
183
184ax = axes[0]
185ax.plot(res_fedfunds4.smoothed_marginal_probabilities[0])
186ax.set(title="Smoothed probability of a low-interest rate regime")
187
188ax = axes[1]
189ax.plot(res_fedfunds4.smoothed_marginal_probabilities[1])
190ax.set(title="Smoothed probability of a medium-interest rate regime")
191
192ax = axes[2]
193ax.plot(res_fedfunds4.smoothed_marginal_probabilities[2])
194ax.set(title="Smoothed probability of a high-interest rate regime")
195
196fig.tight_layout()
197
198# ### Switching variances
199#
200# We can also accommodate switching variances. In particular, we consider
201# the model
202#
203# $$
204# y_t = \mu_{S_t} + y_{t-1} \beta_{S_t} + \varepsilon_t \quad
205# \varepsilon_t \sim N(0, \sigma_{S_t}^2)
206# $$
207#
208# We use maximum likelihood to estimate the parameters of this model:
209# $p_{00}, p_{10}, \mu_0, \mu_1, \beta_0, \beta_1, \sigma_0^2, \sigma_1^2$.
210#
211# The application is to absolute returns on stocks, where the data can be
212# found at https://www.stata-press.com/data/r14/snp500.
213
214# Get the federal funds rate data
215from statsmodels.tsa.regime_switching.tests.test_markov_regression import areturns
216
217dta_areturns = pd.Series(areturns,
218                         index=pd.date_range("2004-05-04",
219                                             "2014-5-03",
220                                             freq="W"))
221
222# Plot the data
223dta_areturns.plot(title="Absolute returns, S&P500", figsize=(12, 3))
224
225# Fit the model
226mod_areturns = sm.tsa.MarkovRegression(
227    dta_areturns.iloc[1:],
228    k_regimes=2,
229    exog=dta_areturns.iloc[:-1],
230    switching_variance=True,
231)
232res_areturns = mod_areturns.fit()
233
234res_areturns.summary()
235
236# The first regime is a low-variance regime and the second regime is a
237# high-variance regime. Below we plot the probabilities of being in the low-
238# variance regime. Between 2008 and 2012 there does not appear to be a clear
239# indication of one regime guiding the economy.
240
241res_areturns.smoothed_marginal_probabilities[0].plot(
242    title="Probability of being in a low-variance regime", figsize=(12, 3))
243