Applied Statistics and Data Mining

markovchainパッケージは、離散状態離散時間マルコフ過程(通称は離散マルコフ過程、あるいはマルコフ連鎖)をRで扱うのに便利な機能を提供する。以下の説明は、https://cran.r-project.org/web/packages/markovchain/vignettes/an_introduction_to_markovchain_package.pdfに基づいている。

• An introduction to Markov chains using R
• Markov chain analysis with R
• Crash introduction
• An introduction

マルコフ過程に基づいて現象の性質を明らかにしたり、未来の時点のおける状態の確率分布を推定することなどを、マルコフ解析という。

マルコフ解析に必要なものは、現在の状態もしくは初期の状態と、対象の状態遷移行列である。

• 状態
• 状態遷移行列

天気の移り変わりを例に説明を進める。

weatherStates <- c("sunny", "cloudy", "rain")

byRow <- TRUE

 weatherMatrix <- matrix(data = c(0.70, 0.2, 0.1,
 0.3, 0.4, 0.3,
 0.2, 0.45, 0.35), byrow = byRow, nrow = 3,
 dimnames = list(weatherStates, weatherStates))

mcWeather <- new("markovchain", states = weatherStates, byrow = byRow,
 transitionMatrix = weatherMatrix, name = "Weather")

mcList <- new("markovchainList", markovchains = list(mcWeather, defaultMc),
 name = "A list of Markov chains")

上で定義したマルコフ連鎖を確認しておく。

weatherStates
weatherMatrix
mcWeather

さて、初期状態は曇りcloudy。

initialState = c(0, 1, 0)

二日後の天気の確率分布を求めてみる。

after2Days <- initialState * (mcWeather * mcWeather)

求めた確率分布を確認する。

after2Days

７日後の天気の条件付き確率分布を推定する。

after7Days <- initialState * (mcWeather ^ 7)

求めた確率分布を確認する。

after7Days

こんなに長い桁は必要ないので、四捨五入して小数点以下第３位まで表示してみる。

round(after7Days, 3)

状態を列ベクトルで表すなら、すべてを転地すればいい。

initialState <- c(0, 1, 0)
after2Days <- (t(mcWeather) * t(mcWeather)) * initialState
after7Days <- (t(mcWeather) ^ 7) * initialState
after2Days

fvals<-function(mchain,initialstate,n) {
out<-data.frame()
names(initialstate)<-names(mchain)
for (i in 0:n)
{
iteration<-initialstate*mchain^(i)
out<-rbind(out,iteration)
}
out<-cbind(out, i=seq(0,n))
out<-out[,c(4,1:3)]
return(out)
}
fvals(mchain=mcWeather,initialstate=c(90,5,5),n=4)

states(mcWeather)

names(mcWeather)

dim(mcWeather)

 name(mcWeather)

markovchain:::sort(mcWeather)

transitionProbability(mcWeather, "cloudy", "rain")

mcWeather[2,3]

print(mcWeather)

show(mcWeather)

plot(mcWeather, package="diagram",box.size = 0.04)

mcDf <- as(mcWeather, "data.frame")
mcNew <- as(mcDf, "markovchain")
mcDf

mcIgraph <- as(mcWeather, "igraph")

require(msm)
Q <- rbind ( c(0, 0.25, 0, 0.25),
c(0.166, 0, 0.166, 0.166),
c(0, 0.25, 0, 0.25),
c(0, 0, 0, 0) )
cavmsm <- msm(state ~ years, subject = PTNUM, data = cav, qmatrix = Q, death = 4)
msmMc <- as(cavmsm, "markovchain")
msmMc

myMatr<-matrix(c(.1,.8,.1,.2,.6,.2,.3,.4,.3), byrow=TRUE, ncol=3)
myMc<-as(myMatr, "markovchain")
myMc

conditionalDistribution(mcWeather, "sunny")

steadyStates(mcWeather)

absorbingStates(mcWeather)

summary(mcWeather)

transientStates(mcWeather)

is.accessible(object = mcWeather, from = "sunny", to = "cloudy")

period(mcWeather)

firstPassage(object = mcWeather, state = "sunny",n = 10)