如何将矩阵的每一行除以R中的向量元素
Posted
技术标签:
【中文标题】如何将矩阵的每一行除以R中的向量元素【英文标题】：How to divide each row of a matrix by elements of a vector in R 【发布时间】：2022-01-12 05:51:23 【问题描述】：我想将矩阵的每一行除以一个固定向量。例如
mat<-matrix(1,ncol=2,nrow=2,TRUE)
dev<-c(5,10)
给mat/dev 将每个列除以dev。
[,1] [,2]
[1,] 0.2 0.2
[2,] 0.1 0.1
但是，我希望得到这个结果，即执行 row-wise 操作：
rbind(mat[1,]/dev, mat[2,]/dev)
[,1] [,2]
[1,] 0.2 0.1
[2,] 0.2 0.1
是否有明确的命令可以到达那里？
【问题讨论】：
请务必注意，mat/dev 只会按照您在length(dev) == nrow(mat) 中显示的那样划分每一列。这是因为 R 将其矩阵信息存储在 column major order 中。
将向量称为vec 会更清楚，就像矩阵是mat 一样，但现在为时已晚。
【参考方案1】：
不仅比最快的基本方法更快（使用矩阵乘法，请参阅@G. Grothendieck's answer），而且更具可读性的是使用rray package，它允许使用其%b/% 运算符进行类似numpy 的广播：
mat %b/% matrix(dev, nrow = 1)
#install.packages("rray")
library(rray)
set.seed(84789)
mat <- matrix(runif(1e6), nrow = 1e5)
dev <- runif(10)
bench::mark(rray = ret <- mat %b/% matrix(dev, nrow = 1); attr(ret, "dimnames") <- NULL; ret,
matmult = mat %*% diag(1/dev)
)
#> # A tibble: 2 x 6
#> expression min median `itr/sec` mem_alloc `gc/sec`
#> <bch:expr> <bch:tm> <bch:tm> <dbl> <bch:byt> <dbl>
#> 1 rray 8.24ms 8.82ms 108. 7.67MB 46.4
#> 2 matmult 11.17ms 12.01ms 77.6 7.66MB 13.9
【讨论】：
【参考方案2】：这里有几种按代码长度递增的顺序：
t(t(mat) / dev)
mat / dev[col(mat)] # @DavidArenburg & @akrun
mat %*% diag(1 / dev)
sweep(mat, 2, dev, "/")
t(apply(mat, 1, "/", dev))
plyr::aaply(mat, 1, "/", dev)
mat / rep(dev, each = nrow(mat))
mat / t(replace(t(mat), TRUE, dev))
mapply("/", as.data.frame(mat), dev) # added later
mat / matrix(dev, nrow(mat), ncol(mat), byrow = TRUE) # added later
do.call(rbind, lapply(as.data.frame(t(mat)), "/", dev))
mat2 <- mat; for(i in seq_len(nrow(mat2))) mat2[i, ] <- mat2[i, ] / dev
数据框
如果mat 是数据框并产生数据框结果，则所有以mat / 开头的解决方案也适用。 sweep 解决方案和最后一个，即mat2 解决方案也是如此。 mapply 解决方案适用于 data.frames，但会生成一个矩阵。
矢量
如果mat 是纯向量而不是矩阵，则其中任何一个都返回单列矩阵
t(t(mat) / dev)
mat / t(replace(t(mat), TRUE, dev))
这个返回一个向量：
plyr::aaply(mat, 1, "/", dev)
其他人给出错误、警告或不是想要的答案。
基准
代码的简洁明了可能比速度更重要，但为了完整起见，这里有一些使用 10 次重复和 100 次重复的基准。
library(microbenchmark)
library(plyr)
set.seed(84789)
mat<-matrix(runif(1e6),nrow=1e5)
dev<-runif(10)
microbenchmark(times=10L,
"1" = t(t(mat) / dev),
"2" = mat %*% diag(1/dev),
"3" = sweep(mat, 2, dev, "/"),
"4" = t(apply(mat, 1, "/", dev)),
"5" = mat / rep(dev, each = nrow(mat)),
"6" = mat / t(replace(t(mat), TRUE, dev)),
"7" = aaply(mat, 1, "/", dev),
"8" = do.call(rbind, lapply(as.data.frame(t(mat)), "/", dev)),
"9" = mat2 <- mat; for(i in seq_len(nrow(mat2))) mat2[i, ] <- mat2[i, ] / dev,
"10" = mat/dev[col(mat)])
给予：
Unit: milliseconds
expr min lq mean median uq max neval
1 7.957253 8.136799 44.13317 8.370418 8.597972 366.24246 10
2 4.678240 4.693771 10.11320 4.708153 4.720309 58.79537 10
3 15.594488 15.691104 16.38740 15.843637 16.559956 19.98246 10
4 96.616547 104.743737 124.94650 117.272493 134.852009 177.96882 10
5 17.631848 17.654821 18.98646 18.295586 20.120382 21.30338 10
6 19.097557 19.365944 27.78814 20.126037 43.322090 48.76881 10
7 8279.428898 8496.131747 8631.02530 8644.798642 8741.748155 9194.66980 10
8 509.528218 524.251103 570.81573 545.627522 568.929481 821.17562 10
9 161.240680 177.282664 188.30452 186.235811 193.250346 242.45495 10
10 7.713448 7.815545 11.86550 7.965811 8.807754 45.87518 10
对所有耗时
microbenchmark(times=100L,
"1" = t(t(mat) / dev),
"2" = mat %*% diag(1/dev),
"3" = sweep(mat, 2, dev, "/"),
"5" = mat / rep(dev, each = nrow(mat)),
"6" = mat / t(replace(t(mat), TRUE, dev)),
"10" = mat/dev[col(mat)])
给予：
Unit: milliseconds
expr min lq mean median uq max neval
1 8.010749 8.188459 13.972445 8.560578 10.197650 299.80328 100
2 4.672902 4.734321 5.802965 4.769501 4.985402 20.89999 100
3 15.224121 15.428518 18.707554 15.836116 17.064866 42.54882 100
5 17.625347 17.678850 21.464804 17.847698 18.209404 303.27342 100
6 19.158946 19.361413 22.907115 19.772479 21.142961 38.77585 100
10 7.754911 7.939305 9.971388 8.010871 8.324860 25.65829 100
所以在这两个测试中#2（使用diag）是最快的。原因可能在于它几乎直接吸引了 BLAS，而 #1 则依赖于更昂贵的 t。
【讨论】：
我希望前两个选项中的一个最快。 不是最快但非常明确：scale(mat, center = FALSE, scale = dev)
@flodel，注意scale 在内部使用sweep。
@tomka，它为问题中的示例提供了正确答案，但我认为其意图是 mat 可以是具有任意元素的通用矩阵，在这种情况下，它给出了错误的答案，一般。
@Alnair。实际上它确实适用于一列矩阵，但在您的代码中 mat[, -1] 不是一列矩阵。如果您想捕捉这种边缘情况，请使用mat[, -1, drop = FALSE]。请参阅 R 常见问题解答 7.5 -- cran.r-project.org/doc/FAQ/…【参考方案3】：
您正在寻找应用于行的apply 函数：
t(apply(mat, 1, function(x) x/dev))
【讨论】：
谢谢——好的，对于如此简单的操作来说，这似乎相当复杂。这是最简单/最短/最简洁的方式吗？ @tomka:m.FUN.m..v <- function(FUN) function(m, v) t(FUN(t(m), v)); '%m/v%' <- m.FUN.m..v('/'); M9 <- matrix(1:9, ncol=3); M9 %m/v% 1:3 足够简短吗？然后你几乎可以免费做m.FUN.m..v('+')等。 （感谢@g-grothendieck 和/或@anton 提供双重转置提示（+1）；对引号感到抱歉——如何标记作为代码一部分的反引号？）以上是关于如何将矩阵的每一行除以R中的向量元素的主要内容，如果未能解决你的问题，请参考以下文章