我有一组DNA序列(字符串(，我以成对的方式相互比较。每次比较都提供了序列之间相似性的确切数量(有多少核苷酸是相同的(，并用于填充较低的对角矩阵。现在，我想在这个矩阵中找到 8 个序列的子集(所有可能的 8 个序列组(，它们之间的相似性最少(这 8 个序列组中的成对相似性应该尽可能低(，但我不知道如何继续......
任何使用R(首选(或Python的帮助将不胜感激！
下面是我的矩阵示例： 这里的主要思想是找到 n 个序列的子集(例如，2 个序列(，这些子集在它们之间共享最少的相似性。 我的原始矩阵是 61X61。
seq1 seq2 seq3 seq4
seq1 NA NA NA NA
seq2 1 NA NA NA
seq3 2 5 NA NA
seq4 3 2 6 NA
在此示例中，相似性最小的 n=2 的子集是 (seq1，seq2(，相似性 = 1。n=3 的子集将是 (seq1，seq2， seq4(，因为在这种情况下，它们的成对相似性之和是最低的(seq1，seq2=1， seq1，seq4=3， seq2，seq4=2; sum = 6(。(我一直使用成对交互作用的最小总和作为目标，但如果无法达到，我很乐意建立一个截止点，例如：子集中的任何成对交互作用都不应> 20(。
不确定我是否完全理解这项任务，我可能过于简单化了，但这是一个尝试。
# some test data
seqs <- matrix(nrow = 10, ncol=10)
x <- length(seqs[lower.tri(seqs)])
seqs[lower.tri(seqs)] <- sample.int(n = 5, size = x, replace = TRUE)
nms <- paste("seq", 1:10, sep="")
rownames(seqs) <- colnames(seqs) <- nms
# all combinations of 4 sequences
all_4 <- combn(x = nms, 4, simplify = FALSE)
names(all_4) <- paste("mat", 1:length(all_4), sep="_")
# a function to subset the matrix to a smaller one
submat <- function(mat, cols) {
mat[cols, cols]
}
mats_4 <- lapply(all_4, function(x) submat(seqs, x))
# similarity per smaller matrix
mats_4_dist <- sapply(mats_4, sum, na.rm=TRUE)
# index of those matrices with similarity < 20
mats_4_lt20_ind <- mats_4_dist < 20
# extract those matrices
mats_4_lt20 <- mats_4[mats_4_lt20_ind]
# alternatively, find the matrices with the minimal sum
mats_4_min <- mats_4[which.min(mats_4_dist)]
编辑：我没有用61x61矩阵和8x8子矩阵测试这种方法。但是我在发布后尝试了它，它肯定会遇到内存问题。即
> combn(61, 8)
Error in matrix(r, nrow = len.r, ncol = count) :
invalid 'ncol' value (too large or NA)
In addition: Warning message:
In combn(61, 8) : NAs introduced by coercion to integer range
这是 python 中的一个实现。请注意，61 选择 8 将达到大约 30 亿，因此检查每个可能的组合，就像我在这里所做的那样，将需要一段时间。
from itertools import combinations
# dataframe stored as df
# assuming columns and indices have same names
subsets_of_columns = combinations(df.columns, 8)
lowest = None
subset = None
for s in subsets_of_columns:
arr = df.loc[s, s].fillna(0).values
if lowest is None:
lowest = arr.sum()
subset = s
else:
if arr.sum() < lowest:
subset = list(s)
lowest = arr.sum()
print(subset, lowest)