2017-05-03 29 views
0

我写了一些代码,有一个逻辑增长组件(即随着N接近'承载能力'增长速度较慢,直到它达到'承载能力'停止增长) 。但是,当我在R中运行它似乎并没有工作。一些人口最终大于承载能力。我看了数学和一切OK。所以我认为问题在于dN/dt只是针对每个人群计算一次。有谁知道如何解决这个问题?desolve微分方程不工作

任何帮助将不胜感激!

实施例以下代码:

library('optimbase') 
library('deSolve') 
library('tidyverse') 
K = 150000     #carrying capacity   
deaths = 0.2567534  #death rate   
treesize = 0.23523  #resource size   
K_mat = K*ones(10, 1) #matrix of Ks   
death_mat = deaths*ones(10, 1) #matrix of deathrates 
tree_mat = treesize*ones(11, 1) #matrix of resources 
for_mat <- matrix(rbinom(11 * 10, 1, 0.2), ncol = 11, nrow = 10) #connection 
#matrix of foraging 
parameters <- c(for_mat, tree_mat, death_mat, K_mat) #outline parameters 
N <- runif(10,0,K) 
state <- N #starting state 

nestchange <- function(t, state, parameters){ 
     with(as.list(c(state, parameters)),{ 
    r = for_mat %*% tree_mat - death_mat 
    dNdt = r*N - r*N*(N/K_mat) 
    list(c(dNdt)) 
    }) 
} 
times <- seq(0,100) 

out <- ode (y = state, 
      times = times, 
      func = nestchange, 
      parms = parameters) 

results <- as.data.frame(out) 
results <- gather(results, 'nest', 'N', 2:11) 

ggplot(data=results, 
     aes(x=time, y=N, colour=nest)) + 
    geom_line() + 
    theme_bw() 

回答

1

如果您的功能实际上是

nestchange <- function(t, state, parameters){ 
     with(as.list(c(state, parameters)),{ 
    r <- for_mat %*% tree_mat - death_mat 
    dNdt <- r*state - r*state*(state/K_mat) 
    list(c(dNdt)) 
    }) 
} 

作为ODE求解器传递state到函数在每个时间步骤,然而该函数使用可变用于ODE解算器不会更新的计算,而不是由ODE解算器更新。

+0

完美!这解决了问题!非常感谢! –

+1

很高兴它的工作! – Lyngbakr