Р: Может ли матрица заполняться в цикле путем индексации строк и столбцов?

Я биолог и новичок в R, и я учусь создавать простую популяционную модель.

Итак, у меня есть популяционная матрица ("поп") из 30 возрастных классов женщин (1:4 - не заводчики, 5:30 - заводчики), которые будут моделироваться в течение 100 лет.

pop <- matrix(0,30,100)

Затем я заполняю эту матрицу 3 молодыми взрослыми женщинами.

pop[5, 1] <- 3

Затем я хочу выполнить это в течение 100 лет со стохастичностью, чтобы увидеть, как это население делает со временем. (Я не заполнил их, но они вам не нужны, все они имеют разные вероятности выживания по сравнению с половозрелыми взрослыми.)

for (t in 1:100)  {  # Edited y to t, typing error!
    pop[1,t+1] <- rbinom(1,colSums(pop[5:30, t]), b/2) 
    pop[5, t+1] <- rbinom(1, pop[4, t], s2)
    pop[6, t+1] <- rbinom(1, pop[5, t], s2)
    .....
    pop[30, t+1] <- rbinom(1, pop[29, t], s2)
}

Итак, мой вопрос: есть ли способ заполнить эту матрицу без необходимости явно писать 30 строк кода? Поскольку строки 5 - 30 все будут одинаковыми, но даже после 5 часов (буквально) веб-поиска и чтения R вручную, я не могу найти способ индексировать строки, что, по-видимому, и здесь необходимо.

Любое понимание здесь приветствуется, в том числе другой способ моделирования этой популяции.

1 ответ

Решение

Поскольку rbinom векторизация, вы должны быть в состоянии заменить

pop[1,t+1] <- rbinom(1,colSums(pop[5:30, t]), b/2) 
pop[5, t+1] <- rbinom(1, pop[4, t], s2)
pop[6, t+1] <- rbinom(1, pop[5, t], s2)
.....
pop[30, t+1] <- rbinom(1, pop[29, t], s2)

с

pop[,t+1] <- rbinom(27,size=c(colSums(pop[5:30,t]),pop[4:29,t]),
                       prob=c(b/2,rep(s2,26)))

или что-то в этом роде (26 и 27 могут быть не правы, я, возможно, не учел)

Тем не менее, я действительно думаю, что вы должны сделать что-то вроде этого:

pop[1,t+1]  <- rpois(1,b/2*sum(pop[5:30,t]))
pop[-1,t+1] <- rbinom(29,size=pop[1:29,t]),
                       prob=c(rep(s1,4),rep(s2,26)))

Это предполагает, что люди в последнем возрастном классе падают с края света и умирают... Я использую rpois здесь, чтобы разрешить воспроизведение Пуассона, а не бином, хотя вы могли бы использовать rbinom(1,size=sum(pop[5:30,t])/2,prob=p) если вы действительно хотите настаивать на том, чтобы у самок было максимум 1 потомство (с вероятностью p)

Другие вопросы по тегам