python - 计算 numpy 数组中成对列的组合

Question

我有一个具有一定列数的矩阵，其中仅包含数字0和1，我想计算[0, 0]、[0, 1]、[1, 0]和[1]的数量, 1] 在每对列中。

例如，如果我有一个四列的矩阵，我想计算第一列和第二列中 00、11、01 和 11 的数量，将最终结果附加到列表中，然后循环第三列和第四列并将答案附加到列表中。

输入示例:

array([[0, 1, 1, 0],
       [1, 0, 1, 0],
       [0, 1, 0, 1],
       [0, 0, 1, 1],
       [1, 1, 0, 0]])

我的预期输出是:

array([[1, 1],
       [2, 1],
       [1, 2],
       [1, 1]])

说明:

前两列有一次 [0, 0]。后两列也有一次 [0, 0]。前两列有 [0, 1] 两次，后两列有 [0, 1] 一次...等等。

<小时/>

这是我最近的尝试，它似乎有效。希望得到反馈。

# for each pair of columns calculate haplotype frequencies
# haplotypes:
# h1 = 11
# h2 = 10
# h3 = 01
# h4 = 00
# takes as input a pair of columns
def calc_haplotype_freq(matrix):
    h1_frequencies = []
    h2_frequencies = []
    h3_frequencies = []
    h4_frequencies = []
    colIndex1 = 0
    colIndex2 = 1
    for i in range(0, 2): # number of columns divided by 2
        h1 = 0
        h2 = 0
        h3 = 0
        h4 = 0
        column_1 = matrix[:, colIndex1]
        column_2 = matrix[:, colIndex2]
        for row in range(0, matrix.shape[0]):
            if (column_1[row, 0] == 1).any() & (column_2[row, 0] == 1).any():
                h1 += 1
            elif (column_1[row, 0] == 1).any() & (column_2[row, 0] == 0).any():
                h2 += 1
            elif (column_1[row, 0] == 0).any() & (column_2[row, 0] == 1).any():
                h3 += 1
            elif (column_1[row, 0] == 0).any() & (column_2[row, 0] == 0).any():
                h4 += 1
        colIndex1 += 2
        colIndex2 += 2
        h1_frequencies.append(h1)
        h2_frequencies.append(h2)
        h3_frequencies.append(h3)
        h4_frequencies.append(h4)
    print("H1 Frequencies (11): ", h1_frequencies)
    print("H2 Frequencies (10): ", h2_frequencies)
    print("H3 Frequencies (01): ", h3_frequencies)
    print("H4 Frequencies (00): ", h4_frequencies)

对于上面的示例输入，这给出:

----------
H1 Frequencies (11):  [1, 1]
H2 Frequencies (10):  [1, 2]
H3 Frequencies (01):  [2, 1]
H4 Frequencies (00):  [1, 1]
----------

哪个是正确的，但是有更好的方法吗？如何从函数返回这些结果以进行进一步处理？

Answer 1

从这里开始 -

x
array([[0, 1, 1, 0],
       [1, 0, 1, 0],
       [0, 1, 0, 1],
       [0, 0, 1, 1],
       [1, 1, 0, 0]])

将数组分成 2 列组并将它们连接起来:

y = x.T
z = np.concatenate([y[i:i + 2] for i in range(0, y.shape[0], 2)], 1).T

现在，执行广播比较和求和:

(z[:, None] == [[0, 0], [0, 1], [1, 0], [1, 1]]).all(2).sum(0)
array([2, 3, 3, 2])

<小时/>

如果您想要每列对计数，那么您可以执行以下操作:

def calc_haplotype_freq(x):
    counts = []
    for i in range(0, x.shape[1], 2):
        counts.append(
             (x[:, None, i:i + 2] == [[0, 0], [0, 1], [1, 0], [1, 1]]).all(2).sum(0)
        )

    return np.column_stack(counts)

calc_haplotype_freq(x)
array([[1, 1],
       [2, 1],
       [1, 2],
       [1, 1]])