R中对角线之间矩阵元素的有效总和
Sum of matrix elements between diagonals efficiently in R
我有 n*n
矩阵形式的数据,我想对其进行一些计算(例如 sum
),其元素位于对角线之间(不包括对角线)。
例如这个矩阵:
[,1] [,2] [,3] [,4] [,5]
[1,] 2 0 1 4 3
[2,] 5 3 6 0 4
[3,] 3 5 2 3 1
[4,] 2 1 5 3 2
[5,] 1 4 3 4 1
sum
(对角元素之间)的结果为:
# left slice 5+3+2+5 = 15
# bottom slice 4+3+4+5 = 16
# right slice 4+1+2+3 = 10
# top slice 0+1+4+6 = 11
# dput(m)
m <- structure(c(2, 5, 3, 2, 1, 0, 3, 5, 1, 4, 1, 6, 2, 5, 3, 4, 0,
3, 3, 4, 3, 4, 1, 2, 1), .Dim = c(5L, 5L))
如何有效地做到这一点?
以下是获取 "top slice" 的方法:
sum(m[lower.tri(m)[nrow(m):1,] & upper.tri(m)])
#[1] 11
形象化:
lower.tri(m)[nrow(m):1,] & upper.tri(m)
# [,1] [,2] [,3] [,4] [,5]
#[1,] FALSE TRUE TRUE TRUE FALSE
#[2,] FALSE FALSE TRUE FALSE FALSE
#[3,] FALSE FALSE FALSE FALSE FALSE
#[4,] FALSE FALSE FALSE FALSE FALSE
#[5,] FALSE FALSE FALSE FALSE FALSE
以下是计算所有 4 个切片的方法:
up <- upper.tri(m)
lo <- lower.tri(m)
n <- nrow(m)
# top
sum(m[lo[n:1,] & up])
# left
sum(m[lo[n:1,] & lo])
# right
sum(m[up[n:1,] & up])
# bottom
sum(m[up[n:1,] & lo])
sum(sapply(1:dim(m)[[2L]], function(i) sum(m[c(-i,-(dim(m)[[1L]]-i+1)),i])))
这是逐列进行的,对于每一列,取出对角线元素并对其余元素求和。然后总结这些部分结果。
我相信这会很快,因为我们逐列进行,R 中的矩阵是逐列存储的(即 CPU 缓存友好)。我们也不必生成大的索引向量,只需为每列生成两个索引(取出的索引)的向量。
编辑:我再次仔细阅读了这个问题。可以更新代码,为 sapply: 中的每个区域的每个元素生成列表四个值。这个想法保持不变,对于大矩阵,如果你逐列而不是在列之间来回跳转,它会很快。
我有 n*n
矩阵形式的数据,我想对其进行一些计算(例如 sum
),其元素位于对角线之间(不包括对角线)。
例如这个矩阵:
[,1] [,2] [,3] [,4] [,5]
[1,] 2 0 1 4 3
[2,] 5 3 6 0 4
[3,] 3 5 2 3 1
[4,] 2 1 5 3 2
[5,] 1 4 3 4 1
sum
(对角元素之间)的结果为:
# left slice 5+3+2+5 = 15
# bottom slice 4+3+4+5 = 16
# right slice 4+1+2+3 = 10
# top slice 0+1+4+6 = 11
# dput(m)
m <- structure(c(2, 5, 3, 2, 1, 0, 3, 5, 1, 4, 1, 6, 2, 5, 3, 4, 0,
3, 3, 4, 3, 4, 1, 2, 1), .Dim = c(5L, 5L))
如何有效地做到这一点?
以下是获取 "top slice" 的方法:
sum(m[lower.tri(m)[nrow(m):1,] & upper.tri(m)])
#[1] 11
形象化:
lower.tri(m)[nrow(m):1,] & upper.tri(m)
# [,1] [,2] [,3] [,4] [,5]
#[1,] FALSE TRUE TRUE TRUE FALSE
#[2,] FALSE FALSE TRUE FALSE FALSE
#[3,] FALSE FALSE FALSE FALSE FALSE
#[4,] FALSE FALSE FALSE FALSE FALSE
#[5,] FALSE FALSE FALSE FALSE FALSE
以下是计算所有 4 个切片的方法:
up <- upper.tri(m)
lo <- lower.tri(m)
n <- nrow(m)
# top
sum(m[lo[n:1,] & up])
# left
sum(m[lo[n:1,] & lo])
# right
sum(m[up[n:1,] & up])
# bottom
sum(m[up[n:1,] & lo])
sum(sapply(1:dim(m)[[2L]], function(i) sum(m[c(-i,-(dim(m)[[1L]]-i+1)),i])))
这是逐列进行的,对于每一列,取出对角线元素并对其余元素求和。然后总结这些部分结果。
我相信这会很快,因为我们逐列进行,R 中的矩阵是逐列存储的(即 CPU 缓存友好)。我们也不必生成大的索引向量,只需为每列生成两个索引(取出的索引)的向量。
编辑:我再次仔细阅读了这个问题。可以更新代码,为 sapply: 中的每个区域的每个元素生成列表四个值。这个想法保持不变,对于大矩阵,如果你逐列而不是在列之间来回跳转,它会很快。