r - R中的高效矩阵运算

以下代码是我项目的核心，不幸的是，考虑到我的问题的规模，目前它太慢了。有没有更有效的方法来达到同样的结果？

nbassets <- 80
nbrisksource <- 100
nbsimul <- 300000
set.seed(100)
#generate random number for each 100 source of risk in many simulations
random <- matrix(runif(nbsimul*nbrisksource)+0.9,nrow=nbsimul,ncol=nbrisksource)
# random vulnerability to each source of risk for each of 120 assets
EL_decomp <- matrix(runif(nbassets*nbrisksource),nrow=nbassets,ncol=nbrisksource)
#initiate matrix to store asset returns
asset_ret <- matrix(NA, nrow=nbsimul,ncol=nbassets)


ptm <- proc.time()
#loop through each asset
 for (i in 1:nbassets){
  #determine if the asset has been impacted by any source of risk, if yes return is -1, otherwise 0
  asset_ret[,i] <- apply(matrix(EL_decomp[i,], nrow=nbsimul,ncol=nbrisksource,byrow=TRUE) < random,1,all)-1
}
print(proc.time() - ptm)

ptm <- proc.time()

最佳答案

事情可以大大改善。下面是新旧代码对比:

nbassets     <- 80
nbrisksource    <- 100
nbsimul <- 300000
set.seed(100)
random <- matrix(runif(nbsimul*nbrisksource)+0.9, nrow=nbsimul,ncol=nbrisksource)

EL_decomp     <- matrix(runif(nbassets    *nbrisksource), nrow=nbassets,     ncol=nbrisksource)
asset_ret1  <- matrix(NA, nrow=nbsimul, ncol=nbassets)
asset_ret2  <- matrix(NA, nrow=nbsimul, ncol=nbassets)

ptm <- proc.time()
for (i in 1:nbassets){
  #determine if the asset has been impacted by any source of risk, if yes return is -1, otherwise 0
  asset_ret1[,i] <- apply(matrix(EL_decomp[i,],nrow=nbsimul,ncol=nbrisksource,byrow=TRUE) < random,1,all)-1
}
print(head(asset_ret1))
print(proc.time() - ptm) #182s on my old mac

#improved version
ptm <- proc.time()
randomt <- t(random)
asset_ret2 <- apply(EL_decomp, 1, function(x) (colSums(x < randomt) == nbrisksource))- 1L
print(head(asset_ret2))
print(proc.time() - ptm) #14s
print(identical(asset_ret1,asset_ret2))

关于r - R中的高效矩阵运算，我们在Stack Overflow上找到一个类似的问题：https://stackoverflow.com/questions/33837486/