QR分解及其应用

《矩阵分析与应用》专题报告——QR分解及应用
学生姓名：卢楠、胡河群、朱浩
2015年11月25日
目录
1 引言 (3)
2 QR分解 (4)
2.1QR分解的性质 (4)
2.2 QR分解算法 (5)
2.2.1 采用修正Gram-Schmidt法的QR分解 (5)
2.2.2 Householder QR分解 (6)
2.2.3 采用Givens旋转的QR分解 (8)
3 QR分解在参数估计中的应用 (9)
3．1 基于QR分解的参数估计问题 (9)
3. 2基于Householder变换的快速时变参数估计 (12)
3. 3基于Givens旋转的时变参数估计 (14)
4 QR分解在通信系统中的应用 (16)
4.1 基于QR分解的稳健干扰对齐算法 (16)
4.2基于QR分解的MIMO置信传播检测器 (19)
总结 (21)
参考文献 (22)
1 引言
矩阵分解是指将一个矩阵表示为结构简单或具有特殊性质的若干矩阵之积或之和，大体上可以分为满秩分解、QR分解和奇异值分解。

矩阵分解在矩阵分析中占有很重要的地位，常用来解决各种复杂的问题。

而QR分解是工程中应用最为广泛的一类矩阵分解。

QR分解是目前求一般矩阵全部特征值的最有效并广泛应用的方法，一般矩阵先经过正交相似变换成为Hessenberg矩阵，然后再应用QR分解求特征值和特征向量。

它是将矩阵分解成一个正交矩阵Q与上三角矩阵R，所以称为QR分解。

参数估计是在已知系统模型结构时，用系统的输入与输出数据计算系统模型参数的过程。

它在系统辨识和无线通信领域有着广泛的应用。

18世纪末德国数学家C.F.高斯首先提出参数估计的方法，他用最小二乘法计算天体运行的轨道。

20世纪60年代，随着电子计算机的普及，参数估计有了迅猛的发展。

参数估计有很多方法，如矩估计、极大似然法、一致最小方差无偏估计、最小风险估计、同变估计、最小二乘法、贝叶斯估计、极小极大熵法等。

其中最基本的是最小二乘法和极大似然法。

本文将重点介绍QR分解及其在参数估计和通信系统中的应用。

2 QR 分解
2.1QR 分解的性质
定理2.1.1（QR 分解）若m n R ⨯∈A ，且m n ≥，则存在列正交矩阵m n R ⨯∈Q 和上三角矩阵m n R ⨯∈R 使得A =QR 。

当m n =时，Q 是正交矩阵。

如果A 是非奇异的n n ⨯ 矩阵，则R 的所有对角线元素均为正，并且在这种情况下Q 和R 二者是唯一的。

若A 是复矩阵，则Q 和R 取复值。

注意到T T T A A =(QR)(QR)=R R ,因此可以得出结论：T G =R 是T A A 的下三角Cholelskey 因子。

由于这个原因，在关于估计的文献中，矩阵R 常称为平方根滤波器（算子）。

下面的引理称为矩阵分解引理，它在矩阵的QR 分解的应用中是一个很有结果。

引理2.2.1 若A 和B 是任意两个m n ⨯ 矩阵，则 H H A A =B B (2.1.1) 当且仅当存在一个m m ⨯ 酉矩阵Q ，使得
QA =B
(2.1.2)
证明 充分性证明：若QA =B ，并且Q 是酉矩阵，则
H H H H B B =A Q QA =A A 。

必要性证明：令A 和B 的奇异值分解分别为 ∑H
A A A A =U V ∑H
B B B B =U V
式中，A U 和B U 均为m m ⨯ 酉矩阵；A V 和B V 都是n n ⨯ 酉矩阵；而m n ⨯ 矩阵∑A 和∑B 分别包含了矩阵A 和B 的非负奇异值。

由于
∑∑H H H
A A A A A A =V V ∑∑H H H
B B B B B B =V V
若H H A A =B B ,则有A B V =V 和∑=∑A B 。

定义矩阵
H
B A Q =U U
易知 ∑∑H H H H
B A B A A A A B B B QA =U U A =U U U V =U V =B
这就证明了引理的必要条件[10]。

2.2 QR 分解算法
2.2.1 采用修正Gram-Schmidt 法的QR 分解
矩阵A 的QR 分解可以利用Gram-Schmidt 正交化方法实现。

Gram-Schmidt 正交化方法原本是一种由n 个向量12n a ,a ,...,a 构造互相正交且范数为1的向量
12n q ,q ,...,q 的方法。

将向量1a 标准正交化的结果取作1q ，即
1111R R ⎫=⎪
⎬=⎪⎭1
11a q q
(2.2.1)
然后，从2a 中除去与1a 平行的向量，再进行标准正交化，并将结果取作2q ，则有
1222121222()R R R R R ⎫
=⎪
=-⎬⎪=-⎭
H 12
21221q a a q q a q (2.2.2)
进而，又从3a 除去与1a 和2a 平行的两个分量，再进行标准正交化，并使用该结果作3q ，即有
1323331323132333()R R R R R R R R ⎫
=⎪=⎪
⎬=--⎪⎪=--⎭
H 13
H
233123312q a q a a q q q a q q
(2.2.3)
如此继续，则对于(2)k k n ≤≤q 有
111
331,1j k 1R ()H jk j k k kk k j jk j k k k j jk j R R R R -=-=⎫
⎪
=≤≤-⎪
⎪⎪
=-⎬⎪
⎪⎪
=-⎪⎭
∑∑q a a q q a q
(2.2.4)
容易验证，i q 是标准正交基，即满足
i j ij δ=H q q
(2.2.5)
其中，ij δ 为Kronecker δ 函数。

如果令m n ⨯ 矩阵A 的列向量,,,12n a a ...a ，则以1,2,,n q q ...q 为列向量的矩阵Q 与A 之间有下列关系：
A =QR (2.2.6)
又由于i q 组成标准正交基，所以
H
n Q Q =I
将A 与Q 重写在同一矩阵，应用以上Gram-Schmidt 正交化的方法叫做经典Gram-Schmidt 正交化法[6]。

2.2.2 Householder QR 分解
Householder 变换可以实现任意m n ⨯ 矩阵A 的QR 分解，其原理是使用变维向量的Householder 变换，使得该向量除第一个元素外，其他元素皆为0。

根据Householder 变换的相关知识，欲使一个p 维向量12,,...,x T
p x x ⎡⎤=⎣⎦
x 的
第1个元素后面的所有元素变为0，则p 维的Householder 向量应取
=
w (2.2.7)
式中
1
11
, x x x ββ=-=-
x x (2.2.8)
假定m n ⨯ 矩阵A 的列分块形式为 1,2,,...m n n ⨯⎡⎤=⎣⎦A a a a
首先令111,21,,1...m a a a ⎡⎤==⎣⎦T
x a ，并取p m = ，则按照式(2.2.7)和式(2.2.8)，可以计算得到1m =u w 。

此时，
(1)(1)(1)
1111112 ,,,n ⎡⎤=-→==⎣⎦T H I u u A H A a a ...a
(2.2.9)
变换后，矩阵1A 的第1列(1)
1a 的第一个元素等于()
12
22
211211...a m a a +++ ，而该列的其他元素全部为0。

第二步针对矩阵1A 的第2列(1)
2a ，令1p m =- 和
(1)(1)(1)
22322,,...,m a a a ⎡⎤=⎣
⎦T
x 又可按照式(2.2.7)和式(2.2.8)求出（m-1）维向量1m -w 。

此时，取210m -⎡⎤
=⎢⎥⎣⎦u w ，又可得到
(1)(2)(2)
2222212112,,,n ⎡⎤=-→===⎣⎦T H I u u A H A H H A a a ...a (2.2.10)
变换后，矩阵2A 的第1列与1A 的第1列相同，而第2列(1)
1a 的第一个元素等于(1)
12
a ，第二个元素等于(
)
12
222(1)(1)(1)22
32
2
...a
m a a
+++ ，而该列的其他元素全部为
0。

类似地，又可针对矩阵2A 的第3列设计Householder 变换矩阵3H ，使得2A
的第一、二个元素保持不变，其他元素组成的m-2维向量(2)(2)(2)
33433,,...,m a a a ⎡⎤=⎣⎦T
x 变
换为除第一个元素外的全部元素变为0。

假定矩阵A 经过k-1次Householder 变换后，已变成(1)k -A ，即
(1)(2)111(1)(1)(1)
12... ,,...,, k=2,3,...
k k k k k k k n -------==⎡⎤=⎣⎦A H A H H A
a a a
并且其前k-1列具有以下变换结果: (1)
(1)(1)
1,...,,0,...,0, j=1,2,...k-1k k k j
j jj ---⎡⎤=⎣⎦T
a
a a
因此，第k 次Householder 变换的目的就是保持前k-1列不变，实现(1)k -A 列第
k 列的下述变换：
(1)(),,(1)1,(1),00k k k k k k k k k k k m k a a a a --+-⎡⎤⎡⎤⎢⎥⎢⎥⎢⎥⎢⎥=⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎣
⎦⎣⎦H
这相当于对矩阵(1)k -A 进行Householder 变换(1)
k k -H A 时取
1k k k -⎡⎤=⎢⎥⎣
⎦I
0H 0H n 次Householder 变换后，即可实现QR 分解。

2.2.3 采用Givens 旋转的QR 分解
Givens 旋转也可以用来计算QR 分解。

这里以43⨯ 矩阵为例，说明Givens QR
分解的思想：
(3,4)(2,3)(1,2)(3,4)(2,3)(3,4)000000000000
00
0⨯⨯⨯⨯
⨯⨯⊗
⨯⨯⨯⨯⨯⎡⎤⎡⎤⎡⎤⎡⎤
⎢⎥⎢⎥⎢⎥⎢⎥⨯⨯⨯⊗⨯⨯⊗⨯⨯⨯⨯⎢⎥⎢⎥⎢⎥⎢⎥−−−→−−−→−−−→−−−→⎢⎥
⎢⎥
⎢⎥
⎢⎥⊗⨯⨯⊗⨯⨯⨯⨯⊗⨯⎢⎥⎢⎥⎢⎥⎢⎥⊗⨯⨯⨯⨯⨯⨯⊗⨯⎣⎦⎣⎦⎣⎦⎣⎦
⨯⨯⨯⨯
⨯⨯⨯
⨯⎡⎤⎡⎤⎢⎥⎢⎥⊗⨯⨯⨯⎢⎥⎢⎥−−−→−−−→⎢⎥
⎢⎥
⊗⨯⊗⎢⎥⎢⎥⨯⊗⎣⎦⎣⎦0000
00⨯⎡⎤⎢⎥⨯⨯⎢⎥⎢⎥
⨯⎢⎥⎣⎦
其中，⊗ 表示用Givens 旋转进行变化你的元素。

从上述说明中易得出结论：如果令j G 代表约化过程中的第j 次Givens 旋转，则T Q A =R 是上三角矩阵，其中t t-11Q =G G LG ，而t 是总的旋转次数。

3 QR 分解在参数估计中的应用
3．1 基于QR 分解的参数估计问题
现在以系统辨识为例，说明如何利用矩阵的QR 分解进行系统参数的递推估计。

令系统在k 时刻的输入为()x k ，系统输出的观测值由卷积方程
()()()()()()0p
T k i i y k x k e k x k i e k e k θθ=*+=-+=+∑x θ (3.1.1)
给出，其中，* 表示离散卷积，()e k 代表k 时刻的观测误差，且
()()()01,1,,,,
,T
k T
p x k x k x k p θθθ=--⎡⎤⎣⎦
⎡⎤=⎣⎦
x θ (3.1.2)
若将1,2,,k n = 的所有观测数据组成一向量，则
n n n =+y A θe
(3.1.3)
式中，()()()1,2,
,T
n y y y n =⎡⎤⎣⎦
y ，()()()1,2,
,T
n e e e n =⎡⎤⎣⎦
e ，
()()()1,2,
,T
n x x x n =⎡⎤⎣⎦A 。

系统辨识问题的提法是：已知系统输入()x k 和输出观测值()y k ，其中，
1,2,
,k n = ，估计系统参数向量θ [7]。

在时变系统的辨识中，则要求在已估计
n 时刻的系统参数向量n θ 的情况下，使用增加的()()1,1x n y n ++ 值，通过简
单的运算，递推出1n + 时刻的系统参数向量1n θ+ 。

n 时刻的系统辨识问题可以简化为最小二乘问题
22
min n
n n n
-θA θy (3.1.4)
求解，并且其解由“法方程”
T T n n n n n xx n n ==A A θA y R θr 或
(3.1.5)
确定。

式中，T xx n n =R A A 代表系统输入()x k 的协方差矩阵，T
n n n =r A y 。

直接求解式(3.1.5)的方法叫做协方差方法。

例如，先计算协方差矩阵xx R
的Cholesky 分解T xx GG =R ,然后利用回带法解三角矩阵1
T n n -=G θG R 直接得到n θ 。

然而，由于T
xx n n =R A A 的条件数是n A 的条件数的平方，因此，直
接计算式(3.1.5)的得到的解有可能是严重病态的（即条件数很大），即使n A 本身的条件数并不大，不是严重病态的。

在系统参数向量θ 的自适应递推辨识中，标准的递推最小二乘()RLS 法和T U DU 分解法都是针对协方差矩阵xx R 进行更新的。

虽然T U DU 分解（其中，U 为上三角矩阵，D 为对角矩阵）在数值上比较稳定，但是这些递推辨识方
法也同样存在条件数变大的毛病。

相比之下，n A 的QR 分解可以保持原问题的条件数不变。

不妨令
n T n n ⎡⎤
=⎢⎥⎣⎦
R Q A O
(3.1.6)
式中，n Q 是n n ⨯ 正交矩阵，n R 是()()11p p +⨯+ 上三角矩阵，而O 为
()()11n p p --⨯+ 维零矩阵。

由于正交变换可以保持被变换向量的Euclidean 长度或范数不变，所以式(3.1.4)的最小二乘问题可等价写作 22
min n
T T
n n n n n
-θQ A θQ y (3.1.7)
或 222
2
min n
n n n
n
-+θR θy y (3.1.8)
式中
n T
n n n ⎡⎤
=⎢⎥⎢⎥⎣⎦
y Q y y
(3.1.9)
且n y 为()11p +⨯ 向量，n y 为()11n p -⨯ 向量，它们可以从T
n n Q y 直接分块
得到。

一旦获得了n y ，即可由n n n =R θy 得到θ 。

解此方程需要()1/2n n + 次计算，并且最小残差值等于22
n
y 。

假定增加了两个已知数值()1x n + 和()1y n + ，我们来讨论如何更新系统参数的估计，即使用已估计的参数向量n θ 和简单的运算，得到1n + 时刻的新估计1n +θ 。

为了减少过去数据数据对参数估计的影响，对数据()x k 和()y k 采取指数加权，即1n + 时刻的数据矩阵和观测数据向量分别取作
()111,1n n n n T n y n λλ+++⎡⎤⎡⎤
==⎢⎥⎢⎥+⎣⎦⎣⎦
y A A y x (3.1.10)
式中，01λ≤< 称为遗忘因子，且()()()11,,,T
n x n x n x n p +=+-⎡⎤⎣⎦x 。

于是，
可以写出式(3.1.4)在()1n + 时刻的形式为
()2111
2
2
12
arg min arg min 1n n n n n T n y n λλ++++=-⎡⎤
⎡⎤=-⎢⎥⎢⎥+⎣⎦⎣⎦
θ
θ
θA θy y A θx (3.1.11)
乍一看，上式似乎没有什么特别吸引人之处，其实不然。

这是因为，如同下面的引理所述，式(3.1.11)的极小化变量等价为下述式的极小化变量，而后者非常适合于递推更新。

引理3.1.2r 若n n n ⎡⎤=⎢⎥⎣⎦R A Q O ，n T
n n n ⎡⎤=⎢⎥⎢⎥⎣⎦
y Q y y ，其中，n Q 是正交矩阵，n R 是
上三角矩阵，则式(3.1.11)的极小化变量等同于下式的极小化变量：
()2111
2
2
12
arg min arg min 1n n n n n T n y n λλ++++=-⎡⎤
⎡⎤
=-⎢⎥⎢⎥+⎣⎦⎣⎦
θ
θ
θA θy y A θx (3.1.12)
证明见[]。

如果将式(3.1.12)的极小化变量记作1n +θ ，则以上讨论可总结为1n +θ 的自适应递推估计算法如下。

算法1（系统参数的自适应估计算法）
步骤1 对矩阵1n T n λ+⎡⎤=⎢⎥⎣⎦R R x 进行QR 分解，得
1111n n T T n n T n λ++++⎡⎤⎡⎤==⎢⎥⎢⎥⎣
⎦⎣⎦R R Q R Q x O (3.1.13)
式中，1n +Q 是()()11n n +⨯+ 正交矩阵，1n +R 为()()11p p +⨯+ 上三角矩阵，且O 是()()1n p p -⨯+ 零矩阵。

步骤二 进行分块运算
111
1n T
n n n ++++⎡⎤=⎢⎥⎢⎥⎣⎦
y Q y
y (3.1.14)
其中，1n +y 为()11p +⨯ 向量，1n +y 为()1n p -⨯ 向量。

步骤三 切结三角矩阵方程111n n n +++=R θy 得到1n +θ 。

3.2基于Householder 变换的快速时变参数估计
考察()1n p ⨯+ 矩阵
1112
1,121222,1,1,2,1p p n n n n p a a a a a a a a a +++⎡⎤⎢⎥
⎢⎥=⎢⎥⎢
⎥⎣⎦
A 的Householder QR 分解，即
***
11121,1**
22
2,1*1,100
00000
0p p n n p p a a a a a a ++++⎡⎤⎢⎥⎢⎥
⎢⎥⎢
⎥=⎢⎥⎢⎥
⎢⎥⎢⎥⎢⎥⎣
⎦
H
A (3.2.1)
显然，只需要进行p 次Householder 变换即可。

换言之，为了得到上述QR 分解，应该选择n H 为p 个Householder 变换矩阵之积，即
()()()11n n n n p p =-H H H H
(3.2.2)
式中
()/,1
,2,,T n j j j j j p σ=-=H I u u (3.2.3)
是对矩阵()
()()()121j n
n n n n j =-A H H H A 第j 列向量()()
()12,,
,T
j j j j j nj a a a ⎡⎤⎣⎦ 进行的
Householder 变换矩阵，其参数选择方法为
()()()()
()
(),1,2,
,0,sgn ,,j j j j j jj j j j jj jj j j ij a j p i j u i a a i j a i j ασααα⎫
⎪⎪⎪=⎪⎪⎪
=+=⎬⎪
⎧⎪<⎪⎪⎪⎪
=+=⎨⎪⎪⎪>⎪⎩⎭
(3.2.4)
其中 ()
()1j j
T n
n j j +=-A A u q
(3.2.5)
并且
()/j
T T j j n j σ=q u A
(3.2.6)
递推的Householder QR 分解算法如下：
基于QR 分解的自适应参数估计算法一般由两个分开的过程组成：（1）递推更新
QR 分解T =Q A R 中的上三角矩阵R ；（2）用回代法求解三角矩阵方程。

由于直接的回代需要()2
O m 次运算（m 为数据长度）。

因此，即便Householder 变
换再快速，整个自适应算法也至少需要()2
O m
次运算。

文献[]将上述快速
Householder QR 分解算法和求解三角矩阵方程的回代法综合起来考虑，提出了只具有()O m 复杂度的快速自适应算法。

3.3基于Givens 旋转的时变参数估计
现在考虑另外一种递推方法，递推求解n θ 的变化量n δ ，而不是直接递推求
1n +θ 本身。

换句话说，令
1n n n +=+θθδ
(3.3.1)
问题是如何更新n δ 。

假定正交矩阵Q 为已知，它满足
11n n T n λ++⎡⎤⎡⎤=⎢⎥⎢⎥⎣
⎦⎣⎦R R Q x O (3.3.2)
由式(3.1.11)，式(3.3.1)和式(3.3.2)易知，n δ 是下式的极小化变量：
()11argmin 1n n n n n n n T T n n y n λλλ++⎧⎫⎡⎤⎡⎤⎡⎤⎪⎪
=--⎨⎬⎢⎥⎢⎥⎢⎥+⎣⎦⎣⎦⎪⎪⎣⎦⎩
⎭δR R y δQ δQ θx x (3.3.3)
此式又可化简为
()1arg min 1n n n n u n +⎡⎤⎡⎤
=-⎢⎥⎢⎥+⎣⎦⎣⎦
δ0R δδQ 0 (3.3.4)
式中，()()111T
n n u n y n ++=+-x θ 。

因此，n δ 可以从三角矩阵方程
11n n n ++=R δy
(3.3.5)
解出，其中，1k +y 满足
(
)()111n u n r n +⎡⎤⎡⎤=⎢⎥⎢⎥++⎣⎦⎣⎦0y Q (3.3.6)
为了求出满足式(3.3.2)的Q ，可以使用Givens 平面旋转进行清零，将式(3.3.2)中的行向量1T
n +x 的全部元素变成零。

由于Q 必须左乘式(3.3.6)，所以对增广的矩阵
()1
1n
T n u n λ+⎡⎤⎢⎥+⎣⎦
R 0
x
(3.3.7)
执行所需要的清零。

综合以上分析，在每一步递推更新中需要的步骤如下[5]：
（1）计算预测误差11T
k k k y φθ++- ；
（2）形成式(3.3.7)中的()()11n n +⨯+ 矩阵；
（3）利用一系列Givens 旋转将上述矩阵最底一行的左边n 个元素扫除为零； （4）解上三角矩阵方程(3.3.5)得到k δ 。

利用Givens 旋转求解方程=A θy 的递推最小二乘算法的程序见文献[]。

该算法中，同时对矩阵A 和向量y 应用Givens 旋转，因此无需存储正交矩阵Q 。

4 QR 分解在通信系统中的应用
4.1 基于QR 分解的稳健干扰对齐算法
考虑K 用户MIMO 干扰信道，每个发送端的天线数为t N ，每个接收端的天线数为r N ，每个用户对应的自由度为[]12,,
,k d d d ，此处的自由度代表每个用
户能使用的独立数据流个数。

为了让系统自由度达到最大值，即
t min(,N )2r K N ，那么每个发送端所提供的信号空间的维数应该相等，故此处不妨设12K d d d d ==
= ，并假设在同一时刻同一频率上的各个发送接收对之
间的信道是平坦衰落的，且信道系数独立同分布。

在一个特定的时频资源上，接收端i 的接收信号可以表示为
1,K
i ii i i ji j j i j j i
=≠=+
+∑
y H Ws H W s n
(4.1.1)
其中维数为r t N N ⨯ 的ii H 和ji H 分别是发送端i 和j 到接收端i 的信道矩阵。

i W 和j W 分别是发送端i 和j 对应接收端i 和j 的预编码矩阵，且满足
i H i i d =W W I ，j H
j j d =W W I 。

维数为1i d ⨯ 的i s 是接收端i 的下行数据矢量信号，且满足功率约束()H i i E P i ⎡⎤=⎣⎦s s 。

维数为1r N ⨯ 的i n 是均值为0，方差为1的加性高斯白噪声噪声，且r H
i i N E ⎡⎤=⎣⎦n n I 。

干扰对齐往往要求完美的CSI ，但在实际通信系统中，发送端得到CSI 常常是有误差的。

为了构建稳健的干扰对齐算法，此处引入信道误差变量
ji ji ji =-E H H ，ji H 表示真实的信道矩阵，ji H 表示具有误差的信道矩阵，并且假设ji E 的元素服从均值为0，方差为2
e σ 的循环对称复高斯分布（GSCG ），即
满足2
ji ji e Nr E σ⎡⎤=⎣⎦
H E E I 。

故式(4.1.1)变为
1,()()K
i ii ii i ii ji ji j j i j j i
=≠=-+
-+∑
y H E Ws H E W s n
(4.1.2)
此时整个系统的联合接收信号可以表示为
1111121
122212
222121121
11111222222212 k k k K K k k kk k k k
k kk K K K ⎡⎤⎡⎤⎡⎤⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥==⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥
⎣
⎦⎣⎦⎣⎦⎡⎤⎡⎤⎡⎤⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥-+⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎣⎦⎣⎦⎣⎦y W s H H H y W s H H H Y y W s H H H E E E W s n E E E W s n E E E W s n
(4.1.3)
对得到的误差联合信道矩阵H 进行QR 分解有
1121
111
21122212
22212k k k k kk k
k
kk ⎡⎤⎡⎤
⎢⎥⎢⎥⎢⎥
⎢⎥===⎢⎥⎢⎥⎢⎥⎢⎥⎣⎦
⎣⎦R R R H H H R R H H H H Q QR R H H H (4.1.4)
其中，Q 是维数为r r KN KN ⨯ 的酉矩阵，R 是维数为r t KN KN ⨯ 的上三角矩阵。

因为Q 是酉矩阵，根据矩阵理论可知R 和H 有相同的统计特性，所以定义R 为系统的误差等效联合信道矩阵。

根据式(4.1.4)，式(4.1.3)可以改写为
112111122
2221121111112
22222212 k k kk K K k k k k
kk K K K ⎡⎤⎡⎤
⎢⎥⎢⎥⎢
⎥⎢⎥=⎢⎥⎢⎥⎢⎥⎢⎥⎣⎦
⎣
⎦⎡⎤⎡⎤⎡⎤
⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥-+⎢⎥⎢⎥⎢⎥⎢
⎥⎢⎥⎢⎥⎣⎦⎣⎦
⎣⎦R R R W s R R W s Y Q R W s E E E W s n E E E W s n E E E W s n
(4.1.5)
这是考虑联合接收，对式(4.1.5)的联合信号接收信号进行左乘1-Q 的预处理，得到如下的联合接收信号：
11
2111121
1111122
21222
222221121121
1122221 k k k k kk k k
kk K K K K k k kk K --⎡⎤⎡⎤⎡⎤⎡⎤
⎢⎥⎢⎥⎢⎥⎢⎥⎢
⎥⎢⎥⎢⎥⎢⎥=-⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎣⎦⎣⎦
⎣
⎦⎣⎦⎡⎤⎡⎤⎢⎥⎢⎥⎢⎥⎢⎥+=⎢⎥⎢⎥⎢⎢⎥⎣⎦⎣⎦R R R E E E W s W s R R E E E W s W s Y Q R E E E W s W s R R R n R R n Q
R n 1122111112112221222212 K K k k K K K k
k kk ⎡⎤⎢⎥⎢⎥⎢⎥⎥⎢⎥⎣⎦
⎡⎤⎡⎤⎡⎤⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥-+⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥
⎣⎦⎣⎦⎣⎦W s W s W s W s n E E E W s n E E E W s n E E E (4.1.6)
因为1-Q 是酉矩阵，于是ij E 和ij E 有相同的统计特性，同理j n 和j n 有相同的统计特性。

通过式(4.1.6)，在接收端i 经过干扰抑制矩阵i U 处理后，接收端i 的接收信号为
1122(1)1122120,0,,,, ,,H i i ii i i Ki K K H H i i i
Ki i i K K +⎡⎤⎢⎥
⎢⎥⎡⎤=⎣⎦⎢⎥
⎢⎥⎣⎦
⎡⎤⎢⎥
⎢⎥⎡⎤-+⎣⎦⎢⎥⎢⎥⎣⎦
W s W s y U R R R W s W s W s U E E E U n W s
(4.1.7)
通过稳健干扰对齐算法得到最优干扰对齐矩阵opt
i U 和opt
i W ，具体算法流程总结为[8]：
（1） 初始化i W ，这里可以随机选择均值为0，方差为1的矩阵
{},1,2,
,,i H i i i d i K ∀∈=W W W I 。

（2） 计算出{},1,2,,i i K ∀∈U ，并且单位化i U （3） 计算出{},1,2,
,i i K ∀∈W ，并且单位化i W
（4） 重复步骤（2）和（3），直到收敛。

4.2基于QR 分解的MIMO 置信传播检测器
在一个N 个发射天线和M 个接收天线的MIMO 系统中，[]1,,T
N s s s = 是
1N ⨯ 传输信号向量。

系统的输入输出关系可以写为
=+y Hs n
(4.2.1)
其中，[]1,,T
m y y =y 是1M ⨯ 接收信号向量，n 是1M ⨯ 噪声向量，n 的元素是零均值、方差2σ 的独立同分布（i.i.d ）复高斯随机变量。

H 是M N ⨯ 的
MIMO 信道矩阵，()rank N =H 。

由QR 分解，MIMO 信道矩阵H 可以写为0T
T T ⎡⎤=⎣⎦H Q R R ,其中Q 是M M ⨯
酉矩阵，0R 是()M N N -⨯ 非零矩阵，R 是N N ⨯ 上三角矩阵
1,11,2
1,12,2
2,2,000
N N N N N r r r r
r r ⎡⎤⎡⎤⎢⎥⎢⎥⎢
⎥⎢⎥==⎢⎥⎢⎥⎢⎥⎢⎥⎣⎦⎣
⎦
r r R r (4.2.2)
由于H
M M ⨯=Q Q I , 0R 是非零矩阵， 0T H N N ⨯⎡⎤=⎣⎦I R Q n R 。

式(4.2.1)可以写
为
0T H
N N ⨯⎡⎤==+⎣⎦y I R Q y Rs n
(4.2.3)
这样，R 相比于由H 得出的全连接的偶图就是一个包含更少的边数和环数的的偶图。

如下图所示。

图中，4N M == ,1Q = 。

在第一次迭代前，全部0
nm α 初始化为0，其中1n NQ ≤≤ 且()1m g n ≤≤。

在第
l 次迭代时，每一个l
mn
β 可以由最大对数近似得到
()()11
10
21
::1,21::1,max max N m n k N m n k NQ
l
l m mn
km
m a A b f f m b k n
NQ
l m km
m a A b f f m b k n
y y β
α
α
-+=-+=-∈==≠-∈==≠⎧⎫⎪⎪--+⎨⎬⎪⎪⎩
⎭
⎧⎫⎪
⎪
---+⎨⎬⎪⎪⎩
⎭
∑
∑
r a r a (4.2.4)
其中，n b 和k b 是a 中的第()n f m - 比特和第()k f m - 比特对应的s 中的第
n 比特和第k 比特，1
l km α- 表示第1l - 次迭代时第k 节点发出、第m 节点检
测的消息。

计算出l
mn β 后，有
()
1,g n l
l mn
kn k k m
α
β=≠=
∑
(4.2.5)
其中1n NQ ≤≤ 且()1m g n ≤≤ 。

第l 次迭代后n b 的软输出l
n L 是
()
1g n l l
n
kn k L β==∑ (4.2.6)
其中，1n NQ ≤≤ 。

这样，给出的QR BP - 检测器就由上三角信道矩阵R 得出的偶图来进行操作。

上述QR BP - 的计算复杂度主要由式(4.2.3)中的线性变换和式(4.2.4)中l
mn β
的计算决定。

线性变换的复杂度开销主要来自于QR 分解，，它的复杂度是
()2O MN [3]。

另外，l mn β 的计算的大部分复杂度开销来自2
m m y -r a 的计算。

文献[1]给出了ν 次迭代后该算法的计算复杂度是212M
mQ m O MN m ν=⎫
⎛+⋅⎪ ⎝⎭
∑ .所以
给出的QR BP - 检测器的计算复杂度小于文献[1]中标准BP 检测器的复杂度
()2NQ O MN ν⋅ 。

总结
通过此次文献调研，我们获益匪浅。

首先，我们了解了矩阵分解和矩阵变换，学习了QR分解、Householder变换、Givens旋转。

其次，我们通过文献调研，了解了QR分解算法在实际中的应用，如参数估计，实际的MIMO系统等。

最后，我们也从中学到了团队分工和协作的重要性。

参考文献
[1]Hu, J., and Duman, T.M.: ‘Graph-based detection algorithms for layered
space-time architectures’, IEEE J. Sel. Areas Commun., 2008, 26, (2),pp.
269–280
[2]Sangjoon P. and Sooyong C., "QR decomposition aided belief propagation
detector for MIMO systems", IET Journals & Magazines Electronics
Letters,2015, 51,(11),pp.873 - 874
[3]Kim, T.: ‘Low-complexity sorted QR decomposition for MIMO systems based on
pairwise column symmetrization’, IEEE Trans. Wirel. Commun., 2014, 13, (3), pp. 1388–1396
[4]Bobrow J.E., Murray W., An algorithm for RLS identification of parameters that
vary quickly with time,IEEE Trans Automatic Control,993,38:351-354
[5]Jianwen Gao,Quasi-Orthogonal Space-time Block Code with Givens Rotation for
OFDM System,IEEE Letters,2013,pp.64—70
[6]Jian Wu and Shu fang ,QR decomposition and gram Schmidt orthogonalization
based low-complexity multi-user MIMO precoding,IEEE conference,2014，pp.1296—1304
[7]Syed A.Raza,Direct and parallel QR based subspace decomposition methods for
system identification,IEEE conference 2014，pp.471—480
[8]谢显中，一种基于QR分解的稳健干扰对齐算法，2015（08）1957—1964
[9]陈慧，基于稀疏矩阵表示的MIMO雷达多目标定位2015（06）260—267
[10]张贤达，矩阵分析与应用。