我需要计算 n × n 矩阵中每个非对角元素的平均值。下三角和上三角是多余的。这是我目前使用的代码:
A <- replicate(500, rnorm(500))
sapply(1:(nrow(A)-1), function(x) mean(A[row(A) == (col(A) - x)]))
这似乎有效,但不适用于较大的矩阵。我所拥有的并不大,大约 2-5000^2,但即使有 1000^2,它也需要比我想要的更长的时间:
A <- replicate(1000, rnorm(1000))
system.time(sapply(1:(nrow(A)-1), function(x) mean(A[row(A) == (col(A) - x)])))
> user system elapsed
> 26.662 4.846 31.494
有没有更聪明的方法来做到这一点?
编辑 为了澄清,我想独立地计算每个对角线的平均值,例如为了:
1 2 3 4
1 2 3 4
1 2 3 4
1 2 3 4
我想:
mean(c(1,2,3))
mean(c(1,2))
mean(1)
最佳答案
只需使用线性寻址直接提取对角线,您就可以显着加快速度:superdiag
这里从 A 中提取第 i 个超对角线(i=1 是主对角线)
superdiag <- function(A,i) {
n<-nrow(A);
len<-n-i+1;
r <- 1:len;
c <- i:n;
indices<-(c-1)*n+r;
A[indices]
}
superdiagmeans <- function(A) {
sapply(2:nrow(A), function(i){mean(superdiag(A,i))})
}
在 1K 方阵上运行它可以得到约 800 倍的加速:
> A <- replicate(1000, rnorm(1000))
> system.time(sapply(1:(nrow(A)-1), function(x) mean(A[row(A) == (col(A) - x)])))
user system elapsed
26.464 3.345 29.793
> system.time(superdiagmeans(A))
user system elapsed
0.033 0.006 0.039
这将为您提供与原始顺序相同的结果。
关于r - 在大矩阵中计算非对角线平均值的更快方法,我们在Stack Overflow上找到一个类似的问题: https://stackoverflow.com/questions/13914845/