How to remove loop from the following model function?-CodePudding

I'm rewriting some code, and I am currently creating a small population model. I have re-created the current model function below from a book, it's a simple population model based on a few parameters. I've left them at default and returned the data frame. Everything works well. However, I was wondering whether I could somehow exclude the loop from the function.

I know R is great because of vectorized calculation, but I'm not sure in this case whether it would be possible. I thought of using something like lead/lag to do it, but would this work? Perhaps not as things need to be calculated sequentially?



# Nt numbers at start of time t
# Ct = removed at the end of time t
# Nt0 = numbers at time 0
# r = intrinsic rate of population growth
# K = carrying capacity


mod_fun = function (r = 0.5, K = 1000, N0 = 50, Ct = 0, Yrs = 10, p = 1) 
{
  # sets years to year value plus 1
  yr1 <- Yrs   1
  # creates sequence of length years from year 1 to Yrs value  !
  years <- seq(1, yr1, 1)
  # uses years length to create a vector of length Yrs   1
  pop <- numeric(yr1)
  # sets population at time 0
  pop[1] <- N0
  
  # creates a loop that calculates  model for each year after first year 
  for (i in 2:yr1) {
    # sets starting value of population for step to one calculated previous step
    # thus Nt is always the previous step pop size
    Nt <- pop[i - 1]
    
    pop[i] <- max((Nt   (r * Nt/p) * (1 - (Nt/K)^p) - 
                     Ct), 0)
  }
  
  # sets pop2 to original pop length
  pop2 <- pop[2:yr1]
  
  # binds together years (sequence from 1 to length Yrs), 
  # pop which is created in loop and is the population at the start of step t
  # pop2 which is the population at the end of step t
  out <- cbind(year = years, nt = pop, nt1 = c(pop2, NA))
  
  # sets row names to 
  rownames(out) <- years
  
  out <- out[-yr1, ]
  
  #returns data.frame
  return(out)
}

  
result = mod_fun()

This is what the output looks like. Basically rowwise starting from row 1 given the starting population of 50 the loop calculates nt1 then sets next nt row to lag(nt1) and then things continue in a similar fashion.

result
#>    year       nt      nt1
#> 1     1  50.0000  73.7500
#> 2     2  73.7500 107.9055
#> 3     3 107.9055 156.0364
#> 4     4 156.0364 221.8809
#> 5     5 221.8809 308.2058
#> 6     6 308.2058 414.8133
#> 7     7 414.8133 536.1849
#> 8     8 536.1849 660.5303
#> 9     9 660.5303 772.6453
#> 10   10 772.6453 860.4776

Created on 2022-04-24 by the reprex package (v2.0.1)

CodePudding user response：

mod_fun = function (r = 0.5, K = 1000, N0 = 50, Ct = 0, Yrs = 10, p = 1) 
{
  years <- seq_len(Yrs)
  pop <- Reduce(function(Nt, y)max((Nt   (r * Nt/p) * (1 - (Nt/K)^p) - Ct), 0),
         years, init = N0, accumulate = TRUE)
  
  data.frame(year = years, nt = head(pop,-1), nt1 = pop[-1])
 
}

   year       nt      nt1
1     1  50.0000  73.7500
2     2  73.7500 107.9055
3     3 107.9055 156.0364
4     4 156.0364 221.8809
5     5 221.8809 308.2058
6     6 308.2058 414.8133
7     7 414.8133 536.1849
8     8 536.1849 660.5303
9     9 660.5303 772.6453
10   10 772.6453 860.4776