LDPC码及其译码算法（全网脉络最清晰）_ldpc译码 spa

  在编码理论第十一章：LDPC码原理学习及其习题分析中当时为便于复习记了个简易的blog。另外起一篇长文，详细谈谈LDPC和SPA算法，以及SPA算法后续的一些改进。

面向人群：具有编码理论背景但通原知识不牢而且对LDPC码知识了解不多的人，大佬们刷到时建议跳过，新进展google学术。

1. LDPC码及其译码算法

参考文献：
[1] Ryan, William. (2004). An introduction to LDPC codes. 10.1201/9780203490310.ch36.
[2] Shokrollahi, Amin. (2004). LDPC Codes: an Introduction. 10.1007/978-3-0348-7865-4_5.
[3] Gallager, Robert G. “Principles of Digital Communication”.

1.1 LDPC码和SPA译码想要做什么

  LDPC码由Gallager于1960年首次提出，由Tanner在1981年推广了LDPC码并引入了称为Tanner图的图形表示。LDPC想要解决的问题是：高速数据无线传输时的误码和噪声问题，其设计初衷在不同的信道条件下通过不同的码率，希望能以较少的冗余比特达到尽可能高的纠错性能，且译码延迟尽可能低。
  根据Claude E. Shannon信息论的定义和证明就是遍历完所有可能的信道输入分布后，找到那个可以使互信息最大的那个分布，互信息最大时就是信道容量。第二定理 (有噪信道编码定理) 给出，当信道的信息传输率R不超过信道容量C时(R < C)，采用合适的信道编码方法可以实现任意高的传输可靠性，但若信息传输率超过了信道容量(R > C)，就不可能实现可靠的传输。

It is possible to transmit information over this channel reliably (with probability of error → 0) if and only if:
$$R\stackrel{\text{ def }}{=}\frac{\#\text{ information bits }}{\text{ channel use }}<C$$信息传输率R即单位时间内由数据通路传输数据的数量，可以是每秒钟传输二进制码元的个数（比特率bit/s），也可以是每秒钟传输信号码元的个数（波特率Bd），但本文都指比特率。

  通俗点想，如果一个信道信息传输率R<C，编码符号数量足够多时，总能找到M个编码符号，组成一个码用对应的译码规则正确译码（译码错误概率任意小）。比如带限AWGN信道容量的无约束容量：$$C_{\text{awgn }}(P,B)=B{\log }_2\left(1+\frac{P}{N_{0}B}\right)\text{ bits }/s$$带宽限制为B Hz（正负两边的），均接收功率限制为P瓦特，$P=E_{b}R$是每个信息位的平均信号能量和数据传输率R的乘积。加性白高斯噪声的功率谱密度为$\frac{N_0}{2}$。

1. 想要知道频谱利用率和功率效率的香农极限，可以参考信道容量和香农极限
2. 想要知道AWGN信道之外的信道容量计算，可以参考信道容量

  虽然有噪信道编码定理理论上可以实现信道容量，但实际上需要高效的编码和译码算法才能用于通信。因为$R=\frac{k}{n}$将k个bit编码为n个bit时意味着，码表大小为$2^k$，虽然采用最大似然译码器时错误概率会在块长度趋向无穷大时迅速趋近于零，但这更说明需要高效的编码和译码算法。
  可以在多项式时间内进行编码的线性码——LDPC出现了。但最大似然译码问题的复杂性对于线性码来说问题很大，所以也引出了多项式时间的次优译码算法SPA。

1.2 LDPC码的编码表示

  LDPC编码可以用矩阵表示也可以用图表示。矩阵和图结构的稀疏性是LDPC码的算法效率的关键特性。

1.2.1 LDPC码的矩阵表示和相关术语

  LDPC码可以概括到非二进制字母表，这里暂时只考虑二进制LDPC码。LDPC码属于一类线性分组码，同样可以描述为二进制域${\mathbb{F}}_2$上的二进制n元组向量空间${\mathbb{F}}_2^n$中的一个k维子空间$\mathcal{C}$：一组基 B = { g 0 , g 1 , … , g k − 1 } B=\{\mathbf{g}_0, \mathbf{g}_1, \ldots, \mathbf{g}_{k-1}\} B={g0​,g1​,…,gk−1​}，它们生成了$\mathcal{C}$，这样每个$\mathbf{c}\in \mathcal{C}$都可以写成$\mathbf{c}=u_0{\mathbf{g}}_0+u_1{\mathbf{g}}_1+\dots +u_{k-1}{\mathbf{g}}_{k-1}$的形式。
  $\mathbf{c}=\mathbf{u}G$，其中$\mathbf{u}$是$u_i$的向量，$G$是称为生成矩阵的$k\times n$矩阵，其行是${\mathbf{g}}_i$向量。生成矩阵$G$的零空间是${\mathcal{C}}^{\perp }$，它是一个$(n-k)$维的向量空间，其中的向量$\mathbf{x}$满足$\mathbf{x}{G}^T=0$。${\mathcal{C}}^{\perp }$的基是 B ⊥ = { h 0 , h 1 , … , h n − k − 1 } B^{\perp}=\{\mathbf{h}_0, \mathbf{h}_1, \ldots, \mathbf{h}_{n-k-1}\} B⊥={h0​,h1​,…,hn−k−1​}，因此对于每个$\mathbf{c}\in \mathcal{C}$，${\mathbf{ch}}_i^T=0$对于所有$i$成立，或者更紧凑地，$\mathbf{c}{H}^T=0$。$H$被称为奇偶校验矩阵，它执行接收字的$m=n-k$个单独的奇偶校验检查。
  LDPC的名字来源就是因为其奇偶校验矩阵$H$中1的密度很低。规则LDPC码是一种线性分组码，其奇偶校验矩阵$H$的每列都包含精确的$w_c$个1 ，每行都包含精确的$w_r=w_c(n/m)$个1，其中$w_c<<m$、$w_c<<n$。码率$R=k/n$与这些参数之间的关系是$R=1-w_c/w_r$（前提是$H$满秩的）。如果$H$的1的密度较低，但每列或每行中的1的数量不是固定的，那么该码是不规则LDPC码。

1.2.2 LDPC码的图表示和相关术语

  Tanner图是LDPC码的图形表示，其中包括变量节点（v-nodes）和校验节点（c-nodes）。Tanner图的构建规则是：如果奇偶校验矩阵$H$中的元素$h_{ji}$为1，那么检查节点$j$与变量节点$i$相连。

PS：变量节点可以理解为bit或者说源符号节点，校验节点可以称为函数节点。

  举例一个$(10,5)$线性分组码，$w_r=w_c(=2)\frac{n}{m}=4$ ，其奇偶校验矩阵$H$包含了特定的1的分布： H = [ 1 1 1 1 0 0 0 0 0 0 1 0 0 0 1 1 1 0 0 0 0 1 0 0 1 0 0 1 1 0 0 0 1 0 0 1 0 1 0 1 0 0 0 1 0 0 1 0 1 1 ] H=\left[

$$\begin{array}{llllllllll} 1 & 1 & 1 & 1 & 0 & 0 & 0 & 0 & 0 & 0 \\ 1 & 0 & 0 & 0 & 1 & 1 & 1 & 0 & 0 & 0 \\ 0 & 1 & 0 & 0 & 1 & 0 & 0 & 1 & 1 & 0 \\ 0 & 0 & 1 & 0 & 0 & 1 & 0 & 1 & 0 & 1 \\ 0 & 0 & 0 & 1 & 0 & 0 & 1 & 0 & 1 & 1 \end{array}$$ \right] H= ​11000​10100​10010​10001​01100​01010​01001​00110​00101​00011​ ​
相应的Tanner图如图中所示，其中$v$节点$c_0,c_1,c_2$和$c_3$与$c$节点$f_0$相连，对应于$H$的第0行（所以H的第$i$行，表示了第$i$个校验节点$f_i$和所有变量节点的链接关系，同理，H的第$i$列，也表明了第$i$个变量节点和所有校验节点的链接关系）。Tanner图的规则性表示是指每个$v$节点具有两个边连接，每个$c$节点具有四个边连接。这与变量节点的度(或者说H矩阵的列权重)$w_c=d_v=2$和检查节点的度（或者说H矩阵的行权重）$w_r=d_c=4$相吻合，同时满足了$m{w}_r=n{w}_c$。

重点：Tanner图中的长度为$\nu$的环（或循环）是由$\nu$条边构成的路径，最终回到自身。示例中图中的六条粗边表示了一个长度为6的环。Tanner图的环的长度$\gamma$是图中的最小循环长度。 在二分图中 最短可能的循环明显是长度为4的循环（无论如何构造），这些循环在$H$矩阵中体现为位于子矩阵的角上的四个1。 短循环会降低LDPC码所使用的迭代译码算法的性能 。 这一事实将在迭代译码算法中讨论。

补充说明两个点：

不规则LDPC怎么去表示节点的分布：
  对于不规则LDPC码，参数$w_c$和$w_r$是列和行数的函数，通常不采用此表示法。相反，通常采用$\lambda (x)$和$\rho (x)$表示$v$节点和$c$节点的度分布多项式。 在多项式$\lambda (x)$中，${\lambda }_d$表示与度为$d$的$v$节点相连的边的比例，$d_v$表示最大的$v$节点度:$$\lambda (x)=\sum _{d=1}^{d_v}{\lambda }_d{x}^{d-1},$$
类似地，在多项式$\rho (x)$中，${\rho }_d$表示与度为$d$的$c$节点相连的边的比例，$d_c$表示最大的$c$节点度:$$\rho (x)=\sum _{d=1}^{d_c}{\rho }_d{x}^{d-1},$$

对举例中的规则码，$w_c=d_v=2$和$w_r=d_c=4$，我们有$\lambda (x)=x$和$\rho (x)=x^3$。

  在说完LDPC编码的表示后应该明晰一点，怎么设计LDPC编码，一个评判标准应该是怎么降低短循环；怎么设计LDPC译码，一个评判标志是应该怎么降低短循环干扰。
  本文后续统一记check节点为c节点，节点集合为$\mathbf{f}$，variable节点为v节点，节点集合为$\mathbf{c}$，并且在tanner图中，c节点在上方，v节点在下方。

1.3 怎么设计LDPC码

  梳理一下参考文献中的，最新进展请自行调研。列出来码是极其有意义的，属于是指出了LDPC码的可能研究方向，很多人沿着他们的方向做。这里没有标注每一个码的设计年代是因为没有意义，很多人已经在他们的方向上进行了改进。

最早的LDPC码：Gallager码，这里说其在1981年由Tanner推广后的形式[3]。

正则LDPC码，其奇偶校验矩阵$H$具有特定的结构。假设有两个任意整数$\mu >1$和$w_r>1$，其$H$矩阵由多个子矩阵$H_d$组成，其中每个子矩阵$H_d$都是$\mu \times \mu w_r$的矩阵。每个$H_d$的行权重为$w_r$，列权重为1。

重点特征为$H_1$具有特定的形式：

1. 对于$i=1,2,\dots ,\mu$，第$i$行包含了所有$w_r$个1且位置在列$(i-1)w_r+1$到$i{w}_r$之间。
2. 其他子矩阵$H_d$只是$H_1$的列置换。

$H$矩阵是正则的，具有维数$\mu w_c\times \mu w_r$，行权重为$w_r$，列权重为$w_c$。$H$矩阵中不一定不存在长度为4的循环，但可以通过计算机设计来避免它们。Gallager证明了这类编码的集合在$w_c\ge 3$和$w_r>w_c$的条件下具有出色的距离特性，而且这些编码具有低复杂度的编码器。

Irregular LDPC 码（luby）：假设使用消息传递解码，针对信道条件设计过的度分布多项式$\lambda (x)$和$\rho (x)$参数化的编码集合中的典型编码能够在更差的信道条件下进行可靠通信，而不在集合外部的编码则不能。

PS1：最差情况下的信道条件被称为解码阈值，而$\lambda (x)$和$\rho (x)$的优化通过密度演化算法和优化算法的组合找到的。

---

PS2：密度演化是指解码器的Tanner图中传递的各种量的概率密度函数（pdf）的演化。给定$\lambda (x)$ - $\rho (x)$对的解码阈值是通过评估编码比特的计算对数似然比的pdf来确定的。单独的优化算法是在$\lambda (x)$ - $\rho (x)$对上进行优化。

---

PS3： 一般来说，通过密度演化设计算法最适用于速率不太高的编码（$R\lesssim 3/4$）和长度不太短的编码（$n\gtrsim 5000$）。 原因在于密度演化设计算法假设$n\to \infty$（因此，$m\to \infty$），因此对于非常长的编码而言，最优的$\lambda (x)$ - $\rho (x)$对将不适用于中等长度和短编码。这些$\lambda (x)$ - $\rho (x)$对应用于中等长度和短编码会导致较高的错误率底线。针对短块对密度演进做出其他调整的除外哈

Finite Geometry Codes（包含了cyclic finite geometry codes和quasi-cyclic finite geometry codes）：基于有限几何的传统技术设计的正则LDPC码，该方法产生了循环和准循环的LDPC码。循环有限几何码往往具有相对较大的最小码字距离，准循环码的最小码字距离往往较小。使用这些技术设计的短LDPC码（$n$在200比特左右）通常优于使用伪随机$H$矩阵设计的短LDPC码。

缺点

1. 在解码中使用的奇偶校验矩阵是$n\times n$而不是$(n-k)\times n$。（可以选择$n\times n$矩阵的$(n-k)\times n$子矩阵进行解码，但性能损失通常不可忽略。）
2. $w_r$和$w_c$的值相对较大，这是不希望的，因为迭代消息传递解码器的复杂度与这些值成正比。
3. 在长度和速率的选择方面没有灵活性，只能通过编码的缩短和削减来处理这个问题。

RA, IRA, and eIRA Codes：一种兼具串行Turbo码和LDPC码特征的编码类型。

RA：用户比特被重复（通常是2次或3次），然后经过排列，最后通过一个累加器（差分编码器）进行编码。这些编码已经证明能够接近容量极限，但它们的缺点是它们通常具有较低的速率（速率$1/2$或更低）。
IRA：一些比特重复的次数多于其他比特，形成了不规则的重复累积。通常是非系统编码，尽管它们可以转换成系统编码形式，但这会大大降低速率
eIRA：系统编码，并允许低速率和高速率。

数组码（Array Codes）：采用一种独特的校验矩阵 H 结构，简化了编码。以低误差率底部而闻名，可以设计用于低码率和高码率。然而，由于 H 矩阵的特定格式，只能实现特定的码长和码率。对于应用中需要较短码长的情况，这种编码方案可以是不错的选择。

Combinatorial LDPC Codes：利用组合数学原理来设计具有特定约束条件的 LDPC 码。对于应用中需要较短码长的情况，这种编码方案可以是不错的选择。

PS：这里列的LPDC发展方向是不全面的，其他方向请自行调研。

2. SPA译码及其之后其他译码算法

  SPA（Sum-Product Algorithm）和积算法是Gallager在1960年提出的迭代解码算法，迭代译码算法也叫belief propagation algorithm (BP algorithm), and the message passing algorithm (MPA)。

PS: BP算法是SPA算法在译码中的具体应用。在LDPC码的上下文中，两者通常可以互换使用。SPA 包括“概率域”版本（计算 APP）及其对数域版本log-SPA（计算 LLR）以及某些近似值，BP一般指其对数域版本。

2.0 前置知识

后验概率：a posteriori probability (APP)

给定

1. 传输码字$\mathbf{c}=\left[\begin{array}{llll}{c}_0& c_1& \dots & c_{n-1}\end{array}\right]$
2. 接受码字$\mathbf{y}=\left[\begin{array}{llll}{y}_0& y_1& \dots & y_{n-1}\end{array}\right]$

定义

1. 第i个bit的后验概率/第i位的后验概率 : $\mathrm{Pr}\left(c_i=1\mid \mathbf{y}\right)$
2. the APP ratio/the likelihood ratio, LR : $l\left(c_i\right)\triangleq \frac{\mathrm{Pr}\left(c_i=0\mid \mathbf{y}\right)}{\mathrm{Pr}\left(c_i=1\mid \mathbf{y}\right)}.$
3. the log-APP ratio/the log-likelihood ratio (LLR) : $L\left(c_i\right)\triangleq \log \left(\frac{\mathrm{Pr}\left(c_i=0\mid \mathbf{y}\right)}{\mathrm{Pr}\left(c_i=1\mid \mathbf{y}\right)}\right)$

其中此处和后续内容均假定自然对数。

迭代译码算法目的：在已知接受码字的情况下，计算第i位的概率分布情况
  计算$\mathrm{Pr}\left(c_i=1\mid \mathbf{y}\right)、l\left(c_i\right)$ 或 $L\left(c_i\right)$
迭代译码算法思想：基于图的消息传递。

2.1 原始消息迭代算法/BP算法

  2.1.1先说迭代过程，2.1.2说迭代过程中的消息传递信息的是什么

  下文中，输出信息被分解为内信息 (intrinsic information) 和内信息 (extrinsic information)， 来自外部的信息叫外部信息，来自节点内部的信息叫做内信息，实际上外信息和外部信息由于没用特殊映射，在LR和LLR中就是线性加减关系，可以认为就是等同。

2.1.1 迭代过程

  在1.2.2中，变量节点v的集合 { c 0 , c 1 , . . . , c n − 1 } \{c_0,c_1,...,c_{n-1}\} {c0​,c1​,...,cn−1​}可以理解为编码后的bit或者说源符号节点，校验节点c的集合 { f 0 , f 1 , . . . , f n − k } \{f_0,f_1,...,f_{n-k}\} {f0​,f1​,...,fn−k​}可以称为函数节点。但是在消息迭代中，将v节点代表一种类型的处理器，c 节点代表另一种类型的处理器，边代表消息路径，会更好理解：

1. 前半次迭代（H第i列，是变量节点$c_0$和所有check节点的关系）：
     单个v节点处理其输入消息，输入消息就是在该bit位，接受到的，并将其生成的输出消息传递到相邻的 c 节点，比如：H第0列前三位为1时，v节点$c_0$相邻的c节点就是$f_0,f_1,f_2$。记从$c_0$传递给$f_2$的信息为$m_{\uparrow 02}$。
   
     所以，前半次迭代过程，就是v节点$c_0$，处理其能获得的（除了要传向的c节点$f_2$）的所有消息后得到的外部信息$m_{\uparrow 0j}$，然后将其传递给c节点$f_j$的过程。
     当然，既然变量节点有外部信息，那么肯定有内部信息和完整信息的概率。完整信息，肯定是$\mathrm{Pr}(c_0=b\mid$ input messages ) , b ∈ { 0 , 1 } ), b \in\{0,1\} ),b∈{0,1}的LR或者LLR形式，表征了$c_0$的可能状态，但玩意儿计算起来很麻烦，所以要用信息迭代的方式去近似计算这个值；外部信息是从排除了要传递信息的校验节点$f_2$外后的相邻校验节点c的集合 { f 0 , f 1 } \{f_0,f_1\} {f0​,f1​}中接收到的信息，由这些校验节点c各自根据其自身内部信息计算一个值并传递给$c_0$，而这个经过c计算其内部信息的过程，本质上是向变量节点$c_0$提供了关于其他变量节点的状态信息; 变量节点的内部信息，就是由于上一轮外部信息传递后，所更新的信息，表征了当前情况下变量节点关于其对应比特的概率。
     变量节点怎么处理得到信息这个问题马上会在2.1.2里讲。先继续把迭代过程整明白。
2. 后半次迭代
   
     现在看H矩阵的第i行，表示了c节点$f_i$和所有v节点的连接关系。在这半次迭代，就是每一个c节点，都会去处理输出信息，然后把消息传递给相邻的v节点。图4给出了$m_{\downarrow 0j}$，表示由c节点$f_0$ 传给$c_4$的外部信息。类比v节点的完整信息：$\mathrm{Pr}(c_0=b\mid$ input messages ) , b ∈ { 0 , 1 } ), b \in\{0,1\} ),b∈{0,1}。这里给出c节点的完整信息：$\mathrm{Pr}(\text{check equation }{f}_0\text{ is satisfied}\mid$ input messages ) , b ∈ { 0 , 1 } ), b \in\{0,1\} ),b∈{0,1}， 表征了$f_0$的可能状态（是否满足校验方程），c节点$f_0$要传递给第j个v节点的外部信息：$m_{\downarrow 0j}$，表征了$f_0$向其他相邻节点提供了关于其他变量节点的状态信息的过程。
     所以后半次迭代就是c节点$f_0$处理外部信息后传递给第j个v节点的外信息：$m_{\downarrow 0j}$的过程。

  经过一定数量的迭代或者达到迭代停止条件后，译码器就会根据LLR计算每一个变量节点$c_i$的bit判决情况，最普适最常见的迭代停止准则就是$\hat{\mathbf{c}}{\mathbf{H}}^T=0$，$\hat{\mathbf{c}}$上的尖尖表示tentatively解码码字，也就是暂定估计值。

2.2.1 补：消息的迭代过程的先决假设

1. 假设传递的消息在整个解码过程中在统计上是独立的。→降低解码复杂性，消息传递时不用考虑多个消息之间的影响了

   • 假设成立条件，当接受bit位 $y_i$ 独立时，如果 Tanner 图不存在环，则该独立性假设成立；如果存在周长为 $\gamma$ 的环，独立性假设仅在第 $\gamma /2$ 次迭代之前成立，之后消息开始在图的各个周期中自行循环。
   • 模拟表明，避免长度为四的周期，消息传递算法通常非常有效，且部分长度为四周期的配置是无害的。H. Tang, J. Xu, S. Lin, and K. Abdel-Ghaffar, “Codes on finite geometries,” submitted to IEEE Trans. Inf. Theory. 2002.
   • 有害性可以通过the message-passing schedule/flooding schedule大量缓解。H. Xiao and A. Banihashemi, “Message-passing schedules for decoding LDPC codes,” submitted to IEEE Trans. Comm., 2003.

2.2 传递什么样的消息——论Log-Domain SPA Decoder

引入文献中的符号表示：

• V j = { v V_j=\left\{\mathrm{v}\right. Vj​={v-nodes connected to c-node $f_j\}$，c节点$f_j$的邻接节点集合
• V j \ i = { v V_j \backslash i=\left\{\mathrm{v}\right. Vj​\i={v-nodes connected to c-node ， f j } \ { \left.f_j\right\} \backslash\left\{\right. fj​}\{ v-node $c_i\}$
• C i = { C_i=\left\{\right. Ci​={ c-nodes connected to v-node $c_i\}$
• C i \ j = { C_i \backslash j=\left\{\right. Ci​\j={ c-nodes connected to v-node c i } \ { \left.c_i\right\} \backslash\left\{\right. ci​}\{ c-node $f_j\}$
• M v ( ∼ i ) = { M_v(\sim i)=\left\{\right. Mv​(∼i)={ messages from all v-nodes except node $c_i\}$，表示v节点中，除了$c_i$以外，其他所有v节点能提供的外部信息
• M c ( ∼ j ) = { M_c(\sim j)=\left\{\right. Mc​(∼j)={ messages from all c-nodes except node $f_j\}$，下标表示信息由谁提供，~表示排除了哪个check节点
• $P_i=\mathrm{Pr}\left(c_i=1\mid y_i\right)$，已知接受符号$y_i$的值（注意经过调制解调后在波形上采样后得到的值是一个数），对应bit位值为1的概率
• $S_i=$ event that the check equations involving $c_i$ are satisfied，$S_i$是一个事件，表征在变量节点$c_i$纳入时校验方程的满足情况
• $q_{ij}(b)=\mathrm{Pr}\left(c_i=b\mid S_i,y_i,M_c(\sim j)\right)$, where b ∈ { 0 , 1 } b \in\{0,1\} b∈{0,1}. For the APP algorithm presently under consideration, $m_{\uparrow ij}=q_{ij}(b)$; for the LR algorithm, $m_{\uparrow ij}=q_{ij}(0)/q_{ij}(1)$; and for the LLR algorithm, $m_{\uparrow ij}=\log \left[q_{ij}(0)/q_{ij}(1)\right]$.变量节点的后验概率，根据消息的形式不同需要对$q_{ij}(b)$做适应操作。上检查节点下变量节点，所以这里的下表ij，表示从变量节点$c_i$到检查节点$f_j$
• $r_{ji}(b)=\mathrm{Pr}($ check equation $f_j$ is satisfied $\mid c_i=b,M_v(\sim i))$, where b ∈ { 0 , 1 } b \in\{0,1\} b∈{0,1}. For the APP algorithm presently under consideration, $m_{\downarrow ji}=r_{ji}(b)$; for the LR algorithm, $m_{\downarrow ji}=$ $r_{ji}(0)/r_{ji}(1)$; and for the LLR algorithm, $m_{\downarrow ji}=\log \left[r_{ji}(0)/r_{ji}(1)\right]$.检查节点的后验概率，根据消息的形式不同需要对$r_{ji}(b)$做适应操作

注意，消息 $q_{ij}(b)$ 和$r_{ji}(b)$ 虽然在这里被解释为概率，但本身还是随机变量 (random variables, rv’s)，因为它们是关于随机变量 $y_i$ 和其他消息的函数。

2.2.1 如何把变量节点的后验概率$q_{ij}(b)$用检查节点的后验概率$r_{ji}(b)$表示？

对于$q_{ij}(b)$的生成，用贝叶斯公式两次，将一大坨常数概率值用常数$K_{ij}$表示 q i j ( 0 ) = Pr ⁡ ( c i = 0 ∣ y i , S i , M c ( ∼ j ) ) = ( 1 − P i ) Pr ⁡ ( S i ∣ c i = 0 , y i , M c ( ∼ j ) ) / Pr ⁡ ( S i ) = K i j ( 1 − P i ) ∏ j ′ ∈ C i \ j r j ′ i ( 0 )

$$\begin{aligned} q_{i j}(0) & =\operatorname{Pr}\left(c_i=0 \mid y_i, S_i, M_c(\sim j)\right) \\ & =\left(1-P_i\right) \operatorname{Pr}\left(S_i \mid c_i=0, y_i, M_c(\sim j)\right) / \operatorname{Pr}\left(S_i\right) \\ & =K_{i j}\left(1-P_i\right) \prod_{j^{\prime} \in C_i \backslash j} r_{j^{\prime} i}(0) \end{aligned}$$ qij​(0)​=Pr(ci​=0∣yi​,Si​,Mc​(∼j))=(1−Pi​)Pr(Si​∣ci​=0,yi​,Mc​(∼j))/Pr(Si​)=Kij​(1−Pi​)j′∈Ci​\j∏​rj′i​(0)​
同理：
$$\begin{array}{r}{q}_{ij}(1)=K_{ij}{P}_i\prod _{j^{\prime }\in C_i\backslash j}{r}_{j^{\prime }i}(1).\end{array}$$常数$K_{ij}$满足约束方程$q_{ij}(0)+q_{ij}(1)=1.$
公式（1）合并为$q_{ij}(b)=K_{ij}{P}_i{\prod }_{j^{\prime }\in C_i\backslash j}{r}_{j^{\prime }i}(b).$示意图如下：

2.2.2 如何把检查节点的后验概率$r_{ji}(b)$用变量节点的后验概率$q_{ij}(b)$表示？

定理1：假设有一个长度为$M$的独立二进制数列{$a_i$}，且$Pr(a_i-1)=p_i$，则 { a i } i = 1 M \left\{a_i\right\}_{i=1}^M {ai​}i=1M​包含偶数个1的概率为$\frac{1}{2}+\frac{1}{2}{\prod }_{l=i}^M\left(1-2p_i\right)$。

根据定理1，考虑令$p_i\leftrightarrow q_{ij}(1)$，则有公式（2）：$$r_{ji}(0)=\frac{1}{2}+\frac{1}{2}\prod _{i^{\prime }\in V_j\backslash i}\left(1-2q_{i^{\prime }j}(1)\right)$$该公式的含义为：当$c_i=0$时， { c i ′ : i ′ ∈ V j \ i } \left\{c_{i^{\prime}}: i^{\prime} \in V_j \backslash i\right\} {ci′​:i′∈Vj​\i}bit位上的bit，必须包含偶数个1，才能满足校验方程$f_i$。此外，$r_{ji}(1)=1-r_{ji}(0)$。

公式（2）可归纳为：$$1-2r_{ji}(b)=\prod _{i^{\prime }\in V_j\backslash i}\left(1-2q_{i^{\prime }j}(b)\right).$$

2.2.3 MPA算法的APP’s初始化及Probability-Domain SPA Decoder的算法流程总结

  MPA算法在计算bit位的APP序列时，初始化根据H矩阵完成：
  对于$h_{ij}=1$（对应Tanner中存在连接关系的一对检查节点&变量节点），由于$q_{ij}(b)=\mathrm{Pr}\left(c_i=b\mid y_i\right)$，根据先前定义的符号，即$q_{ij}(1)=P_i$ and KaTeX parse error: Expected '\right', got 'EOF' at end of input: …_{i j}(0)=1-P_i，显然，$y_i$为实际接收到的信道符号，这个概率是和信道类型有关的，所以初始化也需要根据信道情况进行初始化：

• BEC信道， y i ∈ { 0 , 1 , E } y_i \in\{0,1, E\} yi​∈{0,1,E} where $E$ is the erasure symbol，定义擦除概率为
  $\delta =\mathrm{Pr}(y_i=$ $E\mid c_i=b)$，则初始化为： Pr ⁡ ( c i = b ∣ y i ) = { 1  when  y i = b 0  when  y i = b c 1 / 2  when  y i = E \operatorname{Pr}\left(c_i=b \mid y_i\right)=\left\{

  $$\begin{array}{cc}1 & \text { when } y_i=b \\0 & \text { when } y_i=b^c \\1 / 2 & \text { when } y_i=E\end{array}$$ \right. Pr(ci​=b∣yi​)=⎩ ⎨ ⎧​101/2​ when yi​=b when yi​=bc when yi​=E​其中，$b^c$为b的补集complement。

• BSC信道， y i ∈ { 0 , 1 } y_i \in\{0,1\} yi​∈{0,1}，定义错误概率为$\epsilon =\mathrm{Pr}\left(y_i=b^c\mid c_i=b\right)$，则初始化为： Pr ⁡ ( c i = b ∣ y i ) = { 1 − ε  when  y i = b ε  when  y i = b c . \operatorname{Pr}\left(c_i=b \mid y_i\right)=\left\{

  $$\begin{array}{c}1-\varepsilon \text { when } y_i=b \\\varepsilon \quad \text { when } y_i=b^c\end{array}$$ .\right. Pr(ci​=b∣yi​)={1−ε when yi​=bε when yi​=bc​.

• Binary Input Additive White Gaussian Noise Channel (BI-AWGNC)，假定传输时，第$i$个bit值$c_i$被BPSK映射为$x_i=1-2c_i$，即$x_i=+1(-1)$ 对应$c_i=0(1)$。那么第 $i$ 个接收到的样本是 $y_i=x_i+n_i$，其中 $n_i$ 是独立的并且正态分布为 $\eta \left(0,{\sigma }^2\right)$。此时，初始化为：$$\mathrm{Pr}\left(x_i=x\mid y_i\right)={\left[1+\exp \left(-2yx/{\sigma }^2\right)\right]}^{-1}$$其中 x ∈ { ± 1 } x \in\{ \pm 1\} x∈{±1}。

初步总结：

1. 对于 $i=0,1,\dots ,n-1$，设置 $P_i=\mathrm{Pr}\left(c_i=1\mid y_i\right)$， 其中 $y_i$ 是 $i$- th 接收到的信道符号。 然后为所有 $i,j$ 设置 $q_{ij}(0)=1-P_i$ 和 $q_{ij}(1)=P_i$，其中 $h_{ij}=1$。
2. 使用（2.2.2）中公式（2）更新 { r j i ( b ) } \left\{r_{j i}(b)\right\} {rji​(b)}。
3. 使用（2.2.1）中公式（1）更新 { q j i ( b ) } \left\{q_{j i}(b)\right\} {qji​(b)}。 求解常数 $K_{ij}$。
4. 对于 $i=0,1,\dots ,n-1$，计算$$Q_i(0)=K_i\left(1-P_i\right)\prod _{j\in C_i}{r}_{ji}(0)$$和$$Q_i(1)=K_i{P}_i\prod _{j\in C_i}{r}_{ji}(1)。$$
   其中常数 $K_i$ 的选择是为了确保 $Q_i(0)+Q_i(1)=1$。
5. 对于 $i=0,1,\dots ,n-1$，设置$${\hat{c}}_i=\{\begin{array}{l}1\text{ if }{Q}_i(1)>Q_i(0)\\ 0\text{ else }\end{array}.$$如果${\mathbf{c}}^T=0$或者迭代次数等于最大限制，则停止； 否则，转至步骤 2。

为简化上述算法中的计算过程，步骤4应在步骤3之前计算，然后步骤3可用简单除法$q_{ij}(b)=K_{ij}{Q}_i(b)/r_{ji}(b)$代替。

2.2.4 怎么由Probability-Domain SPA Decoder过渡到Log-Domain SPA Decoder？

  为什么要过渡：与概率域 Viterbi 算法一样，概率域 SPA 因为涉及许多概率乘法，会在数值上变得不稳定。因此，和log - Viterbi一样 ，SPA 的对数域版本是首选，将开销大且不稳定的概率乘法转化为开销小的加法。
  怎么实现过渡：

首先定义以下LLRs：
L ( c i ) = log ⁡ ( Pr ⁡ ( c i = 0 ∣ y i ) Pr ⁡ ( c i = 1 ∣ y i ) ) L ( r j i ) = log ⁡ ( r j i ( 0 ) r j i ( 1 ) ) L ( q i j ) = log ⁡ ( q i j ( 0 ) q i j ( 1 ) ) L ( Q i ) = log ⁡ ( Q i ( 0 ) Q i ( 1 ) )

$$\begin{aligned} L\left(c_i\right) & =\log \left(\frac{\operatorname{Pr}\left(c_i=0 \mid y_i\right)}{\operatorname{Pr}\left(c_i=1 \mid y_i\right)}\right) \\ L\left(r_{j i}\right) & =\log \left(\frac{r_{j i}(0)}{r_{j i}(1)}\right) \\ L\left(q_{i j}\right) & =\log \left(\frac{q_{i j}(0)}{q_{i j}(1)}\right) \\ L\left(Q_i\right) & =\log \left(\frac{Q_i(0)}{Q_i(1)}\right) \end{aligned}$$ L(ci​)L(rji​)L(qij​)L(Qi​)​=log(Pr(ci​=1∣yi​)Pr(ci​=0∣yi​)​)=log(rji​(1)rji​(0)​)=log(qij​(1)qij​(0)​)=log(Qi​(1)Qi​(0)​)​

初始化过程中的转化为：
（ B E C ） L ( q i j ) = L ( c i ) = { + ∞ , y i = 0 − ∞ , y i = 1 0 , y i = E （ B S C ） L ( q i j ) = L ( c i ) = ( − 1 ) y i log ⁡ ( 1 − ε ε ) （ B I − A W G N C ） L ( q i j ) = L ( c i ) = 2 y i / σ 2

$$\begin{aligned} &（BEC） L\left(q_{i j}\right)=L\left(c_i\right)=\left\{\begin{array}{c} +\infty, y_i=0 \\ -\infty, y_i=1\\ 0, y_i=E \end{array}\right. \\ & （BSC）L\left(q_{i j}\right)=L\left(c_i\right)=(-1)^{y_i} \log \left(\frac{1-\varepsilon}{\varepsilon}\right) \\ & （BI-AWGNC）L\left(q_{i j}\right)=L\left(c_i\right)=2 y_i / \sigma^2 \end{aligned}$$ ​（BEC）L(qij​)=L(ci​)=⎩ ⎨ ⎧​+∞,yi​=0−∞,yi​=10,yi​=E​（BSC）L(qij​)=L(ci​)=(−1)yi​log(ε1−ε​)（BI−AWGNC）L(qij​)=L(ci​)=2yi​/σ2​

Step1转化过程：
  想要简化对数计算和除法运算，根据$\tanh \left[\frac{1}{2}\log \left(p_0/p_1\right)\right]=p_0-p_1=1-2p_1$，考虑公式（2）可转化为$$\tanh \left(\frac{1}{2}L\left(r_{ji}\right)\right)=\prod _{i^{\prime }\in V_j\backslash i}\tanh \left(\frac{1}{2}L\left(q_{i^{\prime }j}\right)\right)$$
  想要简化tanh函数，首先将$L\left(q_{ij}\right)$分解成符号函数${\alpha }_{ij}=\mathrm{sign}\left[L\left(q_{ij}\right)\right]$*绝对值${\beta }_{ij}=\left|L\left(q_{ij}\right)\right|$的形式：$$L\left(q_{ij}\right)={\alpha }_{ij}{\beta }_{ij}$$进而将公式（2）转化为：$$\tanh \left(\frac{1}{2}L\left(r_{ji}\right)\right)=\prod _{i^{\prime }\in V_j\backslash i}{\alpha }_{i^{\prime }j}\cdot \prod _{i^{\prime }\in V_j\backslash i}\tanh \left(\frac{1}{2}{\beta }_{i^{\prime }j}\right)$$LLR代入后化简得到公式（4）： L ( r j i ) = ∏ i ′ α i ′ j ⋅ 2 tanh ⁡ − 1 ( ∏ i ′ tanh ⁡ ( 1 2 β i ′ j ) ) = ∏ i ′ α i ′ j ⋅ 2 tanh ⁡ − 1 log ⁡ − 1 log ⁡ ( ∏ i ′ tanh ⁡ ( 1 2 β i ′ j ) ) = ∏ i ′ α i ′ j ⋅ 2 tanh ⁡ − 1 log ⁡ − 1 ∑ i ′ log ⁡ ( tanh ⁡ ( 1 2 β i ′ j ) ) = ∏ i ′ ∈ V j \ i α i ′ j ⋅ ϕ ( ∑ i ′ ∈ V j \ i ϕ ( β i ′ j ) )

$$\begin{aligned} L\left(r_{j i}\right) & =\prod_{i^{\prime}} \alpha_{i^{\prime} j} \cdot 2 \tanh ^{-1}\left(\prod_{i^{\prime}} \tanh \left(\frac{1}{2} \beta_{i^{\prime} j}\right)\right) \\ & =\prod_{i^{\prime}} \alpha_{i^{\prime} j} \cdot 2 \tanh ^{-1} \log ^{-1} \log \left(\prod_{i^{\prime}} \tanh \left(\frac{1}{2} \beta_{i^{\prime} j}\right)\right) \\ & =\prod_{i^{\prime}} \alpha_{i^{\prime} j} \cdot 2 \tanh ^{-1} \log ^{-1} \sum_{i^{\prime}} \log \left(\tanh \left(\frac{1}{2} \beta_{i^{\prime} j}\right)\right) \\ & =\prod_{i^{\prime} \in V_j \backslash i} \alpha_{i^{\prime} j} \cdot \phi\left(\sum_{i^{\prime} \in V_j \backslash i} \phi\left(\beta_{i^{\prime} j}\right)\right) \end{aligned}$$ L(rji​)​=i′∏​αi′j​⋅2tanh−1(i′∏​tanh(21​βi′j​))=i′∏​αi′j​⋅2tanh−1log−1log(i′∏​tanh(21​βi′j​))=i′∏​αi′j​⋅2tanh−1log−1i′∑​log(tanh(21​βi′j​))=i′∈Vj​\i∏​αi′j​⋅ϕ ​i′∈Vj​\i∑​ϕ(βi′j​) ​​其中，$\varphi (x)=-\log [\tanh (x/2)]=\log \left(\frac{e^x+1}{e^x-1}\right)$

$\varphi (x)$函数性质：${\varphi }^{-1}(x)=\varphi (x)$ when $x>0$
该特性能够简化计算，或者更直接的，The function may be implemented by a look-up table.

Step2将公式（1）转化为公式（3）：$$L\left(q_{ij}\right)=L\left(c_i\right)+\sum _{j^{\prime }\in C_i\backslash j}L\left(r_{j^{\prime }i}\right).$$

step3转化为公式（5）：$$L\left(Q_i\right)=L\left(c_i\right)+\sum _{j\in C_i}L\left(r_{ji}\right).$$

总结：

1. 对于 $i=0,1,\dots ,n-1$， 所有 $i,j$ 满足 $h_{ij}=1$的，设置设置 $L\left(q_{ij}\right)$。
2. 更新$\left\{L\left(r_{j^i}\right)\right\}$ 用公式（4）.
3. 更新 { L ( q i j ) } \left\{L\left(q_{i j}\right)\right\} {L(qij​)}用公式（3）.
4. 更新 { L ( Q i ) } \left\{L\left(Q_i\right)\right\} {L(Qi​)}用公式（5）.
5. 获得估计值${\hat{c}}_i=\{\begin{array}{l}1\text{ if }L\left(Q_i\right)<0\\ 0\text{ else }\end{array}$进行$\hat{\mathbf{c}}{\mathbf{H}}^T=0$校验。如果校验通过，则输出$\hat{\mathbf{c}} $；否则转2。

2.2.5 最小和算法Min-Sum Decoder和其他变体

  与概率域 SPA 算法相比，对数域 SPA 算法的复杂度更低，数值更稳定。Min-Sum Decoder是复杂度更低的解码器，在性能方面通常只受到一点影响，性能上通常约弱 0.5dB。
  Min-Sum Decoder的思想是信息在节点之间传递时，只需进行最小值的计算：
ϕ ( ∑ i ′ ϕ ( β i ′ j ) ) ≃ ϕ ( ϕ ( min ⁡ i ′ β i ′ j ) ) = min ⁡ i ′ ∈ V j \ i β i ′ j .

$$\begin{aligned} \phi\left(\sum_{i^{\prime}} \phi\left(\beta_{i^{\prime} j}\right)\right) & \simeq \phi\left(\phi\left(\min _{i^{\prime}} \beta_{i^{\prime} j}\right)\right) \\ & =\min _{i^{\prime} \in V_j \backslash i} \beta_{i^{\prime} j} . \end{aligned}$$ ϕ(i′∑​ϕ(βi′j​))​≃ϕ(ϕ(i′min​βi′j​))=i′∈Vj​\imin​βi′j​.​
  此时，step2的更新公式转化为：$$L\left(r_{ji}\right)=\prod _{i^{\prime }\in V_j\backslash i}{\alpha }_{i^{\prime }j}\cdot \underset{i^{\prime }\in V_j\backslash i}{\min }{\beta }_{i^{\prime }j}.$$
  注意：在 BI-AWGNC 情况下，当采用最小和算法时，初始化 $L\left(q_{ij}\right)=2y_i/{\sigma }^2$ 可以替换为 $L\left(q_{ij}\right)=y_i$。 优点是在这种情况下不需要知道噪声功率 ${\sigma }^2$。

min-sum变体1：Normalized Min-Sum Decoding
$$L\left(r_{ji}\right)=（系数\alpha ）\prod _{i^{\prime }\in V_j\backslash i}{\alpha }_{i^{\prime }j}\cdot \underset{i^{\prime }\in V_j\backslash i}{\min }{\beta }_{i^{\prime }j}.$$，$\alpha \in (0,1]$

min-sum变体2：Offset Min-Sum Decoding
$$L\left(r_{ji}\right)=\prod _{i^{\prime }\in V_j\backslash i}{\alpha }_{i^{\prime }j}\cdot (\underset{i^{\prime }\in V_j\backslash i}{\min }{\beta }_{i^{\prime }j}-系数\beta ).$$，其中$\beta \ge 0$

3.关于消息传递算法及其变体算法的性能分析

  在长码LDPC情况下，2304bit+200次迭代次数

短块长度LDPC 码的性能存在译码的有限长度效应，导致 error floors 的出现。