Как исправить простую функцию: получение замены ошибки имеет нулевую длину

Question

Я пытаюсь смоделировать процесс брендинга с отрицательным биномиальным распределением потомства. Когда я запускаю один процесс ветвления, код работает нормально. Когда я обертываю его в функцию и использую функцию «replicate» для имитации многих ветвящихся процессов, он выдает ошибку: «замена имеет нулевую длину»

Я - преобразователь SAS, поэтому относительно новый для функций Rдля помощи! Надеюсь, это простое исправление, и любые советы по улучшению всегда приветствуются. Заранее спасибо!

#######
#Single NB Branching Process with 20 generations
n<-20 #20 generations
r0<-0.9 
k<-0.25

#initialize list of population size at generation n
Z<-1
#Initiate with one index case generation 0
Z[0] <- 1
#Cluster size a generation 1
Z[1] <- rnbinom(Z[0], k,r0)
for (i in 2:n)
{
  if(Z[i-1]==0) {Z[i]=0} else 
  {
    x<-rnbinom(Z[i-1], k,r0) 
    Z[i]<- sum(x)
  }
}
print(Z)

######################
#Wrap in a function and replicate 300 times
nbbp<-function(n, r0, k)
{
  #initialize list of population size at generation n
  Z<-1
  #Initiate with one index case generation 0
  Z[0] <- 1
  #Cluster size a generation 1
  Z[1] <- rnbinom(Z[0], k,r0)
  for (i in 2:N)
  {
    if(Z[i-1]==0) {Z[i]=0} else 
    {
      x<-rnbinom(Z[i-1], k,r0) 
      Z[i]<- sum(x)
    }
  }
}
ds<-replicate(100,nbbp(20,0.9,0.25))
#Returns: Error in Z[1] <- rnbinom(Z[0], k, r0) : replacement has length zero

Answer 1

Функция ниже делает то, о чем спрашивает вопрос, и проще.

nbbp <- function(n, r0, k){
  # initialize return vector
  Z <- integer(n)
  #Initiate with 1 index case generation 1
  Z[1] <- 1
  for (i in seq_len(n)[-1]){
    # No need to check the case Z[i - 1] == 0,
    #   all cases are initialized to 0 and
    #   this one hasn't changed yet
    if(Z[i - 1] > 0) {
      x <- rnbinom(Z[i - 1], k, r0) 
      Z[i] <- sum(x)
    }
  }
  Z
}

set.seed(1234)
ds <- replicate(10, nbbp(20, 0.9, 10))
ds

Answer 2

В вашем автономном коде у вас есть

for (i in 2:n)

, а в функции у вас есть

for (i in 2:N)

(обратите внимание на заглавную N). Я представляю, что N определяется как что-то в вашем рабочем пространстве, но не то, что вы хотите. Если вы измените это на n, это работает?