防止Gillespie SSA隨機模型運行陰性

我使用Gillespie SSA產生了感染（寄生蠕蟲）的隨機模型。該模型使用「GillespieSSA」包（https://cran.r-project.org/web/packages/GillespieSSA/index.html）。簡而言之，代碼模擬了不連續隔間的種羣。隔室之間的運動取決於用戶定義的速率方程。 SSA算法用於計算給定時間步（tau）由每個速率方程產生的事件數，並相應地更新總體，過程重複到給定時間點。問題是，假設事件的數量是泊松分佈（Poisson（rate [i] * tau）），因此當速率爲負時（包括羣體數量變爲負數時）會產生錯誤。防止Gillespie SSA隨機模型運行陰性

# Parameter Values 
sir.parms <- c(deltaHinfinity=0.00299, CHi=0.00586, deltaH0=0.0854, aH=0.5, 
       muH=0.02, SigmaW=0.1, SigmaM =0.8, SigmaL=104, phi=1.15, f = 0.6674, 
       deltaVo=0.0166, CVo=0.0205, alphaVo=0.5968, beta=52, mbeta=7300 ,muV=52, g=0.0096, N=100) 
# Inital Population Values 
sir.x0 <- c(W=20,M=10,L=0.02) 
# Rate Equations 
sir.a <- c("((deltaH0+deltaHinfinity*CHi*mbeta*L)/(1+CHi*mbeta*L))*mbeta*L*N" 
      ,"SigmaW*W*N", "muH*W*N", "((1/2)*phi*f)*W*N", "SigmaM*M*N", "muH*M*N", 
      "(deltaVo/(1+CVo*M))*beta*M*N", "SigmaL*L*N", "muV*L*N", "alphaVo*M*L*N", "(aH/g)*L*N") 
# Population change for even 
sir.nu <- matrix(c(+0.01,0,0, 
        -0.01,0,0, 
        -0.01,0,0, 
        0,+0.01,0, 
        0,-0.01,0, 
        0,-0.01,0, 
        0,0,+0.01/230, 
        0,0,-0.01/230, 
        0,0,-0.01/230, 
        0,0,-0.01/230, 
        0,0,-0.01/32),nrow=3,ncol=11,byrow=FALSE) 
runs <- 10 
set.seed(1) 

# Data Frame of output 
sir.out <- data.frame(time=numeric(),W=numeric(),M=numeric(),L=numeric()) 
# Multiple runs and combining data and SSA methods 
for(i in 1:runs){ 
    sim <- ssa(sir.x0,sir.a,sir.nu,sir.parms, method="ETL", tau=1/12, tf=140, simName="SIR") 
    sim.out <- data.frame(time=sim$data[,1],W=sim$data[,2],M=sim$data[,3],L=sim$data[,4]) 

    sim.out$run <- i 
    sir.out <- rbind(sir.out,sim.out) 
}

因此，速率是計算的和模型更新羣體值對於每個時間步驟中，在數據幀中的數據存儲，然後與先前運行連接在一起。然而，當人口的水平變得非常低時，可能發生事件，使得減少人口的事件的數量大於隔室中的事件的數量。一種方法是使時間步長非常小，但這會大大增加模擬的時間長度。

我的問題是有沒有辦法增加代碼，以便在每個時間步驟創建/計算數據時，負數的任何人口數值都將轉換爲0？

我已經嘗試過解決這個問題，但似乎只能在模擬完成後改變這些值的方法，而負值仍然會導致運行本身出現問題。例如
如果（sir.out $ L < 0）sir.out $ L == 0

任何幫助，將不勝感激

來源

2016-06-10 Rnought

我相信這個問題是在SSA設置方法（「ETL」）功能。 ETL方法最終會產生負數。你可以試試「OTL」的方法，基於高效的步長選擇的tau蛋白跳躍模擬方法 - 其中有一些你可以調整幾個參數，但基本命令是：

ssa(sir.x0,sir.a,sir.nu,sir.parms, method="OTL", tf=140, simName="SIR")

或直接法，不會產生任何負數：

ssa(sir.x0,sir.a,sir.nu,sir.parms, method="D", tf=140, simName="SIR")

來源

2016-09-13 15:53:31 fernandosm

防止Gillespie SSA隨機模型運行陰性

回答

相關問題