Блог

Я хочу ввести в R следующую матрицу (вероятность перехода), но что-то я делаю неправильно, и результат неправильный. Матрица изображена на рисунке :

мои усилия заключаются в следующем :

 n=100
A = matrix(c(0),n 1,n 1);A
A[1,2] = 1/4 
A[1,1] = 3/4
A[2,1] = 1/4 
A[2,2] = 1/2 
A[2,1] = 1/4
for (k in 2:n) {
  A[k,k 1] = 1/4
  A[k,k]   = 1/2 
  A[k,k-1] = 1/6 
  A[k,k-2] = 1/12
  A[n,n]   = 3/4 
  A[n,n-1] = 1/6 
  A[n,n-2] = 1/12}
A


pA = A
for (i in 10000){
  pA = pA %*% A }
pA

 

но результирующие столбцы (первые 3 столбца) должны быть: 0,2087, 0,1652, 0,1307

Какая-нибудь помощь?

Вы могли бы попробовать немного другой подход здесь:

 n <- 100
A <- matrix(0, n   1, n   1)

# create index for the diagonal and the "other" diagonales
A_diag <- matrix(c(1:(n   1), 1:(n   1)), ncol = 2)
A_diag_1 <- matrix(c(1:(n   1), 1   1:(n   1)), ncol = 2)
A_diag_minus_1 <- matrix(c(1:(n   1), - 1   1:(n   1)), ncol = 2)
A_diag_minus_2 <- matrix(c(1:(n   1), - 2   1:(n   1)), ncol = 2)

# populate those diagonales
A[A_diag] <- 1/2
A[A_diag_1[A_diag_1[,2] <= n   1,]] <- 1/4
A[A_diag_minus_1[A_diag_minus_1[,2] >= 1,]] <- 1/6
A[A_diag_minus_2[A_diag_minus_2[,2] >= 1,]] <- 1/12

# create starting values
A[1, 1] <- 3/4
A[2, 1] <- 1/4

# calculate pA
pA <- A
for (i in 1:10000){ pA <- pA %*% A }
 

Это возвращает

 pA[1, 1:3]
#> [1] 0.2087122 0.1651514 0.1306823
 

Для успокоения pA вы могли бы использовать expm пакет, чтобы ускорить процесс:

 library(expm)

pB <- A %^% 10000

pB[1,1:3]
#> [1] 0.2087122 0.1651514 0.1306823