ページ

2016年4月25日月曜日

【R】マルコフ解析の基礎(3)

Getting Started with Markov Chains

前回はマルコフ連鎖の状態空間に含まれる状態を分類してみた。今回はマルコフ連鎖そのものを分解しつつ分類してみる。吸収的マルコフ連鎖、エルゴード的マルコフ連鎖、周期的マルコフ連鎖の3種が重要だ。

いずれもt→∞での挙動に注目した分類となる。吸収的マルコフ連鎖は、いずれかの吸収状態、いわばブラックホールに吸い込まれて終わるパターン。エルゴード的、周期的は吸収状態を持たないが、エルゴード的の場合は状態分布が特定の定常分布に収束していく、周期的は一つの定常分布に収束するのではなく、複数の状態分布の間を周期的にいったりきたりするものだ。

今回はRを使って状態遷移図を描きながら解説していこう。



diagramパッケージで状態遷移図を描く

markovchainパッケージはマルコフ連鎖を分析する専用パッケージで、状態遷移図も描けるのだが、ややしょぼい。そこで、別のパッケージで状態遷移図を描くことにする。Joseph Rickertさんの記事(冒頭のリンク「Getting Started with Markov Chains」)を参考にした。

まずはパッケージを管理者権限でインストールして、パッケージを呼び出し、plotmat関数で表示させる。
# -----------------------------
# Use Library
# -----------------------------
library("diagram")

# -----------------------------
# my_plotmat : plotmatの呼び出し
# -----------------------------
my_plotmat <- function(x){
 plotmat(x,pos = c(1,2), 
        lwd = 1, box.lwd = 2, 
        cex.txt = 0.8, 
        box.size = 0.1, 
        box.type = "circle", 
        box.prop = 0.5,
        box.col = "light yellow",
        arr.length=.1,
        arr.width=.1,
        self.cex = .4,
        self.shifty = -.01,
        self.shiftx = .13,
        main = "")
}

# -----------------------------
# 遷移確率行列の作成
# -----------------------------
Oz<-matrix(c(0.6,0.3,0.1,0.3,0.6,0.1,0.2,0.3,0.5), 
  nrow=3, byrow=TRUE)
stateNames=c("Fine","Cloudy","Rain")
row.names(Oz) <- stateNames
colnames(Oz) <- stateNames

# -----------------------------
# 遷移確率から状態遷移図
# -----------------------------
my_plotmat(Oz)
実行すると、冒頭の図のようになる。これは天気の移り変わりをマルコフ連鎖で表現したものだ。数値は「マルコフ解析」の例題を参考にした。 この例自体は、吸収状態を持たない(なので、エルゴード的か周期的かのどっちか)。

 

状態確率の予測

 

遷移確率行列をt乗することでt期先を予測することが出来る。ここではexpmパッケージの機能を使っている。
# -----------------------------
# 遷移確率行列から高次の状態確率予測
# -----------------------------
library(expm)
Oz3 <- Oz %^% 3
u <- c(1/3,1/3,1/3)
u3 <- u %*% Oz3
3期先の状態を予測している。結果はこうなる。
> round(u3,4)
       Fine Cloudy   Rain
[1,] 0.3967  0.426 0.1773

 

マルコフ連鎖オブジェクトの構築

 

markovchainパッケージの機能を使い、次のようにしてオブジェクトを構築し、状態の判定を行う。
# -----------------------------
# Use Library
# -----------------------------
library("markovchain")

# オブジェクト構築
mcWether<-new("markovchain", states=stateNames, transitionMatrix=Oz)

# 3期先の表示
mcWether^3

# 一時的状態の識別
transientStates(mcWether) 

# 吸収状態の識別
absorbingStates(mcWether) 

# 定常状態の状態確率
steadyStates(mcWether)

結果はこうなる。
> # 3期先の表示
> mcWether^3
Unnamed Markov chain^3 
 A  3 - dimensional discrete Markov Chain with following states 
 Fine Cloudy Rain 
 The transition matrix   (by rows)  is defined as follows 
        Fine Cloudy  Rain
Fine   0.427  0.417 0.156
Cloudy 0.400  0.444 0.156
Rain   0.363  0.417 0.220
 
> # 一時的状態の識別
> transientStates(mcWether) 
character(0)
> 
> # 吸収状態の識別
> absorbingStates(mcWether) 
character(0)
> 
> # 定常状態の状態確率
> steadyStates(mcWether)
          Fine    Cloudy      Rain
[1,] 0.4047619 0.4285714 0.1666667
つまり、1)ここには一時的状態は含まれておらず、各状態を何回も訪問すること、2)吸収状態も含まれておらず、ブラックホールのようなものはないこと、が分かる。

 

まとめ

 

マルコフ連鎖の分類までたどり着かなかった。今回は遷移確率行列と状態遷移図のTipsだけだったが、急がずに解説していきたい。