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

Posted on by

描述目标函数

​ 之前的文章问题引入：自相关引入了自相关问题以及其经典生成算法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