R中的代码更快

Question

仅供参考：自从我的第一版以来，我已经对其进行了重大编辑。这种模拟从14小时减少到14分钟。

我是编程新手，但我做了一个模拟试图跟随生物体中的无性复制，并量化父母和子女生物之间染色体数量的差异。模拟运行非常缓慢。大约需要6个小时才能完成。我想知道什么是使模拟运行更快的最佳方法。

这些数字生物具有x个染色体。与大多数生物体不同，染色体彼此独立，因此它们有可能被转移到任何一个女性生物体中。

在这种情况下，染色体在子细胞中的分布遵循二项分布，概率为0.5。

函数sim_repo取一个具有已知染色体数目的数字生物矩阵，并将它们通过12代复制。它复制这些染色体，然后使用rbinom函数随机生成一个数字。然后将该号码分配给女儿单元。由于在无性生殖期间没有染色体丢失，所以其他子细胞接收剩余的染色体。然后这个重复的代数为G。然后从矩阵的每一行中采样一个值。

​​

我我的研究很感兴趣，在最初的祖先生物，在这个模拟的最后时间点的染色体变异中的变化。以下功能代表将细胞转移到新的生活环境中。它采用函数sim_rep的输出并使用它来生成更多代。然后它找出第一个和最后一个矩阵列中的行之间的差异，并找出它们之间的差异。

# The following function is mostly the same as I talked about in the description. 
    # The only difference is I changed some aspects to take into account I am using 
    # matrices and not lists. 
    # The function outputs the difference between the intial variance component between 
    # 'cell lines' with the final variance after t number of transfers 

sim_exp = function(x1, G=12, k=1, t=25, h=1000) { 

    xn <- matrix(NA, nrow(x1), t) 
    x <- x1 
    xn[,1] <- x1 
    for (l in 2:t) { 
     x <- sim_repo(x, G, k, t, h) 
     xn[, l] <- x 
    } 

    colvar <- matrix(apply(xn,2,var),ncol=ncol(xn)) 
    ivar <- colvar[,1] 
    fvar <- colvar[,ncol(xn)] 
    deltavar <- fvar - ivar 
    deltavar 
} 

我需要复制这个模拟^h的次数。因此，我做了以下功能，将调用功能sim_exph次数。

sim_1000 = function(x1, G=12, k=1, t=25, h=1000) { 
    xn <- vector(length=h) 
    for (l in 2:h) { 
     x <- sim_exp(x1, G, k, t, h) 
     xn[l] <- x 
    } 
     xn 
} 

当我打电话与像6倍的值需要花费大约52秒来完成的值sim_exp功能。

x1 <- matrix(data=c(100,100,100,100,100,100),ncol=1) 
system.time(sim_1000(x1,h=1)) 
    user system elapsed 
    1.280 0.105 1.369 

如果我能得到更快然后我可以完成更多的这些模拟和仿真应用选择模型。

我的输入将类似于以下x1，在其自己的行每祖生物基质：

x1 <- matrix(data=c(100,100,100,100,100,100),ncol=1) # a matrix of 6 organisms 

当我运行：

a <- sim_repo(x1, G=12, k=1) 

我的预期产出将是：

a 
    [,1] 
[1,] 137 
[2,] 82 
[3,] 89 
[4,] 135 
[5,] 89 
[6,] 109 

system.time(sim_repo(x1)) 
    user system elapsed 
    1.969 0.059 2.010 

当我调用sim_exp函数时，

b < - sim_exp（X1，G = 12，k = 1时，T = 25）

它调用sim_repo函数G倍和输出：

b 
[1] 18805.47 

当我调用sim_1000功能，我通常会将h设置为1000，但这里我将它设置为2.所以这里sim_1000将调用sim_exp并复制它两次。

c <- sim_1000(x1, G=12, k=1, t=25, h=2) 
c 
[1] 18805.47 18805.47 

Answer 1

正如在评论中提到别人，如果我们在函数sim_repo只能看，并更换行：

dup <- x1 * 2 

dup <- apply(x1, c(1,2),"*",2) 

线

z <- apply(dup,c(1,2),rbinom,n=1,prob=0.5) 

与

z <- matrix(rbinom(n=rep(1,length(dup)),size = as.vector(dup),prob = 0.5),nrow = nrow(dup)) 

和内部的循环与

x1 <- matrix(apply(x1,1,sample,size = 1), ncol=1) 

，我收到了，好了，大的速度增加：

system.time(sim_exp(x1)) 
    user system elapsed 
    0.655 0.017 0.686 
> system.time(sim_expOld(x1)) 
    user system elapsed 
21.445 0.128 21.530 

我验证了它在做同样的事情：

set.seed(123) 
out1 <- sim_exp(x1) 

set.seed(123) 
out2 <- sim_expOld(x1) 

all.equal(out1,out2) 
> TRUE 

而这甚至没有削减进入预先分配，根据您构建代码的方式，如果不完全重新设计，可能实际上很难实现。

而这也甚至没有开始看你是否真的需要所有三种功能...