快速计算大矩阵中的非对角线平均值

我需要计算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平方公尺它花费的时间比我想：

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

是否有这样做的一个更聪明的方式？为了澄清，我想每个对角线的平均值是独立的，例如，为：

我想：

mean(c(1,2,3)) 
mean(c(1,2)) 
mean(1)

来源

2012-12-17 blmoore

你可以得到显著更快只是直接使用线性寻址提取对角线：superdiag这里提取的第i superdiagonal从A（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方阵运行这给了〜800X加速：

> 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

这给你的结果顺序与原文相同。

来源

2012-12-17 14:02:40

很好的使用索引。我投这个票作为公认的答案，因为它说明了指数的强大程度。 –

谢谢你，但是你更清楚@JorisMeys;只有当你必须做一个_lot_和每秒十分之一秒的广告时，这种方法才值得额外复杂化。 –

这非常聪明 - 我必须通过索引生成来了解发生了什么。感谢您的回答 – blmoore

您可以使用以下功能：

diagmean <- function(x){ 
    id <- row(x) - col(x) 
    sol <- tapply(x,id,mean) 
    sol[names(sol)!='0'] 
}

如果我们检查这对您的矩阵，速度增益是巨大的：

> system.time(diagmean(A)) 
    user system elapsed 
    2.58 0.00 2.58 

> system.time(sapply(1:(nrow(A)-1), function(x) mean(A[row(A) == (col(A) - x)]))) 
    user system elapsed 
    38.93 4.01 42.98

请注意，该函数计算上下三角形。你可以使用下面的三角计算例如：

diagmean <- function(A){ 
    id <- row(A) - col(A) 
    id[id>=0] <- NA 
    tapply(A,id,mean) 
}

这会导致另一个速度增益。请注意，该解决方案将比较你的逆转：

> A <- matrix(rep(c(1,2,3,4),4),ncol=4) 

> sapply(1:(nrow(A)-1), function(x) mean(A[row(A) == (col(A) - x)])) 
[1] 2.0 1.5 1.0 

> diagmean(A) 
-3 -2 -1 
1.0 1.5 2.0

来源

2012-12-17 13:44:29

非常好，我的机器上的1k^2矩阵不到1秒。非常感谢 – blmoore

快速计算大矩阵中的非对角线平均值

回答

相关问题