如何为循环的每次迭代随机化矩阵元素？

Question

library(popbio)

babies <- 0.3 
kids <- 0.5 
teens <- 0.75
adults <- 0.98

A <- c(0,0,0,0,teens*0.5,adults*0.8,
         babies,0,0,0,0,0,
         0,kids,0,0,0,0,
         0,0,kids,0,0,0,
         0,0,0,teens,0,0,
         0,0,0,0,teens,adults
)
A <- matrix ((A), ncol=6, byrow = TRUE)

N<-c(10,10,10,10,10,10)
N<-matrix (N, ncol=1)

model <- pop.projection(A,N,iterations=10)
model

更新：这与运行10个不同初始随机值的独立模型是不同的。我希望更新发生在单个模型运行中，该模型本身在pop.projection函数中有10次迭代。

希望你能帮忙。

Answer 1

我知道这个答案已经很晚了，但这里有一个使用表达式的方法。首先，使用表达式创建矩阵。

vr <- list( babies=0.3, kids=0.5, teens=0.75, adults=0.98 )
Ax <- expression( matrix(c(
         0,0,0,0,teens*0.5,adults*0.8,
         babies,0,0,0,0,0,
         0,kids,0,0,0,0,
         0,0,kids,0,0,0,
         0,0,0,teens,0,0,
         0,0,0,0,teens,adults), ncol=6, byrow = TRUE ))
A1 <- eval(Ax, vr)
lambda(A1)
[1] 1.011821

接下来，使用表达式使用nrorm或其他函数创建重要速率。

vr2 <- expression(  list( babies=rnorm(1,0.3,0.1), kids=0.5, teens=0.75, adults=0.98 ))
A2 <- eval(Ax, eval( vr2))
lambda(A2)
[1] 1.014586

将表达式应用于100个矩阵。

x <- sapply(1:100, function(x) lambda(eval(Ax, eval(vr2))))
quantile(x, c(.05,.95))
5%      95% 
0.996523 1.025900 

pop.projection2 <- function (Ax, vr, n, iterations = 20) 
{
    x <- length(n)
    t <- iterations
    stage <- matrix(numeric(x * t), nrow = x)
    pop <- numeric(t)
    change <- numeric(t - 1)
    for (i in 1:t) {
        stage[, i] <- n
        pop[i] <- sum(n)
        if (i > 1) {
            change[i - 1] <- pop[i]/pop[i - 1]
        }
        ## evaluate Ax
        A <- eval(Ax, eval(vr))
        n <- A %*% n
    }
    colnames(stage) <- 0:(t - 1)
    w <- stage[, t]
    pop.proj <- list(lambda = pop[t]/pop[t - 1], stable.stage = w/sum(w), 
        stage.vectors = stage, pop.sizes = pop, pop.changes = change)
    pop.proj
}

n <-c(10,10,10,10,10,10)
pop.projection2(Ax, vr2, n, 10)
$lambda
[1] 0.9874586

$stable.stage
[1] 0.33673579 0.11242588 0.08552367 0.02189786 0.02086656 0.42255023

$stage.vectors
      0         1         2         3          4          5          6          7          8         9
[1,] 10 11.590000 16.375700 19.108186 20.2560223 20.5559445 20.5506251 20.5898222 20.7603581 20.713271
[2,] 10  4.147274  3.332772  4.443311  5.6693931  1.9018887  6.8455597  5.3879202 10.5214540  6.915534
[3,] 10  5.000000  2.073637  1.666386  2.2216556  2.8346965  0.9509443  3.4227799  2.6939601  5.260727
[4,] 10  5.000000  2.500000  1.036819  0.8331931  1.1108278  1.4173483  0.4754722  1.7113899  1.346980
[5,] 10  7.500000  3.750000  1.875000  0.7776139  0.6248948  0.8331209  1.0630112  0.3566041  1.283542
[6,] 10 17.300000 22.579000 24.939920 25.8473716 25.9136346 25.8640330 25.9715930 26.2494195 25.991884

$pop.sizes
 [1] 60.00000 50.53727 50.61111 53.06962 55.60525 52.94189 56.46163 56.91060 62.29319 61.51194

$pop.changes
[1] 0.8422879 1.0014610 1.0485765 1.0477793 0.9521023 1.0664832 1.0079517 1.0945797 0.9874586