在R中使用MSwM包复制哈密尔顿的马尔可夫切换模型。

时间:2022-05-01 15:51:22

I'm trying to estimate the basic Markov Switching Model of Hamilton (1989) as is post in E-views webpage. This model is itself is an exact replication of the existing in RATS.

我正在尝试估计Hamilton(1989)的基本马尔可夫转换模型,在E-views网页中。这个模型本身就是对老鼠存在的精确复制。

This is the time series of the example:

这是例子的时间序列:

gnp <- 
structure(c(2.59316410021381, 2.20217123302681, 0.458275619103479, 
0.968743815568942, -0.241307564718414, 0.896474791426144, 2.05393216767198, 
1.73353647046698, 0.938712869506845, -0.464778333117193, -0.809834082445603, 
-1.39763692441103, -0.398860927649558, 1.1918415768741, 1.4562004729396, 
2.1180822079447, 1.08957867423914, 1.32390272784813, 0.87296368144358, 
-0.197732729861307, 0.45420214345009, 0.0722187603196887, 1.10303634435563, 
0.820974907499614, -0.0579579499110212, 0.584477722838197, -1.56192668045796, 
-2.05041027007508, 0.536371845140342, 2.3367684244086, 2.34014568267516, 
1.23392627573662, 1.88696478737248, -0.459207909351867, 0.84940472194713, 
1.70139850766727, -0.287563102546191, 0.095946277449187, -0.860802907461483, 
1.03447124467041, 1.23685943797014, 1.42004498680119, 2.22410642769683, 
1.3021017302965, 1.0351769691057, 0.925342521818, -0.165599507925585, 
1.3444381723048, 1.37500136316918, 1.73222186043569, 0.716056342342333, 
2.21032138350616, 0.853330335823775, 1.00238777849592, 0.427254413549543, 
2.14368353713136, 1.4378918561536, 1.5795993028646, 2.27469837381376, 
1.95962653201067, 0.2599239932111, 1.01946919515563, 0.490163994319276, 
0.563633789161385, 0.595954621290765, 1.43082852218349, 0.562301244017229, 
1.15388388887095, 1.68722847001462, 0.774382052478202, -0.0964704476805431, 
1.39600141863966, 0.136467982223878, 0.552237133917267, -0.399448716111952, 
-0.61671104590512, -0.0872256083215416, 1.21018349098461, -0.907297546921259, 
2.64916154469762, -0.00806939681695959, 0.511118931407946, -0.00401437145032572, 
2.1682142321342, 1.92586729194597, 1.03504719187207, 1.85897218652101, 
2.32004929969819, 0.255707901889092, -0.0985527428151145, 0.890736834018326, 
-0.55896483237131, 0.283502534230679, -1.31155410054958, -0.882787789285689, 
-1.97454945511993, 1.01275266533046, 1.68264718400186, 1.38271278970291, 
1.86073641586006, 0.444737715592073, 0.414490009766608, 0.992022769383933, 
1.36283572253682, 1.59970527327726, 1.98845814838348, -0.256842316681229, 
0.877869502339381, 3.10956544706826, 0.853244770655281, 1.23337321374495, 
0.0031430232743432, -0.0943336967005583, 0.898833191548979, -0.190366278407953, 
0.997723787687709, -2.39120056095144, 0.0664967330277127, 1.26136016443398, 
1.91637832265846, -0.334802886728505, 0.44207108280265, -1.40664914211265, 
-1.52129894225829, 0.299198686266393, -0.801974492802505, 0.152047924379708, 
0.985850281223592, 2.1303461510993, 1.34397927090998, 1.61550521216825, 
2.70930096486278, 1.24461416484445, 0.508354657516633, 0.148021660957899
), .Tsp = c(1951.25, 1984.75, 4), class = "ts")

I want to use the MSwM package, so I wrote the following code:

我想使用MSwM包,所以我写了以下代码:

library(MSwM) #Load the package  
# Create the model with only an intercept (that after will be switching) 
mod=lm(gnp~1)
# Estimate the Markov Switching Model with only an intercept switching, 
# four lags and two regimes as in Hamilton. 
mod.mswm=msmFit(mod,k=2,p=4,sw=c(T,F,F,F,F,F), control=list(parallel=F)) 
summary(mod.mswm)

I get a result that is very different to obtained in Eviews or RATS:

我得到的结果与在Eviews或RATS中得到的结果非常不同:

Coefficients:
Regime 1 
---------
               Estimate Std. Error t value  Pr(>|t|)    
(Intercept)(S)   0.5747     1.0044  0.5722 0.5671865    
gnp_1            0.3097     0.0903  3.4297 0.0006042 ***
gnp_2            0.1273     0.0900  1.4144 0.1572445    
gnp_3           -0.1213     0.0867 -1.3991 0.1617830    
gnp_4           -0.0892     1.6918 -0.0527 0.9579709    
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 0.98316
Multiple R-squared: 0.1437

Standardized Residuals:
        Min          Q1         Med          Q3         Max 
-1.86974671 -0.37107376  0.03466299  0.39090950  1.67876663 

Regime 2 
---------
               Estimate Std. Error t value  Pr(>|t|)    
(Intercept)(S)   0.5461     1.0044  0.5437 0.5866479    
gnp_1            0.3097     0.0903  3.4297 0.0006042 ***
gnp_2            0.1273     0.0900  1.4144 0.1572445    
gnp_3           -0.1213     0.0867 -1.3991 0.1617830    
gnp_4           -0.0892     1.6918 -0.0527 0.9579709    
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 0.98316
Multiple R-squared: 0.1431

Standardized Residuals:
        Min          Q1         Med          Q3         Max 
-2.51219057 -0.46185366  0.06749067  0.52368275  2.11071358 

Transition probabilities:
          Regime 1  Regime 2
Regime 1 0.3879799 0.3651762
Regime 2 0.6120201 0.6348238

The main difference is obtained in the intercept, because in both regimes a positive value is obtained instead of values in Eviews or RATS. This difference is due to maximization algortihm used (EM in MsWm)? or I have done some mistake in my R-Code?

主要的区别是在截距中得到的,因为在这两种情况下,都得到了一个正值,而不是在Eviews或RATS中得到的值。这种差异是由于algortihm (EM)在MsWm中使用的最大化?或者我在R-Code中犯了一些错误?

Thanks a lot.

非常感谢。

1 个解决方案

#1


4  

The difference that I see is that the model that you are defining contains a switching intercept, while the model of Hamilton (1989) specifies a switching mean instead. That is, your model is:

我看到的不同之处在于,您定义的模型包含一个切换拦截,而Hamilton(1989)的模型则指定了切换的平均值。也就是说,你的模型是:

在R中使用MSwM包复制哈密尔顿的马尔可夫切换模型。

and Hamilton's (1989) model is defined as:

而Hamilton的(1989)模型定义为:

在R中使用MSwM包复制哈密尔顿的马尔可夫切换模型。

In an AR model the parameters alpha and mu will take, in general, different values. This may be somewhat confusing in R as discussed here.

在AR模型中,参数alpha和mu通常取不同的值。在这里讨论的R可能有些令人困惑。

By taking expectations in your model (and omitting for simplicity the switching term S_t) we arrive to the following relationship:

通过在模型中获取期望(并省略了转换项S_t),我们到达了以下关系:

在R中使用MSwM包复制哈密尔顿的马尔可夫切换模型。

Upon this relationship, we could expect to be able to recover the mean. However, in this case the switching intercepts does not lead to the switching means found in Hamilton (1989).

在这段关系中,我们可以期望能够恢复平均值。然而,在这种情况下,切换拦截不会导致在Hamilton(1989)中发现的切换方式。

0.5747 / (1 - sum(c(0.3097, 0.1273, -0.1213, -0.0892)))
#[1] 0.7429864
0.5461 / (1 - sum(c(0.3097, 0.1273, -0.1213, -0.0892)))
#[1] 0.7060116

This mapping can usually be applied, for example, with an AR(4) model:

这种映射通常可以应用,例如,使用AR(4)模型:

fit <- lm(gnp[5:135] ~ 1 + gnp[4:134] + gnp[3:133] + gnp[2:132] + gnp[1:131])
fit
# Coefficients:
# (Intercept)   gnp[4:134]   gnp[3:133]   gnp[2:132]   gnp[1:131]  
#     0.55679      0.30974      0.12726     -0.12126     -0.08923 
#
# the mapping from the intercept to mean leads to a value close to the sample mean
coef(fit)[1]/(1 - sum(coef(fit)[-1]))
# 0.7198458 
mean(gnp)
# 0.7445979
# or close to the mean in an AR(4) model, (labelled as intercept)
arima(gnp, order = c(4,0,0), include.mean = TRUE)
# Coefficients:
#          ar1     ar2      ar3      ar4  intercept
#       0.3188  0.1226  -0.1191  -0.0895     0.7441
# s.e.  0.0860  0.0900   0.0898   0.0872     0.1108

It seems that in this case the model should be defined in terms of the mean in order to get estimates of the switching parameter close to those reported in the reference paper.

在这种情况下,模型应该用平均值来定义,以得到与参考文献中所报告的开关参数接近的估计值。

If the function msmFit allowed as input the result returned by arima, it could be used as follows:

如果该函数msmFit允许输入由arima返回的结果,则可以如下所示:

fit <- arima(gnp, order = c(4,0,0), include.mean = TRUE)
msmFit(fit, k = 2, p = 0, sw = c(T,F,F,F,F,F))

I don't know a straightforward way to define an AR model with mean using lm, which is the output required to use msmFit.

我不知道如何用lm来定义AR模型,这是使用msmFit所需的输出。

I think that this difference in the parameterization of the model is more likely to explain the difference in the results rather than the use of the EM algorithm.

我认为这种模型参数化的差异更有可能解释结果的差异而不是EM算法的使用。

#1


4  

The difference that I see is that the model that you are defining contains a switching intercept, while the model of Hamilton (1989) specifies a switching mean instead. That is, your model is:

我看到的不同之处在于,您定义的模型包含一个切换拦截,而Hamilton(1989)的模型则指定了切换的平均值。也就是说,你的模型是:

在R中使用MSwM包复制哈密尔顿的马尔可夫切换模型。

and Hamilton's (1989) model is defined as:

而Hamilton的(1989)模型定义为:

在R中使用MSwM包复制哈密尔顿的马尔可夫切换模型。

In an AR model the parameters alpha and mu will take, in general, different values. This may be somewhat confusing in R as discussed here.

在AR模型中,参数alpha和mu通常取不同的值。在这里讨论的R可能有些令人困惑。

By taking expectations in your model (and omitting for simplicity the switching term S_t) we arrive to the following relationship:

通过在模型中获取期望(并省略了转换项S_t),我们到达了以下关系:

在R中使用MSwM包复制哈密尔顿的马尔可夫切换模型。

Upon this relationship, we could expect to be able to recover the mean. However, in this case the switching intercepts does not lead to the switching means found in Hamilton (1989).

在这段关系中,我们可以期望能够恢复平均值。然而,在这种情况下,切换拦截不会导致在Hamilton(1989)中发现的切换方式。

0.5747 / (1 - sum(c(0.3097, 0.1273, -0.1213, -0.0892)))
#[1] 0.7429864
0.5461 / (1 - sum(c(0.3097, 0.1273, -0.1213, -0.0892)))
#[1] 0.7060116

This mapping can usually be applied, for example, with an AR(4) model:

这种映射通常可以应用,例如,使用AR(4)模型:

fit <- lm(gnp[5:135] ~ 1 + gnp[4:134] + gnp[3:133] + gnp[2:132] + gnp[1:131])
fit
# Coefficients:
# (Intercept)   gnp[4:134]   gnp[3:133]   gnp[2:132]   gnp[1:131]  
#     0.55679      0.30974      0.12726     -0.12126     -0.08923 
#
# the mapping from the intercept to mean leads to a value close to the sample mean
coef(fit)[1]/(1 - sum(coef(fit)[-1]))
# 0.7198458 
mean(gnp)
# 0.7445979
# or close to the mean in an AR(4) model, (labelled as intercept)
arima(gnp, order = c(4,0,0), include.mean = TRUE)
# Coefficients:
#          ar1     ar2      ar3      ar4  intercept
#       0.3188  0.1226  -0.1191  -0.0895     0.7441
# s.e.  0.0860  0.0900   0.0898   0.0872     0.1108

It seems that in this case the model should be defined in terms of the mean in order to get estimates of the switching parameter close to those reported in the reference paper.

在这种情况下,模型应该用平均值来定义,以得到与参考文献中所报告的开关参数接近的估计值。

If the function msmFit allowed as input the result returned by arima, it could be used as follows:

如果该函数msmFit允许输入由arima返回的结果,则可以如下所示:

fit <- arima(gnp, order = c(4,0,0), include.mean = TRUE)
msmFit(fit, k = 2, p = 0, sw = c(T,F,F,F,F,F))

I don't know a straightforward way to define an AR model with mean using lm, which is the output required to use msmFit.

我不知道如何用lm来定义AR模型,这是使用msmFit所需的输出。

I think that this difference in the parameterization of the model is more likely to explain the difference in the results rather than the use of the EM algorithm.

我认为这种模型参数化的差异更有可能解释结果的差异而不是EM算法的使用。