创建一个 8x8 矩阵,其中元素由二进制字符串标记,条目遵循某种概率分布

Creating a 8x8 matrix with elements labeled by binary strings and entries following some probability distribution

我想创建一个 8x8 矩阵,它提供位通信中的错误概率。矩阵如下所示:

列代表观察到的数量,行代表测量数量。元素 p[i,j] 等于条件概率 p(j|i)。例如,元素 p[0,1] 给出当实际值为 000 时观察字符串 001 的概率,即它测量 p(001|000).

问题:如何在Python中创建这样的矩阵,使得

  1. 翻转次数越多,等效条件概率越小(例如p(100|000)<p(110|000)?
  2. 如何启用“不对称”。即 p(001|000)< p(000|001) 的概率。也就是说,与从 0 到 1 的转换相比,从 1 到 0 的转换概率更高。

当然,每一行的概率之和必须等于1。

总而言之,我想在 Python 中创建一个函数,该函数将整数 n 作为输入(矩阵的大小,或者等效地,其中 2^n 是位串)并输出具有上述指定规则的概率转移矩阵。

难点在于如何实现概率分布来填充单元格。

创建 8x8 数组并填充对角线很简单:

P = np.zeros((8,8))
for i in range(8):
    for j in range(8):
        if i==j:
            P[i,j]=1

同样,用固定数字填充给定行或给定列也是微不足道的。但是,我无法弄清楚(甚至如何开始)按照上述规则填充这样的矩阵,甚至无法准确定义元素必须遵循的分布。

概率论

您需要查看 Bernoulli or Binomial 发行版。矩阵中的每一行都可以用其中之一建模。

要利用 Python 中的这些分布,您可以使用 scipy。特别是,请参阅 scipy 文档的 bernoulli and binom in the Discrete Distributions 部分。

这里有一个小示例,可以帮助您开始。将矩阵中的每个单元格视为两种可能性的乘积(或更准确地说,convolution):n 位从 0 翻转到 1 的概率,以及 m 位的概率从 1 翻转到 0(因为你想要不对称)。然后,您可以使用二项分布对这两种情况进行建模,并使用 pmf 方法提取每个结果的概率,然后将它们相乘:

from scipy.stats import binom

# parameters
n = 3  # number of bits
p_0_to_1 = 0.01  # a 0 flips to a 1 1% of the time
p_1_to_0 = 0.03  # a 1 flips to a zero 3% of the time

# the distributions for both scenarios
dist_0_to_1 = binom(n, p_0_to_1)
dist_1_to_0 = binom(n, p_1_to_0)

# print the results
for i in range(n):
    print(f"Probability of {i} bit flips (0->1): {dist_0_to_1.pmf(i)}")
for i in range(n):
    print(f"Probability of {i} bit flips (1->0): {dist_1_to_0.pmf(i)}")
Probability of 0 bit flips (0->1): 0.7290000000000001
Probability of 1 bit flips (0->1): 0.243
Probability of 2 bit flips (0->1): 0.027

Probability of 0 bit flips (1->0): 0.3429999999999999
Probability of 1 bit flips (1->0): 0.4409999999999998
Probability of 2 bit flips (1->0): 0.18899999999999992

注意每个分布的概率加起来为 1。

申请您的案例

对于您的示例,上面矩阵第二行中的红色单元格都具有相同的发生概率,即 P(one bit flips from 0 to 1, and no bits flip 1 to 0)。因此,该特定场景的概率为 P(one 0->1 flip) * P(zero 1->0 flips):

# show the result of a particular scenario
p_single_0_to_1_flip = dist_0_to_1.pmf(1) * dist_1_to_0.pmf(0)
print("Probability of single 0->1 bit flip: ", p_single_0_to_1_flip)
Probability of single 0->1 bit flip:  0.026835324218999996

以这种方式,您可以通过计算两者之间的0->11->0转换次数,然后遵循相同的过程来计算将任意位串转换为另一个位串的概率。

最简单的方法是首先制作另一个矩阵,它显示了所有的可能性。我用 pandas 让它打印得更漂亮:

# generate a matrix of all possible scenarios
probabilities = np.zeros((n + 1, n + 1))
for i in range(n + 1):
    for j in range(n + 1):
        probabilities[j, i] = dist_0_to_1.pmf(i) * dist_1_to_0.pmf(j)
print(pd.DataFrame(probabilities).round(6))
          0         1         2         3
0  0.885566  0.026835  0.000271  0.000001
1  0.082166  0.002490  0.000025  0.000000
2  0.002541  0.000077  0.000001  0.000000
3  0.000026  0.000001  0.000000  0.000000

此矩阵包含您需要的所有数据,并且可以轻松转换为您想要的 table。列标签表示 0->1 转换次数,行标签表示 1->0 转换次数,单元格值显示每个组合场景的概率。

作为完整性检查,您可以看到所有可能情况的概率为 1:

# see if it all sums to 1 as it should 
print(np.sum(probabilities))
1.0

现在添加一个可以 return 转换计数的辅助函数,并尝试一下:

# make a helper function to count the transitions
def count_transitions(s1, s2):
    """Return the number of 0->1 and 1->0 transitions between two binary
    strings, in that order.
    """
    cnts = [0, 0]
    for i, val in enumerate(s1):
        if s1[i] != s2[i]:
            cnts[int(val)] += 1
    return cnts


# test out the helper function on a test case
i, j = count_transitions("000", "001")
print(probabilities[j, i])
0.026835324218999996

请注意,这与我们上面得到的相同。

综合起来

我在这里使用 pandas,因为我可以像使用 f-string:

一样命名列
# put it all together
labels = [f"{i:03b}" for i in range(0, 2**n)]
result = pd.DataFrame(index=labels, columns=labels)
for rowIndex, row in result.iterrows():
    for columnIndex, value in row.items():
        i, j = count_transitions(rowIndex, columnIndex)
        result.at[rowIndex, columnIndex] = probabilities[j, i]
print(result)
          000       001       010       011       100       101       110       111
000  0.885566  0.026835  0.026835  0.000271  0.026835  0.000271  0.000271  0.000001
001  0.082166  0.885566   0.00249  0.026835   0.00249  0.026835  0.000025  0.000271
010  0.082166   0.00249  0.885566  0.026835   0.00249  0.000025  0.026835  0.000271
011  0.002541  0.082166  0.082166  0.885566  0.000077   0.00249   0.00249  0.026835
100  0.082166   0.00249   0.00249  0.000025  0.885566  0.026835  0.026835  0.000271
101  0.002541  0.082166  0.000077   0.00249  0.082166  0.885566   0.00249  0.026835
110  0.002541  0.000077  0.082166   0.00249  0.082166   0.00249  0.885566  0.026835
111  0.000026  0.002541  0.002541  0.082166  0.002541  0.082166  0.082166  0.885566

一些不错的完整性检查:

  • 所有行和列加起来为 1
  • 没有错误的概率(对角线)是一致的
  • 各种错误的概率也差不多

一个包装精美的函数,没有所有解释性的绒毛

import numpy as np
import pandas as pd
from scipy.stats import binom


def generate_error_matrix(n, p_0_to_1, p_1_to_0):
    """Print a matrix of the probabilities of all permutations of bit
    errors for a given bit length.

    Parameters:
        n - the number of bits
        p_0_to_1 - probability of a bit flipping from 0 to 1
        p_1_to_0 - probability of a bit flipping from 1 to 0
    """

    def count_transitions(s1, s2):
        """Return the number of 0->1 and 1->0 transitions between two 
        binary strings, in that order.
        """
        cnts = [0, 0]
        for i, val in enumerate(s1):
            if s1[i] != s2[i]:
                cnts[int(val)] += 1
        return cnts

    dist_0_to_1 = binom(n, p_0_to_1)
    dist_1_to_0 = binom(n, p_1_to_0)

    probabilities = np.zeros((n + 1, n + 1))
    for i in range(n + 1):
        for j in range(n + 1):
            probabilities[j, i] = dist_0_to_1.pmf(i) * dist_1_to_0.pmf(j)

    labels = [f"{i:03b}" for i in range(0, 2**n)]
    result = pd.DataFrame(index=labels, columns=labels)
    for rowIndex, row in result.iterrows():
        for columnIndex, _ in row.items():
            i, j = count_transitions(rowIndex, columnIndex)
            result.at[rowIndex, columnIndex] = probabilities[j, i]
    return result

还有一些例子:

# A symmetric case 
print(generate_error_matrix(n=3, p_0_to_1=0.01, p_1_to_0=0.01))
          000       001       010       011       100       101       110       111
000   0.94148   0.02853   0.02853  0.000288   0.02853  0.000288  0.000288  0.000001
001   0.02853   0.94148  0.000865   0.02853  0.000865   0.02853  0.000009  0.000288
010   0.02853  0.000865   0.94148   0.02853  0.000865  0.000009   0.02853  0.000288
011  0.000288   0.02853   0.02853   0.94148  0.000009  0.000865  0.000865   0.02853
100   0.02853  0.000865  0.000865  0.000009   0.94148   0.02853   0.02853  0.000288
101  0.000288   0.02853  0.000009  0.000865   0.02853   0.94148  0.000865   0.02853
110  0.000288  0.000009   0.02853  0.000865   0.02853  0.000865   0.94148   0.02853
111  0.000001  0.000288  0.000288   0.02853  0.000288   0.02853   0.02853   0.94148
# A case where zeroes can't flip to ones 
print(generate_error_matrix(n=3, p_0_to_1=0, p_1_to_0=0.05))
          000       001       010       011       100       101       110       111
000  0.857375       0.0       0.0       0.0       0.0       0.0       0.0       0.0
001  0.135375  0.857375       0.0       0.0       0.0       0.0       0.0       0.0
010  0.135375       0.0  0.857375       0.0       0.0       0.0       0.0       0.0
011  0.007125  0.135375  0.135375  0.857375       0.0       0.0       0.0       0.0
100  0.135375       0.0       0.0       0.0  0.857375       0.0       0.0       0.0
101  0.007125  0.135375       0.0       0.0  0.135375  0.857375       0.0       0.0
110  0.007125       0.0  0.135375       0.0  0.135375       0.0  0.857375       0.0
111  0.000125  0.007125  0.007125  0.135375  0.007125  0.135375  0.135375  0.857375

最后一个案例表明某处仍然存在错误,因为第一行应该有 100% 的成功率,但我相信这会给你所需要的。

如果位转换的概率只依赖于原始位值,而与位置无关(即P(xy|ab) == P(yx|ba),那么你可以简单地block-multiply一个转换概率核:

x 是一个 2x2 矩阵,使得 x[i,j] 是在给定真值 i 的情况下观察位 j 的概率。即:

x = [[a, b]
     [c, d]]

2位概率矩阵为:

x2 = [[a, a, b, b],          [[a, b, a, b],
      [a, a, b, b],    *      [c, d, c, d],
      [c, c, d, d],           [a, b, a, b],
      [c, c, d, d]]           [c, d, c, d]]

这样的block-multiplication可以简单的表示为numpy:

def bmul(a, x):
    n = a.shape[0] * x.shape[0]
    return (a[:, None, :, None] * x[None, :, None, :]).reshape(n, n)

示例:

u = .2  # "up": p(1|0)
d = .1  # "down": p(0|1)
x = np.array([[1-u, u], [d, 1-d]])

>>> x
array([[0.8, 0.2],
       [0.1, 0.9]])

x2 = bmul(x, x)
>>> x2
array([[0.64, 0.16, 0.16, 0.04],
       [0.08, 0.72, 0.02, 0.18],
       [0.08, 0.02, 0.72, 0.18],
       [0.01, 0.09, 0.09, 0.81]])

x3 = bmul(x2, x)
>>> x3
array([[0.512, 0.128, 0.128, 0.032, 0.128, 0.032, 0.032, 0.008],
       [0.064, 0.576, 0.016, 0.144, 0.016, 0.144, 0.004, 0.036],
       [0.064, 0.016, 0.576, 0.144, 0.016, 0.004, 0.144, 0.036],
       [0.008, 0.072, 0.072, 0.648, 0.002, 0.018, 0.018, 0.162],
       [0.064, 0.016, 0.016, 0.004, 0.576, 0.144, 0.144, 0.036],
       [0.008, 0.072, 0.002, 0.018, 0.072, 0.648, 0.018, 0.162],
       [0.008, 0.002, 0.072, 0.018, 0.072, 0.018, 0.648, 0.162],
       [0.001, 0.009, 0.009, 0.081, 0.009, 0.081, 0.081, 0.729]])

最后一个值就是您要查找的矩阵。

随机检查:

# P(100|010) is u*d*(1-u), and we should find it in x3[4,2]
>>> u * d * (1-u)
0.016000000000000004

>>> x3[4,2]
0.016000000000000004

有趣的事实:

bmul 是结合性的但不是交换性的。换句话说:

  • bmul(bmul(a, b), c) == bmul(a, bmul(b, c),但是
  • bmul(a, b) != bmul(b, a)