流形优化:设计低相关序列生成算法

描述目标函数

​ 之前的文章问题引入:自相关引入了自相关问题以及其经典生成算法CAN。这次我们将该问题一般化,考虑MIMO系统,即设计一个包含M个低自相关序列的序列。
{\bold x_m}=[x_m(1)…x_m(N)]^T,m=1,…,M
​ 其中|x_m|=1,给定两个序列{\bold x_i},{\bold x_j},这两个序列的周期相关性为:
r_{i,j}(k)=\sum_{n=k+1}^Nx_i(n)x_j^*(n-k),k=0,…,N-1\space \space \space i,j=1,…,M\tag{1}
​ 上式中,当i\neq j时,(1)式表示{\bold x_i}{\bold x_j}之间的互相关;当i= j时,(1)式表示{\bold x_i}{\bold x_j}之间的自相关。我们的目标是设计一个低自相关和低互相关的序列,也就是说我们需要优化以下问题:
\min_{{\bold x_m}_{m=1}^M}\sum_{m=1}^{M}\sum_{n=1-N}^{N-1}|r_{m,m}(n)|^2+\sum_{i=1}^{M}\sum_{{j=1}\\j\neq i}^{M}\sum_{n=1-N}^{N-1}|r_{i,j}(n)|^2\\ s.t.\space|x_m(n)|=1,n=1,…N\space m=1,…,M\tag{2}
​ 设\bold X=[\bold x_1\space \bold x_2 \space …\space \bold x_M]_{N\times M},可将上式写为更紧凑的形式:
\bold X^H \bold X=\begin{bmatrix}\bold x_1^H\\\bold x_2^H \\ \vdots\\ \bold x_M^H \end{bmatrix}[\bold x_1\space \bold x_2 \space …\space \bold x_M]=\begin{bmatrix} \bold x_1^H\bold x_1&\bold x_1^H\bold x_2 &\cdots&\bold x_1^H\bold x_M\\ \bold x_2^H\bold x_1&\bold x_2^H\bold x_2 &\cdots&\bold x_2^H\bold x_M\\ \vdots & \vdots & \ddots & \vdots \\ \bold x_M^H\bold x_1&\bold x_M^H\bold x_2 &\cdots&\bold x_M^H\bold x_M\\ \end{bmatrix}\tag{3}
​ 由于{\bold x_i}为单模序列,即\bold x^H\bold x=N,则(3)式可写成:
\bold X^H \bold X=\begin{bmatrix} N&\bold x_1^H\bold x_2 &\cdots&\bold x_1^H\bold x_M\\ \bold x_2^H\bold x_1&N &\cdots&\bold x_2^H\bold x_M\\ \vdots & \vdots & \ddots & \vdots \\ \bold x_M^H\bold x_1&\bold x_M^H\bold x_2 &\cdots&N\\ \end{bmatrix}\tag{4}
​ 注意到,上式除了对角线外的项均为互相关,r_{i,j}(0)=\bold x_i^H \bold x_j,于是可以将(4)式与互相关联系起来:
\sum_{i=1}^{M}\sum_{{j=1}\\j\neq i}^{M}\sum_{n=1-N}^{N-1}|r_{i,j}(0)|^2=||\bold X^H\bold X-N\bold I_M||_F^{2}
​ 为了用\bold X表示(2)式中的其他相关项,引入大小为N\times N的矩阵\bold S_n,它的定义如下:
\begin{aligned} \bold S_n(i,j)&=\begin{cases} 1,&\text{if } i-j=n\\ 0,&\text{Otherwise} \end{cases}\space i,j=1,\cdots,N\\ &=\begin{bmatrix} 0&\cdots &0&\cdots & 0\\ \vdots&\ddots&\vdots&\ddots&\vdots\\ 1&0&0&\cdots&0\\ \vdots&\ddots&\vdots&\ddots&\vdots\\ 0&0&1&\cdots&0\\ \end{bmatrix}\\ &=\begin{bmatrix} 0&0\\ \bold I_n&0 \end{bmatrix} \end{aligned}
​ 注意到,r_{i,j}(n)=\bold x_i^H \bold S_n\bold x_j,则(2)式中的其他相关项可以用如下式子表示:
2\sum_{n=1}^{N-1}||\bold X^H\bold S_n\bold X||^2_F
​ 令\bold x=vec(\bold X)\in \mathbb C^L,L=M\times N,则(2)式可以表示成如下形式:
\begin{aligned} &\min_{\bold x}||\bold X^H\bold X-N\bold I_M||_F^{2}+2\sum_{n=1}^{N-1}||\bold X^H\bold S_n\bold X||^2_F\\ &s.t.\space|\bold x_l|=1,l=1,…L\tag{5} \end{aligned}
​ 为了让目标函数更灵活,给(4)式各项增加权重系数,如下:
\begin{aligned} &\min_{\bold x}f(\bold x)=\gamma_0||\bold X^H\bold X-N\bold I_M||_F^{2}+2\sum_{n=1}^{N-1}\gamma_n||\bold X^H\bold S_n\bold X||^2_F\\ &s.t.\space\bold x\in S^L\tag{6} \end{aligned}
​ 其中\gamma_n为非负实数,通过更改\gamma_n的值,可以按比例更改时间滞后为n的相关项的值。若对所有的n均取\gamma_n=1,则(6)式与(5)式与(2)式等价。上式中S^L为复数圆流形(complex circle manifold),定义为:
S^L\triangleq{\bold x\in\mathbb C^L:|\bold x_l|=1,l=1,…L}
​ 在优化上述问题之前,我们再在(6)式后加一项,用来控制目标函数的收敛性(后文有证明)。
\begin{aligned} &\min_{\bold x}\bar f(\bold x)=\gamma_0||\bold X^H\bold X-N\bold I_M||_F^{2}+2\sum_{n=1}^{N-1}\gamma_n||\bold X^H\bold S_n\bold X||^2_F+\varsigma \bold x^H\bold x\bold x^H\bold x\\ &s.t.\space\bold x\in S^L\tag{7} \end{aligned}
​ 上式中,\varsigma\geq0。由于x^H\bold x=N\times M=L,\bold x^H\bold x\bold x^H\bold x=L^2是常数,所以加上该项后不会影响原式的解。

流形优化

复数圆流形


​ 在流形优化:从一个简单的例子开始流形优化:初步的理论分析中,介绍了流形优化的步骤和基础理论。接下来运用上述知识解决复数圆流形(complex circle manifold)上的优化。

​ 之前介绍过球流形上的流形优化,球流形与复数原流形十分相似,若把复数原流形中每个元素的实部和虚部拆开将其变为S^L\in\mathbb R^{L \times 2},那么该流形的每一行都是一个2阶的球流形。那么,我们先从复数原流形中的每个元素开始分析,对于向量\bold z\in S^L,\bold w\in \mathbb C^L,令\bar {\bold z}_l,\bar {\bold w}_l为二维向量,该二维向量的第一项分别为\bold z,\bold w的第l项的实部,第二项分别为\bold z,\bold w的第l项的虚部,即
\bar {\bold z}_l=[\text{Re}{\bold z_l} \space \text{Im}{\bold z_l}]^T=[\bold z_{lr}\space\bold z_{li}]^T\in \mathbb R^2\\ \bar {\bold w}_l=[\text{Re}{\bold w_l} \space \text{Im}{\bold w_l}]^T=[\bold w_{lr}\space\bold w_{li}]^T\in \mathbb R^2
​ 对于单位圆流形,线性空间中的向量向单位元流形上的投影为\text{Proj}_{\bar {\bold z}}(\bar {\bold w})
\begin{aligned} \text{Proj}_{\bar {\bold z}_l}(\bar {\bold w_l})&=\bar {\bold w_l}-({\bar {\bold z_l}}^T\bar {\bold w_l}){\bar {\bold z_l}}\\ &=\bar {\bold w_l}-(\bold z_{lr}\bold w_{lr}+\bold z_{li}\bold w_{li}){\bar {\bold z_l}}\\ &=\bar {\bold w_l}-\text{Re}(\bold w_{l}^{\star}\bold z_l){\bar {\bold z_l}}\\ \end{aligned}
​ 对于从单位圆流形上一点出发的向量\bold v,拉回函数可以写成
R(\bold v)=\frac{\bold v}{||\bold v||_2}
​ 由于复数圆流形中每个元素的实部虚部可以构成单位圆流形,因此对于复数圆流形,线性空间中的向量向单位元流形上的投影为\text{Proj}_{ {\bold z}}( {\bold w})
\text{Proj}_{\bar {\bold z}}(\bar {\bold w})= {\bold w}-\text{Re}(\bold w^*\bold z){ {\bold z}}
​ 对于从复数圆流形上一点出发的向量\bold v,拉回函数可以写成
R(\bold v)=\bold w\odot\frac{1}{|\bold w|}

目标函数梯度


​ 有了投影函数和拉回函数,我们便可以在流形上用梯度下降法优化非凸问题。还有一个重要的一步是对目标函数求梯度。对于原问题(7)式,为了方便求导,对其做一些变换。由||\bold A||^2_F=tr(\bold A \bold A^H),可得
\begin{aligned} f(\bold x)&=\gamma_0||\bold X^H\bold X-N\bold I_M||_F^{2}+2\sum_{n=1}^{N-1}\gamma_n||\bold X^H\bold S_n\bold X||^2_F+\varsigma \bold x^H\bold x\bold x^H\bold x\\ &=\gamma _0 \text{tr}(\bold X^H\bold X\bold X^H\bold X+N^2\bold I_M-2N\bold X^H\bold X)\\&\space\space\space+2\sum_{n=1}^{N-1}\gamma_n\text {tr}(\bold X^H\bold S_n\bold X\bold X^H\bold S_n^H\bold X)+\varsigma \bold x^H\bold x\bold x^H\bold x\\ &=\gamma _0 \text{tr}(\bold X^H\bold X\bold X^H\bold X)-\gamma_0N^2M\\&\space\space\space+2\sum_{n=1}^{N-1}\gamma_n\text {tr}(\bold X^H\bold S_n\bold X\bold X^H\bold S_n^H\bold X)+\varsigma \bold x^H\bold x\bold x^H\bold x\tag{8} \end{aligned}
​ 接下来求(8)式的梯度,将(8)式拆开,第一部分为
\text{tr}(\bold X^H\bold X\bold X^H\bold X)=\text{vec}(\bold X)^H\text{vec}(\bold X)\text{vec}(\bold X)^H\text{vec}(\bold X)=\bold x^H\bold x\bold x^H\bold x
​ 则该部分的梯度为:
\begin{aligned} \frac{\partial\text{tr}(\bold X^H\bold X\bold X^H\bold X)}{\partial \bold x}&=\frac{\partial \bold x^H\bold x\bold x^H\bold x}{\partial \bold x}\\ &=\frac{\partial (\bold x^H\bold x)^2}{\partial \bold x^H\bold x}\frac{\partial \bold x^H\bold x}{\partial \bold x}\\ &=4\bold x\bold x^H\bold x \end{aligned}
​ 对于第三部分,梯度与第一部分相同。

​ 令\text {vec}(\bold S_n)=\bold s_n,对于(8)式,第二部分为:
{tr}(\bold X^H\bold S_n\bold X\bold X^H\bold S_n^H\bold X)=\bold x^H\bold s_n\bold x^H\bold x\bold s_n^H\bold x
​ 值得注意的是,上式中变量的顺序有所交换,这是因为\text{tr}(\bold A\bold A^H)=\text{tr}(\bold A^H\bold A)=\text{vec}(\bold A)^H\text{vec}(\bold A)

​ 该部分的梯度为:
\begin{aligned} \frac{\partial{tr}(\bold X^H\bold S_n\bold X\bold X^H\bold S_n^H\bold X)}{\partial \bold x}&=\frac{\partial [(\bold x^H\bold s_n)(\bold x^H\bold x)(\bold s_n^H\bold x)]}{\partial \bold x}\\ \tag{9} \end{aligned}
​ 由于
\partial (\bold{xyz})=(\partial\bold{x})\bold{yz}+\bold x(\partial\bold{y})\bold z+\bold{xy}(\partial\bold{z})\tag{10}
​ 且
\frac{\partial (\bold x^H\bold s_n)}{\partial \bold x}=\frac{\partial (\bold s_n^H\bold x)}{\partial \bold x}=\bold s_n\\ \frac{\partial (\bold x^H\bold x)}{\partial \bold x}=2\bold x\\ \tag{11}
​ 结合(9)(10)(11)式,可以得出
\frac{\partial [(\bold x^H\bold s_n)(\bold x^H\bold x)(\bold s_n^H\bold x)]}{\partial \bold x}=2(\bold s_n\bold x^H\bold x\bold s_n^H\bold x+\bold x\bold s_n\bold x^H\bold s_n^H\bold x)
​ 根据上述分析,可以得出(8)式的梯度:
\nabla_{\bold x}f(\bold x)=4(\gamma_0+\varsigma)\bold x\bold x^H\bold x+4\sum_{n=1}^{N-1}\gamma_n(\bold s_n\bold x^H\bold x\bold s_n^H\bold x+\bold x\bold s_n\bold x^H\bold s_n^H\bold x)

Armijo-Goldstein线搜索


​ 在线性搜索中,精确的寻找搜索补偿通常不是十分有效,Armijo和Goldstein分别提出了不精确一维搜索过程。在梯度下降法中,设第k步的步长为d_k,需要确保
image-20200710111115216
​ 设满足这样条件的\alpha的集合为J,即:
image-20200710111401981
​ 以下图为例,集合J代表的范围为MH。

image-20200710100518275

​ 现为了让搜索范围不靠近J的端点,一个合理的要求是:

image-20200710111530875

​ 为了避免搜索范围太小,再加上另一个要求:

f(x_k+s_k)>f(x_k)+(1+\rho) g_k^Ts_k

在复数圆流形优化过程中,为了优化迭代次数,我们使用Armijo线搜索。

复数圆流形优化的收敛性


实验结果

​ 对于(7)式中的目标函数f(\bold x),取\gamma_n=1,设计参数\varepsilon用来表示拟合程度:
\varepsilon=\frac{f(\bold x)}{MN^2}
​ 将复数圆流形优化算法(CGD)与经典的WeCAN算法作比较:

image-20200710103107304

​ 对于(7)式中的目标函数f(\bold x),取\gamma_n为:
\gamma_n=\begin{cases} 1, &n\in[1,19]\cup[236,255]\\ 0, &n\in[20,235] \end{cases}
​ 与WeCAN的对比结果如下:

image-20200710103416046