【此解法中求各n阶概率转移矩阵的R代码】
library(MASS)
data<-c(0,0,0.5,0.5,1,0,0,0,0,1,0,0,0,0,1,0)
p<-matrix(data,nrow=4,ncol=4,byrow=T)
p2=p%*%p
p3=p2%*%p
p4=p3%*%p
p5=p4%*%p
p6=p5%*%p
p7=p6%*%p
p8=p7%*%p
p9=p8%*%p
p10=p9%*%p
fractions(p10)
p
p2
p3
p4
p5
p6
p7
p8
p9
p10
【运行结果】
> library(MASS)
> data<-c(0,0,0.5,0.5,1,0,0,0,0,1,0,0,0,0,1,0)
> p<-matrix(data,nrow=4,ncol=4,byrow=T)
> p2=p%*%p
> p3=p2%*%p
> p4=p3%*%p
> p5=p4%*%p
> p6=p5%*%p
> p7=p6%*%p
> p8=p7%*%p
> p9=p8%*%p
> p10=p9%*%p
> fractions(p10)
[,1] [,2] [,3] [,4]
[1,] 3/8 3/8 3/16 1/16
[2,] 1/8 3/8 3/8 1/8
[3,] 1/4 1/8 3/8 1/4
[4,] 1/2 1/4 1/8 1/8
>
> p
[,1] [,2] [,3] [,4]
[1,] 0 0 0.5 0.5
[2,] 1 0 0.0 0.0
[3,] 0 1 0.0 0.0
[4,] 0 0 1.0 0.0
> p2
[,1] [,2] [,3] [,4]
[1,] 0 0.5 0.5 0.0
[2,] 0 0.0 0.5 0.5
[3,] 1 0.0 0.0 0.0
[4,] 0 1.0 0.0 0.0
> p3
[,1] [,2] [,3] [,4]
[1,] 0.5 0.5 0.0 0.0
[2,] 0.0 0.5 0.5 0.0
[3,] 0.0 0.0 0.5 0.5
[4,] 1.0 0.0 0.0 0.0
> p4
[,1] [,2] [,3] [,4]
[1,] 0.5 0.0 0.25 0.25
[2,] 0.5 0.5 0.00 0.00
[3,] 0.0 0.5 0.50 0.00
[4,] 0.0 0.0 0.50 0.50
> p5
[,1] [,2] [,3] [,4]
[1,] 0.0 0.25 0.50 0.25
[2,] 0.5 0.00 0.25 0.25
[3,] 0.5 0.50 0.00 0.00
[4,] 0.0 0.50 0.50 0.00
> p6
[,1] [,2] [,3] [,4]
[1,] 0.25 0.50 0.25 0.00
[2,] 0.00 0.25 0.50 0.25
[3,] 0.50 0.00 0.25 0.25
[4,] 0.50 0.50 0.00 0.00
> p7
[,1] [,2] [,3] [,4]
[1,] 0.50 0.25 0.125 0.125
[2,] 0.25 0.50 0.250 0.000
[3,] 0.00 0.25 0.500 0.250
[4,] 0.50 0.00 0.250 0.250
> p8
[,1] [,2] [,3] [,4]
[1,] 0.25 0.125 0.375 0.250
[2,] 0.50 0.250 0.125 0.125
[3,] 0.25 0.500 0.250 0.000
[4,] 0.00 0.250 0.500 0.250
> p9
[,1] [,2] [,3] [,4]
[1,] 0.125 0.375 0.375 0.125
[2,] 0.250 0.125 0.375 0.250
[3,] 0.500 0.250 0.125 0.125
[4,] 0.250 0.500 0.250 0.000
> p10
[,1] [,2] [,3] [,4]
[1,] 0.375 0.375 0.1875 0.0625
[2,] 0.125 0.375 0.3750 0.1250
[3,] 0.250 0.125 0.3750 0.2500
[4,] 0.500 0.250 0.1250 0.1250
【此解法中求解线性方程组的R代码】
library(MASS)
data<-c(1,-1,0,0,0,1,-1,0,0.5,0,-1,1,1,1,1,1)
P<-matrix(data,nrow=4,ncol=4,byrow=T)
b<-matrix(c(0,0,0,1),nrow=4,ncol=1,byrow=T)
x<-solve(P,b)
fractions(x)
【运行结果】
> library(MASS)
> data<-c(1,-1,0,0,0,1,-1,0,0.5,0,-1,1,1,1,1,1)
> P<-matrix(data,nrow=4,ncol=4,byrow=T)
> b<-matrix(c(0,0,0,1),nrow=4,ncol=1,byrow=T)
> x<-solve(P,b)
> fractions(x)
[,1]
[1,] 2/7
[2,] 2/7
[3,] 2/7
[4,] 1/7