模拟马尔可夫链路径 1000 次

Question

我有以下马尔可夫链代码：

simulation.mc=function(i0,P,n.sim){
S=1:nrow(P)
X=rep(0,n.sim)
X[1]=i0
for (n in 2:n.sim){
    X[n]=sample(x=S,size=1,prob=P[X[n-1],])
    }
return(X)
}

P=matrix(
c(
    0,1/2,0,1/2,0,0,0,
    1/2,0,1/2,0,0,0,0,
    0,1,0,0,0,0,0,
    1/3,0,0,0,1/3,1/3,0,
    0,0,0,1,0,0,0,
    0,0,0,1/2,0,0,1/2,
    0,0,0,0,0,0,1
),nrow=7,byrow=T);P

X=simulation.mc(1,P,100)


T=min(which(X==7))

我必须计算到达状态 7 之前的平均步数。

我知道我需要运行至少 1000 个路径样本，统计每个样本的步数然后计算平均值（虽然有些路径不会达到 7 状态）。

我这样做了，但还是不行：

n.sim=100
X[i]=rep(0,n.sim)
 for (i in 1:100)
 { X[i]=simulation.mc(1,P,100)
 }

为什么这不起作用？如何在循环中包含循环以包含计算 os 步数的函数？在此先感谢您的任何建议。

Answer 1

您可以使用 replicate 而不是循环：

replicate(1000, min(which(simulation.mc(1,P,100)==7)))

@JDB 提供了一个使用循环的选项。这里还有几个：

# To save each entire chain in a list
n.sim=100
X = list()
for (i in 1:1000) { 
  X[[i]] = simulation.mc(1,P,n.sim)
}

# To save just the number of steps to get to 7
n.sim=100
X = rep(NA, 1000)  
for (i in 1:1000) { 
  X[i] = min(which(simulation.mc(1,P,n.sim)==7))
}

Answer 2

您似乎在尝试创建一个 100 x 100 的数据框，但您调用的是一个已经创建的向量。

这就是您对 X[i]=rep(0,n.sim) 的调用不起作用的原因。（我在想你的意思可能是 X[1]=rep(0,n.sim)，因为你只在循环中定义了 i。你不能像 [=20] 那样填充向量的列=]...

尝试：

X <- data.frame(matrix(nrow=100, ncol=100))
for (i in 1:100)
   { X[i]=simulation.mc(1,P,100)
}

模拟马尔可夫链路径 1000 次

Simulate a Markov chain path 1000 times

r

markov-chains

stochastic-process