計算從s1：s400到sn中的大馬爾可夫轉換矩陣的概率需要永久

我成功實現了以下任務，但我只是將它用於小型測試數據集。對於我的真實數據集，它只是不斷計算。計算從s1：s400到sn中的大馬爾可夫轉換矩陣的概率需要永久

我在R中有一個400x400的轉換概率矩陣。如果用戶在馬爾可夫步行上轉換，用戶點擊「轉換」。所有用戶的吸收狀態爲「空」。「開始」是我的開始狀態。

兩件事情我需要計算：

命中狀態s_j上隨機遊走在「開始」開始
命中「轉換」上在每個397個其他國家開始隨機遊走

第一個是容易R：

v <- numeric(length = ncol(transitionMatrix1)) 
v[1] <- 1 
i <- 2 
R0 <- v%*%(transitionMatrix1 %^% 1) 
R <- R0 
repeat { 
    R1 <- v%*%(transitionMatrix1 %^% i) 
    R <- rbind(R, R1) 
    if (rowSums(R[nrow(R)-1,] - R1) == 0) { 
    #if (rowSums(R[nrow(R)-1,] - R1) < epsilon) { 
    break 
    } 
    else { 
    i <- i+1 
    } 
} 
visit1 <- colSums(R)

我成功地實施2，但我只有它與一個小矩陣一起工作。它需要一個大的永遠：

w <- i 
C1 <- matrix(nrow = w, ncol = ncol(transitionMatrix1)) 
for (i in 1:ncol(transitionMatrix1)) { 
    x <- numeric(length = ncol(transitionMatrix1)) 
    x[i] <- 1 
    for (j in 1:w) { 
    C1[j,i] <- x%*%(transitionMatrix1 %^% j)[,ncol(transitionMatrix1)-1] 
    } 
} 
convert1 <- colSums(C1)

我應該不使用循環。不幸的是，我沒有成功實現這些操作的矢量化。

來源

2016-04-25 Sebastian Lechner

什麼可以幫助我在同一時間保持兩個矩陣在內存中 - 原來的一個和轉置一個，並相應地使用一個或另一個 –

雖然我不完全明白你想要計算什麼，我認爲你可以得到你想要的表達式，而不需要重複乘上矩陣。有關快速概述，請參閱https://en.wikipedia.org/wiki/Absorbing_Markov_chain並進行完整分析，請參閱http://www.math.pku.edu.cn/teachers/yaoy/Fall2011/Kemeny-Snell1976.pdf。你需要計算基本矩陣，然後你可以獲得各種性能指標。 – Forzaa

感謝@Forzaa，與你鏈接的幫助下，我能找到在MATLAB仿真（https://de.mathworks.com/matlabcentral/newsreader/view_thread/48984）一個解決方案，我轉換爲R：

f.absorb <- function(P) { 
    # input: P an absorbing Markov transition matrix in 'from rows to column' form. Rows are 
    # probability vectors and must sum to 1.0. At least one main diagonal 
    # element=1 output: 
    # tau: Time before absorption from a transient state. 
    # B: P(ending up in absorbing state from a transient state) 
    # tau2: variance matrix for tau; 
    # N: fundamental matrix 
    # N2: covariance matrix for N; 
    # CP: transition matrix in canonical form. 
    # Q: Submatrix of transient states, 
    # R: from transient to absorbing state submatrix, 
    # H: hij the probability process will ever go to transient state j starting 
    # in transient state i (not counting the initial state (K & S, p. 61)) 
    # MR: expected number of times that a Markov process remains in a transient 
    # state once entered (including entering step) 
    # VR: Variance for MR (K & S, Theorem 3.5.6, p. 
    # 61) written by E. Gallagher, Environmental Sciences Program UMASS/Boston 
    # Internet: [email protected] written 9/26/93, revised: 10/26/94 
    # refs: Kemeny, J. G. and J. L. Snell. 1976. Finite Markov chains. 
    # Springer-Verlag, New York, New York, U.S.A. Roberts, F. 1976. Discrete 
    # Mathematical Models with applications to social, biological, and 
    # environmental problems. Prentice-Hall 
    # R adaption from MatLab Code found at 
    # https://de.mathworks.com/matlabcentral/newsreader/view_thread/48984 
    pdim = nrow(P) 
    if (rowSums ((colSums(t(P))-matrix(1,1,pdim)) * (colSums(t(P))-matrix(1,1,pdim))) > 0.1) { 
    stop('Rows of P must sum to 1.0') 
    } 
    dP=diag(P); 
    ri=which(dP==1.0); 
    if (is.null(ri)) { 
    stop('No absorbing states (reenter P or use ergodic.m)') 
    } 
    rdim=length(ri); 
    qi=which(dP!=1.0); 
    qdim=length(qi); 
    I=diag(qdim) 
    Q=P[qi,qi]; 
    N=solve(I-Q); # the fundamental matrix 
    tau=data.matrix(colSums(t(N))); 
    CP=P[c(t(ri),t(qi)),c(t(ri),t(qi))]; # the canonical form of P 
    R=CP[(rdim+1):pdim,1:rdim]; 
    B=N%*%R; 
    Ndg=diag(diag(N)); 
    N2=N%*%(2*Ndg-I)-N*N 
    tau2=(2*N-I)%*%tau-tau*tau; 
    H=solve(Ndg,N-I); 
    dQ=diag(Q); 
    oneq=matrix(1,qdim,1); 
    MR=oneq/(oneq-dQ); 
    VR=(dQ/(oneq-dQ)) * (dQ/(oneq-dQ)); 
    list(tau=tau, N=N, B=B, CP=CP, Q=Q, R=R, N2=N2, H=H, MR=MR, VR=VR) 
}

我正在尋找的值可以發現在矩陣中H

來源

2016-04-27 10:47:27

計算從s1：s400到sn中的大馬爾可夫轉換矩陣的概率需要永久

回答

相關問題