R中马尔可夫链的手动模拟 [英] Manual simulation of Markov Chain in R
问题描述
考虑具有状态空间S = {1, 2}、转移矩阵的马尔可夫链
Consider the Markov chain with state space S = {1, 2}, transition matrix
和初始分布α = (1/2, 1/2).
模拟马尔可夫链的5个步骤(即模拟X0、X1, . . . , X5).重复模拟 100次.使用您的模拟结果解决以下问题.
Simulate 5 steps of the Markov chain (that is, simulate X0, X1, . . . , X5). Repeat the simulation 100 times. Use the results of your simulations to solve the following problems.
- 估计P(X1 = 1|X0 = 1).将您的结果与准确概率进行比较.
- Estimate P(X1 = 1|X0 = 1). Compare your result with the exact probability.
我的解决方案:
# returns Xn
func2 <- function(alpha1, mat1, n1)
{
xn <- alpha1 %*% matrixpower(mat1, n1+1)
return (xn)
}
alpha <- c(0.5, 0.5)
mat <- matrix(c(0.5, 0.5, 0, 1), nrow=2, ncol=2)
n <- 10
for (variable in 1:100)
{
print(func2(alpha, mat, n))
}
<小时>
如果我运行此代码一次或 100 次(如问题陈述中所述)有什么区别?
What is the difference if I run this code once or 100 times (as is said in the problem-statement)?
如何从这里找到条件概率?
How can I find the conditional probability from here on?
推荐答案
让
alpha <- c(1, 1) / 2
mat <- matrix(c(1 / 2, 0, 1 / 2, 1), nrow = 2, ncol = 2) # Different than yours
是初始分布和转移矩阵.您的 func2
仅找到不需要的第 n 步分布,并且不会模拟任何内容.相反,我们可以使用
be the initial distribution and the transition matrix. Your func2
only finds n-th step distribution, which isn't needed, and doesn't simulate anything. Instead we may use
chainSim <- function(alpha, mat, n) {
out <- numeric(n)
out[1] <- sample(1:2, 1, prob = alpha)
for(i in 2:n)
out[i] <- sample(1:2, 1, prob = mat[out[i - 1], ])
out
}
其中 out[1]
仅使用初始分布生成,然后对于后续项,我们使用转移矩阵.
where out[1]
is generated using only the initial distribution and then for subsequent terms we use the transition matrix.
然后我们有
set.seed(1)
# Doing once
chainSim(alpha, mat, 1 + 5)
# [1] 2 2 2 2 2 2
所以链从 2 开始并由于指定的转移概率而卡在那里.
so that the chain initiated at 2 and got stuck there due to the specified transition probabilities.
重复 100 次
# Doing 100 times
sim <- replicate(chainSim(alpha, mat, 1 + 5), n = 100)
rowMeans(sim - 1)
# [1] 0.52 0.78 0.87 0.94 0.99 1.00
最后一行显示了我们进入状态 2 而不是状态 1 的频率.这给出了为什么 100 次重复提供更多信息的一个(非常多的)原因:我们陷入状态 2 只做一次模拟,同时重复100 次我们探索了更多可能的路径.
where the last line shows how often we ended up in state 2 rather than 1. That gives one (out of many) reasons why 100 repetitions are more informative: we got stuck at state 2 doing just a single simulation, while repeating it for 100 times we explored more possible paths.
然后可以找到条件概率
mean(sim[2, sim[1, ] == 1] == 1)
# [1] 0.4583333
而真实概率为 0.5(由转移矩阵的左上角条目给出).
while the true probability is 0.5 (given by the upper left entry of the transition matrix).
这篇关于R中马尔可夫链的手动模拟的文章就介绍到这了,希望我们推荐的答案对大家有所帮助,也希望大家多多支持IT屋!