R矩阵/数组求和积转换

如何解决R矩阵/数组求和积转换

我正在尝试使用两个矩阵进行矩阵求和运算。 sumproduct在一个矩阵中具有行，而在另一个矩阵中具有列。

set.seed(123)
y <- matrix(sample(1:6,6,FALSE),nrow=3,ncol=2) 
x <- matrix(sample(1:5,15,TRUE),nrow=5,ncol=3)

> x
     [,1] [,2] [,3]
[1,]    5    4    3
[2,]    2    3    2
[3,]    2    3    3
[4,]    4    3    1
[5,]    3    5    5

> y
     [,2]
[1,]    3    4
[2,]    6    5
[3,]    2    1

如果结果为[1,1]的像元值为x [1，]和y [，1]的和。 [1,2]中的像元值是x [1，]和y [，2]的和。 [2,1]中的像元值是x [2，]和y [，1]的和。等等。最终结果应该像这样

> result

         [,2] 
    [1,]   45   43 
    [2,]   28   25 
    [3,]   30   26 
    [4,]   32   32 
    [5,]   49   42

我可以使用循环来执行此操作，但是如果有一个函数可以自动执行此操作，那么会容易得多。

解决方法

我建议采用下一种方法。有一个名为crossprod()的函数，但可用于定义的矩阵维数。您可以检查一下。接下来的代码可能会有用：

set.seed(123)
#Data
y <- matrix(sample(1:6,6,FALSE),nrow=3,ncol=2) 
x <- matrix(sample(1:5,15,TRUE),nrow=5,ncol=3)

矩阵：

x
     [,1] [,2] [,3]
[1,]    5    5    1
[2,]    4    3    1
[3,]    1    3    5
[4,]    2    1    3
[5,]    3    4    2

y
     [,2]
[1,]    3    4
[2,]    6    5
[3,]    2    1

代码：

#Code
z <-   t(sapply(1:nrow(x),function(i){
    x[i,] %*% sapply(1:ncol(y),function(j) {y[,j]})}))

输出：

z
     [,]   47   46
[2,]   32   32
[3,]   31   24
[4,]   18   16
[5,]   37   34

这是使用outer对矩阵的索引进行运算的方法。 Vectorize允许我们在outer内合并结果。如果您需要性能，我真的认为使用rcpp相对容易一些。

set.seed(123)
y <- matrix(sample(1:6,ncol=3)

outer(seq_len(nrow(x)),seq_len(ncol(y)),FUN  = Vectorize(function(i,j) sum(x[i,] * y[,j])))
#>      [,2]
#> [1,]   47   46
#> [2,]   32   32
#> [3,]   31   24
#> [4,]   18   16
#> [5,]   37   34

出于好奇，这花了大约5分钟的时间来提出rcpp解决方案：

#include <Rcpp.h>
using namespace Rcpp;

// [[Rcpp::export]]
IntegerMatrix custom_crossprod(IntegerMatrix x,IntegerMatrix y) {
  int nrow = x.rows();
  int ncol = y.cols();
  IntegerMatrix ans(nrow,ncol);
  
  for (int i = 0; i < nrow; i++) {
    IntegerVector x_tmp = x(i,_);
    for (int j = 0; j < ncol; j++) {
      ans(i,j) = sum(x_tmp * y(_,j));
    }
  }
  return(ans);
}

基准标记：：

bench::mark(
  use_outer = 
outer(seq_len(nrow(x)),j]))),custom_crossprod(x,y)
)
# A tibble: 2 x 13
##  expression                 min  median `itr/sec` mem_alloc `gc/sec`
##  <bch:expr>             <bch:t> <bch:t>     <dbl> <bch:byt>    <dbl>
##1 use_outer              179.8us 191.2us     4394.   11.73KB     4.37
##2 custom_crossprod(x,y)   4.9us   6.6us   135309.    2.49KB     0

@ user20650提供了最合适的答案。

似乎x％*％y足够了

set.seed(123)
y <- matrix(sample(1:6,ncol=3)
x %*% y

> x %*% y
     [,]   37   34

感谢所有答复。您的编码使我学习了很多有关R如何工作以及如何为某些功能编写脚本的知识。谢谢