在R中生成特殊矩阵的更有效函数

白胡子

我正在尝试对R函数进行编程,以创建一个称为Ui的特殊矩阵(第5页):http : //joshuachan.org/papers/Chan-Jeliazkov-2009.pdf

到目前为止,我已经开发了此功能,我怀疑可以对其进行改进(即提高效率):

create_Uj <- function(uj) {
  q <- length(uj) 
  if (q == 1) return(0)
  nr <- q
  nc <- q*(q-1)/2
  Uj <- matrix(0, nrow = nr, ncol = nc)
  for (kk in 2:nr) {
    uj_sub <- uj[1:(kk-1)]
    Uj[kk, 1:(kk*(kk-1)/2)] <- c(rep(0, (kk-1)*(kk-2)/2), uj[1:(kk-1)])
  }
  -Uj
}

create_Uj(1:4)

     [,1] [,2] [,3] [,4] [,5] [,6]
[1,]    0    0    0    0    0    0
[2,]   -1    0    0    0    0    0
[3,]    0   -1   -2    0    0    0
[4,]    0    0    0   -1   -2   -3

有没有更好的方法对此进行编码?

注意:我在论文中使用下标j代替i

冰激凌巨嘴鸟

您只需避免创建即可获得不错的加速效果rep(0, (kk-1)*(kk-2)/2)删除此步骤可以大大加快操作速度,这一事实使我认为如果不使用,则不可能更快Rcpp

create_Uj2 <- function(uj) {
  q <- length(uj) 
  if (q == 1) return(0)
  nr <- q
  nc <- q*(q-1)/2
  Uj <- matrix(0, nrow = nr, ncol = nc)
  for (kk in 2:nr) {
    Uj[kk, ((kk-1)*(kk-2)/2 + 1):(kk*(kk-1)/2)] <- uj[1:(kk-1)]
  }
  -Uj
}

all.equal(create_Uj2(1:400), create_Uj(1:400))
# [1] TRUE
microbenchmark(
  create_Uj2(1:400),
  create_Uj(1:400),
  times = 10,
  unit = 'relative'
)
# Unit: relative
#               expr      min       lq     mean   median       uq       max neval
#  create_Uj2(1:400) 1.000000 1.000000 1.000000 1.000000 1.000000 1.0000000    10
#   create_Uj(1:400) 2.070708 1.847242 1.529658 2.028489 1.496275 0.9441935    10

本文收集自互联网,转载请注明来源。

如有侵权,请联系 [email protected] 删除。

编辑于
0

我来说两句

0 条评论
登录 后参与评论

相关文章