如何在 Stata 或 R 中计算"compound"马尔可夫转移矩阵？

Question

所谓"复合"，我的意思是转移矩阵满足马尔可夫性质，即我有两列s_t和s_t+k，它们分别表示每个个体在两个周期t和t+k中的状态。我想要的是找到的矩阵M

s_t+k = M^k * s_t

使得矩阵CCD_ 5满足马尔可夫性质。

我的默认工作语言是Stata，其中像tab、svy:tab或xttran这样的命令可以生成一个周期的转换矩阵，但这些矩阵不一定满足马尔可夫性质。所以我想知道如何用Stata或其他通用语言(如R或Python)来实现我的目标。

PS：这个问题是从一篇研究1960-2010年许多国家GDP_per_cpita转型动态的论文中提出的。比方说，在每几十年之初，我们将所有国家分为5组(从1组：赤贫国家到5组：高收入国家)，因此我们有一个由5个州组成的国家分布。如果我使用markovchain类简单地估计十年到十年的转换矩阵，这很容易。然而，提交人声称(第11页，脚注4)

"根据利用数值优化程序。而不是简单的平均值对于五个转移矩阵(其受到Jensen不等式)，我们估计一个转移矩阵，它可以给我们一个精确的5年持续时间过渡矩阵(1960年进入，2010年退出)通过采取其功率5。">

Answer 1

在R中，可以使用markovchain包来获得满足马尔可夫性质的转移矩阵。您可以使用以下示例代码。。。

library(markovchain)
data(rain)
mysequence<-rain$rain
createSequenceMatrix(mysequence)
myFit<-markovchainFit(data=mysequence,,method="bootstrap",nboot=5, name="Bootstrap Mc")
myFit

myFit是您估计的过渡矩阵。本例使用Alofi降雨量数据集。

Answer 2

矩阵在R中的乘积不是*，而是%*%

我在R中编写了一个简单的函数来解决这个问题。

trans_mat=函数(k，s_t，M){对于(1:k中的i){M=M%*%M}return(M%*%s_t)}

现在，您需要做的是键入k(您想要的周期有多长)、s_t(原始状态)和M(马尔可夫属性)。

s_t+k=trans_mat(k，s_t，M)

Answer 3

markovchain包直接实现任何markovchainobject的power：

require(markovchain)
#creating the MC
myMatr<-matrix(data=c(0.2,0.8,.6,.4),ncol=2,byrow=TRUE)
myMc<-as(myMatr,"markovchain")
#5th power of the MC
myMc5<-myMc^5
myMc5