R function for carrying out hybrid continuous deterministic / discrete stochastic simulations of the logistic population model. http://cnr.lwlss.net/DiscreteStochasticLogistic/ #R #Rstats #discrete #stochastic #logistic

simDSLogisticHybrid.R

simCellsHybrid=function(K,r,N0,NSwitch,detpts=100){
  # Every event produces one cell and consumes one unit of nutrients
  if(NSwitch>N0){
    # Unusually, for this model, we know the number of events a priori
    eventNo=NSwitch-N0
    # So we can just generate all required random numbers (quickly) in one go
    unifs=runif(eventNo)
    clist=(N0+1):NSwitch
    # Time between events
    dts=-log(1-unifs)/(r*clist*(1-clist/K))
    # Absolute times
    ats=cumsum(dts)
    tmax=max(ats)
  }else{
    clist=c()
    ats=c()
    tmax=0
  }
  # Switch to discrete deterministic logistic function
  clistdet=seq(NSwitch+(K-NSwitch)/detpts,K,(K-NSwitch)/detpts)
  tsdet=log((clistdet*(K - NSwitch))/((K - clistdet)*NSwitch))/r
  return(data.frame(t=c(0,ats,tmax+tsdet),c=c(N0,c(clist,clistdet))))
}