2

我不知道我正在尝试构建的“签名矩阵”是否在任何字段中都有适当的预先存在的名称或定义,但以下代码似乎可以在某些玩具矩阵上生成正确的结果。我很难在不引起混淆的情况下解释我到底想要做什么,但如果我提供的代码不足以推断出我想要做什么,我很乐意试一试。

当我使用我的实际数据(大小均为 300 x 20,000 个元素的两个整数矩阵)运行此代码时,它似乎正在工作,但经过数小时后仍无法完成。

我知道迭代可能是这里最大的问题,但我无法弄清楚如何删除它。

编码:

# Load required library
library(Matrix)

# Load in the test data
mut <- matrix(data=c(1,1,1,0,0,0,0,1,0,1,1,0,0,1,1,0,0,0,1,0),
              nrow=5,ncol=4,
              dimnames=list(c("p1","p2","p3","p4","p5"),c("GA","GB","GC","GD")))

oute <- matrix(data=c(1,1,0,1,0,1,0,0,1,1,1,1,1,0,0,1,1,0,0,1),
              nrow=5,ncol=4,
              dimnames=list(c("p1","p2","p3","p4","p5"),c("GQ","GW","GE","GR")))

patOutMatrix <- Matrix(data=oute,sparse=TRUE)
patMutMatrix <- Matrix(data=mut,sparse=TRUE)

transposePatMutMatrix <- t(patMutMatrix)

# Build the empty matrix (with row and col names)
sigMatrix <- Matrix(0,nrow=ncol(patMutMatrix), ncol=ncol(patOutMatrix),sparse=TRUE)
rownames(sigMatrix) <- colnames(patMutMatrix)
colnames(sigMatrix) <- colnames(patOutMatrix)

# Populate sigMatrix
for (mgene in rownames(transposePatMutMatrix))
{
  a <- patOutMatrix[which(transposePatMutMatrix[mgene, ] == 1, arr.ind = T), ]

  # Using an IF here to get around a problem with colSums() not working on single rows
  sigMatrix[mgene,] <- if (dim(as.matrix(a))[2] == 1) {
    a
  } else {
    colSums(patOutMatrix[which(transposePatMutMatrix[mgene, ] == 1, arr.ind = T), ])
  }
}

有谁知道我如何在这里更改任何内容以使其执行得更快?

4

1 回答 1

0

如何向量化计算

看起来你只是在那里计算一个矩阵产品。所以这样写:

sigMatrix <- t(patMutMatrix) %*% patOutMatrix

或(更难阅读但性能更好):

sigMatrix <- crossprod(patMutMatrix, patOutMatrix)

代码的作用

假设您的矩阵只有 0 和 1 个条目,代码执行以下操作:从第一个矩阵中的一列和第二个矩阵中的一列中获取每个可能的组合,它计算在这两列中都有 1 的行数。对于您的示例数据,这给出了以下结果:

> patMutMatrix   > patOutMatrix   > sigMatrix
   GA GB GC GD      GQ GW GE GR      GQ GW GE GR
p1  1  .  1  .   p1  1  1  1  1   GA  2  1  3  2
p2  1  .  .  .   p2  1  .  1  1   GB  .  1  1  1
p3  1  1  .  .   p3  .  .  1  .   GC  2  3  1  2
p4  .  .  1  1   p4  1  1  .  .   GD  1  1  .  .
p5  .  1  1  .   p5  .  1  .  1

如果矩阵不限于零和一,我的代码将执行与您的代码不同的操作,至于第一个矩阵,您将除 1 之外的任何值视为 0。

其他要学习的东西

避免大小写区分

  # Using an IF here to get around a problem with colSums() not working on single rows

您可以通过传递drop = FALSE给子集操作来避免这种情况,如下所示:

patOutMatrix[which(transposePatMutMatrix[mgene, ] == 1, arr.ind = T), , drop = FALSE]

这样,子集运算的结果将始终是一个矩阵,即使该矩阵只有一行。

避免按名称访问

for (mgene in rownames(transposePatMutMatrix))

通常最好迭代索引,而不是名称,因为按名称选择内容会导致额外的查找将该名称转换回索引。所以我宁愿做这个

for (mgene in 1:nrow(transposePatMutMatrix))
于 2012-10-17T15:30:42.320 回答