我想计算矩阵中每一行的方差。对于以下矩阵A
[,1] [,2] [,3]
[1,] 1 5 9
[2,] 5 6 10
[3,] 50 7 11
[4,] 4 8 12
我想得到
[1] 16.0000 7.0000 564.3333 16.0000
我知道我可以通过apply(A,1,var) 实现这一目标，但有更快或更好的方法吗？从八度开始，我可以用var(A,0,2) 做到这一点，但我不知道如何使用R 中var() 函数的Y 参数。
编辑：典型块的实际数据集大约有 100 行和 500 列。不过，数据总量约为 50GB。
【问题讨论】：
标签： r
您可以使用rowSums 和rowMeans 对行（或列）上的var 进行矢量化处理
RowVar <- function(x, ...) {
rowSums((x - rowMeans(x, ...))^2, ...)/(dim(x)[2] - 1)
}
RowVar(A)
#[1] 16.0000 7.0000 564.3333 16.0000
使用@Richards 数据，产生
microbenchmark(apply(m, 1, var), RowVar(m))
## Unit: milliseconds
## expr min lq median uq max neval
## apply(m, 1, var) 343.369091 400.924652 424.991017 478.097573 746.483601 100
## RowVar(m) 1.766668 1.916543 2.010471 2.412872 4.834471 100
您还可以创建一个更通用的函数，该函数将接收类似于 apply 的语法，但仍将保持矢量化（列方式方差会较慢，因为需要先转置矩阵）
MatVar <- function(x, dim = 1, ...) {
if(dim == 1){
rowSums((x - rowMeans(x, ...))^2, ...)/(dim(x)[2] - 1)
} else if (dim == 2) {
rowSums((t(x) - colMeans(x, ...))^2, ...)/(dim(x)[1] - 1)
} else stop("Please enter valid dimension")
}
MatVar(A, 1)
## [1] 16.0000 7.0000 564.3333 16.0000
MatVar(A, 2)
V1 V2 V3
## 547.333333 1.666667 1.666667
【讨论】：
RowVar(A[, 1:3) 将仅在前 3 列上运行它。另外，看看matrixStats 包，它们几乎涵盖了那里的所有内容。也可以直接从GH下载这些函数，见here
这是apply() 有用的主要原因之一。它旨在对数组或矩阵的边缘进行操作。
set.seed(100)
m <- matrix(sample(1e5L), 1e4L)
library(microbenchmark)
microbenchmark(apply(m, 1, var))
# Unit: milliseconds
# expr min lq median uq max neval
# apply(m, 1, var) 270.3746 283.9009 292.2933 298.1297 343.9531 100
进行 10,000 次计算需要 300 毫秒吗？
【讨论】：