基于张量分解的DS-CDMA盲数据检测研究 Research on Blind Data Detection for DS-CDMA Based on Tensor Decomposition

doi:10.12677/MOS.2023.123275

Modeling and Simulation
Vol. 12 No. 03 ( 2023 ), Article ID: 66450 , 18 pages
10.12677/MOS.2023.123275

基于张量分解的DS-CDMA盲数据检测研究

周婧，戴泽安，张志超

●How to Cite this Article

南京信息工程大学数学与统计学院，江苏南京

收稿日期：2023年4月3日；录用日期：2023年5月24日；发布日期：2023年5月31日

摘要

为了节约建造天线的成本，本文借助张量分解技术，在DS-CDMA系统中建立PARAFAC模型，固定用户数量(6)、扩频码数(4)和符号快照的块长度数(50)，改变信噪比和天线的数量，通过蒙特卡洛方法模拟生成多维阵列信号，执行COMFAC算法进行盲数据检测来求出误码率，进而在相同的信噪比条件下将误码率与天线数量进行数据拟合，最后结合拟合函数曲线和拟合方程得出最优天线数。研究结果表明拟合优度都达到90%以上，高信噪比的最优天线数少，稳定为5根；低信噪比的最优天线数多，并且不稳定。

关键词

张量分解，PARAFAC模型，COMFAC算法，DS-CDMA，蒙特卡洛

Research on Blind Data Detection for DS-CDMA Based on Tensor Decomposition

Jing Zhou, Ze’an Dai, Zhichao Zhang

School of Mathematics and Statistics, Nanjing University of Information Science and Technology, Nanjing Jiangsu

Received: Apr. 3^rd, 2023; accepted: May 24^th, 2023; published: May 31^st, 2023

ABSTRACT

In order to save the cost of antenna construction, the PARAFAC model is established in DS-CDMA system by means of tensor decomposition technology. Then the number of users (6), spread spectrum code (4) and block length of symbolic snapshot (50) are fixed. Also the signal-to-noise ratio and the number of antennas are changed. Multi-dimensional array signal is simulated by Monte Carlo method. COMFAC algorithm is used to detect blind data to calculate bit error rate. Then, under the same SNR condition, the bit error rate and the number of antennas were fitted. Finally, the optimal number of antennas was obtained by combining the fitting function curve and fitting equation. The results show that the goodness of fit is more than 90%, and the number of optimal antennas with high signal-to-noise ratio is less, and the stable number is 5. The optimal antennas with low signal-to-noise ratio are numerous and unstable.

Keywords:Tensor Decomposition, PARAFAC Model, COMFAC Algorithm, DS-CDMA, Monte Carlo

This work is licensed under the Creative Commons Attribution International License (CC BY 4.0).

http://creativecommons.org/licenses/by/4.0/

1. 引言

张量是描述多维向量和对偶空间的笛卡尔上多重映射关系的几何体 [1] ，是进行多维阵列信号处理的有力工具 [2] 。随着科技的飞速发展，扩展频谱通信技术是现代通信技术的热点之一。因此，国内外学者对张量分解技术在直接序列码分多址(Direct Sequence-Code Division Multiple Access, DS-CDMA)盲数据检测中的应用进行大量的研究。

Sidiropoulos N D等人建立了平行因子分析(Parallel Factor Analysis, PARAFAC)模型与DS-CDMA系统盲多用户分离均衡检测之间的联系，基于低秩三向阵列分解和三线性交替最小二乘回归的唯一性，提出了一种确定性盲PARAFAC接收机，研究结果表明：该接收机具有性能接近非盲最小均方误差的优点 [3] 。Sørensen M等人介绍几种正交约束条件下的PARAFAC的数值算法，这就为不相关信号的接收提供了基础，经过研究结果发现PARAFAC接收机的性能最好 [4] 。Fernandes C A R等人提出了用于上行多用户协同分集系统的新型盲接收机，同时还引入CANDECOMP/PARAFAC (CP)模型，最后通过计算机模拟结果发现各项性能显著提高 [5] 。上述研究表明PARAFAC模型在盲多用户分离检测中的优越性。

在PARAFAC模型的基础上，刘旭和许宗泽研究DS-CDMA系统的盲检测算法，并且提出一种基于正交约束的PARAFAC的DS-CDMA盲接收机，通过研究结果发现具有正交约束的三线性最小二乘算法在处理DS-CDMA比传统的三线性最小二乘算法有更低的克拉美罗界。此外，正交约束–平行因子分析接收机改善了误码率的性能，收敛速度也大大提高 [6] 。Bro R等人提出了复平行因子分析(Complex Parallel Factor Analysis, COMFAC)算法，该算法用于三线性并行因子分析模型拟合到DS-CDMA信号中，同时基于对若干辅助子程序的组合和改进，研究结果表明：算法达到缩短计算时间的目的 [7] 。De Lathauwer L和Castaing J提出了一种基于张量的DS-CDMA信号在天线阵列上的盲分离机制，研究结果表明：相对较高的用户数量，新技术比ALS (Alternate Least Squares)算法计算盲接收机更加准确，并且计算要求更低 [8] 。曾卉露提出基于复平行因子多故障盲分离方法，克服了传统算法在分解的时候运算时间比较长的困难，通过研究发现了方法的有效性 [9] 。于桂晨提出一种基于耦合随机投影和耦合张量分解的大规模张量分解算法，通过对比研究结果表明大规模张量分解算法比传统的算法具有更高的精度；而相对传统张量分解算法虽然有性能损失，但是获得了更少的空间占用和更快的运行效率 [10] 。上述研究表明基于张量分解的算法相对于传统的算法而言，具有提高收敛速度、缩短计算时间等显著的优点。

基于张量分解的技术，蒙特卡洛思想在模拟仿真阶段起着关键性的作用。张珽研究了蒙特卡洛方法在通信系统中采用不同的调制方式下受高斯白噪声信道影响的误码率，并且通过研究证明蒙特卡洛方法的正确性 [11] 。曾璐和谢晓尧通过MATLAB研究不同噪声在直接序列扩频通信系统中误码率与信噪比之间的关系，结果表明：数字信号在复高斯白噪声信道产生的误码率比在高斯白噪声信道中产生的误码率略大 [12] 。

综上所述，国内外相关专家、学者对张量分解技术在DS-CDMA盲数据检测中的应用进行了非常深入的研究。本文就是在DS-CDMA系统中建立PARAFAC模型进行分解，通过蒙特卡洛方法模拟生成多维阵列信号，借助COMFAC算法对模型求解得到误码率的大小，最后在相同的信噪比条件下，将误码率与天线数进行数据拟合来研究在各个信噪比条件下最优天线的个数，在仿真实验中发现规律来达到节约资源的目的。

2. 张量的理论基础

本小节主要介绍有关张量的数学理论知识，包含张量的数学表示和基本的运算。张量是一个定义在一些向量空间和一些对偶空间的笛卡尔积上的多重线性映射 [13] ，是多维数组的简称。本文研究的重点是张量在通信领域中的应用问题，因此，本小节主要以三阶张量为主展开介绍。

2.1. 张量的表示

1) 三阶张量

三阶张量通常用花体符号表示，记作 $χ \in ℂ^{L \times M \times N}$ ，它的元素表示为 $X_{l m n} ， l \in L, m \in M, n \in N$ 。 $L 、 M 、 N$ 表示三个不同的维度，结合不同的实际问题就是代表三个不同的物理意义。

2) 张量纤维

三阶张量的三路阵列称为张量纤维，即改变一个下标，固定其他下标不变得到每一路阵列。三阶张量可以用不同的模式来表达。模式-1表示三阶张量的行纤维，记作 $X_{l : n}$ ；模式-2表示三阶张量的列纤维，记作 $X_{: m n}$ ；模式-3表示三阶张量的管纤维，记作 $X_{l m :}$ 。对于高阶张量，模式-n表示N阶张量有N种不同的纤维。

3) 张量切片

张量切片是指将三阶张量按照指定单个索引将产生一个特定方向的矩阵。模式-1表示水平切片，记作 $X_{l : :}$ ；模式-2表示侧面切片，记作 $X_{: m :}$ ；模式-3表示正向切片，记作 $X_{: : n}$ 。

4) 张量的矩阵化

张量的矩阵展开是对高阶张量进行降维的过程。对张量 $χ \in ℂ^{L \times M \times N}$ 进行n模展开，记作 $m a t_{n} (χ)$ 或者 $X_{(n)}$ 。按照Kiers水平展开方式 [14] 进行矩阵化，可以得到模式-1、模式-2、模式-3的水平展开分别为如下式所示：

$X_{(1)} = X^{L \times M N} = [X_{: : 1}, \dots, X_{: : N}]$ (1)

$X_{(2)} = X^{M \times L N} = [X_{1 : :}, \dots, X_{L : :}]$ (2)

$X_{(3)} = X^{N \times L M} = [X_{: 1 :}, \dots, X_{: M :}]$ (3)

如图1所示是具体描述了三阶张量 $χ \in ℂ^{L \times M \times N}$ 的Kiers水平展开模式。

2.2. 张量的运算

1) 张量的外积

对于N阶张量 $A \in ℂ^{I_{1} \times I_{2} \times \dots \times I_{N}}$ 和M阶张量 $B \in ℂ^{J_{1} \times J_{2} \times \dots \times J_{M}}$ ，上述二者外积定义为 $C = A \circ B$ ，其中 $C \in ℂ^{I_{1} \times I_{2} \times \dots \times I_{N} \times J_{1} \times J_{2} \times \dots \times J_{M}}$ ，具体如下所示：

$c_{i_{1}, i_{2}, \dots, i_{N}, j_{1}, j_{2}, \dots, j_{M}} = a_{i_{1}, i_{2}, \dots, i_{N}} \cdot b_{j_{1}, j_{2}, \dots, j_{M}}$ (4)

Figure 1. Kiers horizontal expansion of third order tensors

图1. 三阶张量的Kiers水平展开图

2) 张量的内积

对于N阶张量 $χ_{1} \in ℂ^{I_{1} \times I_{2} \times \dots \times I_{N}}$ 和张量 $χ_{2} \in ℂ^{I_{1} \times I_{2} \times \dots \times I_{N}}$ ，它们的内积定义如下所示：

$〈 χ_{1}, χ_{2} 〉 \overset{d e f}{=} 〈 v e c (χ_{1}), v e c (χ_{2}) 〉 = (v e c {(χ_{1})}^{H} v e c (χ_{2})) = \sum_{i_{1} = 1}^{I_{1}} \sum_{i_{2} = 1}^{I_{2}} \dots \sum_{i_{n} = 1}^{I_{N}} a^{*}_{i_{1}, i_{2}, \dots, i_{n}} b_{i_{1}, i_{2}, \dots, i_{n}}$ (5)

3) 张量的Frobenius范数

定义张量 $χ$ 的Frobenius范数为：

${‖ χ ‖}_{F} = \sqrt{〈 χ, χ 〉}$ (6)

4) Kronecker积

对于任意二维矩阵 $A \in ℂ^{I \times J}$ 和矩阵 $B \in ℂ^{M \times N}$ 的Kronecker积是一个 $I M \times J N$ 的矩阵，记作 $A \otimes B$ ，具体运算过程如下所示：

$A \otimes B = [\begin{matrix} a_{11} B & a_{12} B & \dots & a_{1 n} B \\ a_{21} B & a_{12} B & \dots & a_{2 n} B \\ ⋮ & ⋮ & ⋱ & ⋮ \\ a_{m 1} B & a_{m 2} B & \dots & a_{m n} B \end{matrix}] \in ℂ^{I M \times J N}$ (7)

5) Khatri-Rao积

任意两相同列数的矩阵 $A \in ℂ^{K \times J}$ 和 $B \in ℂ^{I \times J}$ 的Khatri-Rao积记作 $A ⊙ B$ ，具体如下所示：

$A ⊙ B = [a_{1} \otimes b_{1}, a_{2} \otimes b_{2}, \dots, a_{j} \otimes b_{j}] \in ℂ^{K I \times J}$ (8)

其中Khatri-Rao积还有如下的运算性质 [15] ：

$A ⊙ B = [\begin{array}{l} B D_{1} (A) \\ B D_{2} (A) \\ ⋮ \\ B D_{j} (A) \end{array}]$ (9)

$(A \otimes B) (C ⊙ D) = A C ⊙ B D$ (10)

其中 $D_{j} (A)$ 表示 $J \times J$ 的对角矩阵，其中的元素表示矩阵 $A$ 中的第j行中的元素。

6) Kruskal秩

任意一个矩阵 $A \in ℂ^{I \times J}$ 的Kruskal秩，记为k-秩或者 $k_{A}$ ，它表示矩阵的任意r个列向量都是线性无关的最大整数 [16] 。矩阵的k-秩小于等于矩阵的秩。如果由列向量构成的矩阵 $A \in ℂ^{K \times J}$ ，并且这些列向量是从绝对连续分布中独立提取出的元素，矩阵不仅满秩而且满足k-秩，即 $r a n k (A) = k_{A} = \min (I, J)$ [16] 。

7) 秩一张量

如果一个三阶张量 $χ \in ℂ^{L \times M \times N}$ 当且仅当它能表示成三个向量的外积的形式，即

$χ = a_{1} \circ a_{2} \circ a_{3}$ (11)

称此张量为秩-1张量。对于N阶张量，以此类推。

3. 模型的建立

3.1. PARAFAC模型(Parallel Factor Analysis)

PARAFAC (Parallel Factor Analysis)模型中文解释为平行因子分解法，又被称作CANDECOMP、PARAFAC或者CANDECOMP/PARAFAC (CP)。它实际上是对高维张量进行拆分的方法，其核心思想是用有限个秩1张量的和来近似某一个张量，使其误差达到最小的目的。PARAFAC模型在统计学、信号处理中得到广泛应用。

任意一个 $I \times J \times K$ 三阶张量 $χ$ ，其中的每一个元素记为 $x_{i j k}$ ， $i = 1, 2, \dots, I$ ， $j = 1, 2, \dots, J$ ， $k = 1, 2, \dots, K$ 。 $χ_{i} (J \times K)$ 、 $χ_{j} (I \times K)$ 和 $χ_{k} (I \times J)$ 表示三阶张量在三个不同的方向上的切片，其中元素的对应关系如下所示：

$χ_{i} (j, k) = χ_{j} (i, k) = χ_{k} (i, j) = x_{i j k}$ (12)

得到PARAFAC模型如下式所示：

$x_{i j k} = \sum_{f = 1}^{F} a_{i f} b_{j f} c_{k f} + e_{i j k}$ (13)

其中 $a_{i f}$ 、 $b_{j f}$ 和 $c_{k f}$ 都是要估计的参数，而 $e_{i j k}$ 残差项， $F$ 表示成分的求和数量。 $χ_{k}$ 表示三阶张量按照第k个模式展开得到 $I \times J$ 的矩阵，那么可以得到PARAFAC模型：

$χ_{k} = A D_{k} B^{T} + E_{k}$ (14)

其中矩阵 $A$ 表示含有元素 $a_{i f}$ 的 $I \times F$ 阶矩阵，同理矩阵 $B$ 是 $J \times F$ 阶矩阵，但是 $D_{k}$ 表示对角矩阵，包含 $K \times F$ 阶矩阵 $C$ 的第k行，对角元素为 $c_{k f}$ 。为了在最小二乘中拟合模型，应该最小化损失函数，如下所示：

$σ = \sum_{k = 1}^{K} {‖ χ_{k} - A D_{k} B^{T} ‖}_{F}^{2}$ (15)

3.2. 模型的唯一性分析

PARAFAC模型具有模糊性问题，即不确定性问题，具体就是排列模糊和尺度模糊问题 [3] 。对于上一小节中的三阶张量 $χ$ 和因子矩阵 $A$ 、 $B$ 和 $C$ ，排列模糊问题 [2] 具体如下式所示：

$χ = 〚 A, B, C 〛 = 〚 A Π, B Π, C Π 〛$ (16)

其中 $〚 \cdot 〛$ 表示多线性乘积， $Π$ 表示 $F \times F$ 列交换矩阵， $A = [a_{1}, a_{2}, \dots, a {}_{F}]$ ， $B = [b_{1}, b_{2}, \dots, b {}_{F}]$ ， $C = [c_{1}, c_{2}, \dots, c {}_{F}]$ 。

尺度模糊问题 [2] 为：当满足 $α_{f} β_{f} γ_{f} = 1, \forall f = 1, \dots, F$ ，即

$χ = \sum_{f = 1}^{F} (α_{f} a_{f}) \circ (β_{f} b_{f}) \circ (γ_{f} c_{f}) = \sum_{f = 1}^{F} a_{f} \circ b_{f} \circ c_{f}$ (17)

针对PARAFAC模型的排列模糊和尺度模糊问题，Kruskal在1977年提出了PARAFAC模型唯一性分解的充分条件 [17] ：

$k_{A} + k_{B} + k_{C} \geq 2 F + 2$ (18)

其中 $F$ 表示三个因子矩阵的列向量数。满足式(17)的时候，张量分解唯一，此时只存在排列模糊和尺度模糊问题。排列模糊和尺度模糊表示如下所示：

$\bar{A} = A Π Δ_{1}, \bar{B} = B Π Δ_{2}, \bar{C} = C Π Δ_{3}$ (19)

其中 $Π \in ℂ^{F \times F}$ 为列交换矩阵， $Δ_{i} \in ℂ^{F \times F} (i = 1, 2, 3)$ 表示对角模糊矩阵，并且满足 $Δ_{1} Δ_{2} Δ_{3} = I_{F}$ 。

在2.1中介绍如果由列向量构成的矩阵 $A \in ℂ^{I \times J}$ ，并且这些列向量是从绝对连续分布中独立提取出的元素，因此矩阵不仅满秩而且满足k-秩。矩阵 $A$ 、 $B$ 和 $C$ 的列向量是从绝对连续分布中独立提取的元素，式(17)就可以写成：

$\min (I, F) + \min (J, F) + \min (K, F) \geq 2 F + 2$ (20)

由式(19)可知：PARAFAC模型可以通过调整因子矩阵的行数来使其分解唯一。在实际应用的场景中就可以改变系统的相关参数来使模型分解唯一。

4. 模型的求解

PARAFAC模型最传统的求解方式是三线性交替最小二乘算法(Trilinear Alternating Least Square, TALS)。它的基本原理是首先更新一个矩阵，然后对其余的矩阵进行最小二乘算法进行更新，如此迭代更新矩阵直到收敛，最终算法结束。由于TALS算法的缺陷是要花费大量的时间来计算，因此本文主要借助一种提高TALS收敛速度的算法来对PARAFAC模型进行求解，称为COMFAC (Complex Parallel Factor Analysis)算法 [7] 。该算法的原理框架如图2所示。

Figure 2. COMFAC algorithm flowchart

图2. COMFAC算法流程图

首先是对数据阵列进行压缩。数据阵为三阶张量 $χ \in ℂ^{I \times J \times K}$ ，分别在这三个维度上新定义三个压缩矩阵 $U (I \times P)$ 、 $V (J \times Q)$ 、 $Z (K \times R)$ ，值得注意的是 $P \leq I, Q \leq J, R \leq K$ 。因此，三阶张量 $χ$ 被压缩成 $G (P \times Q \times R)$ 。在压缩过程中使用的是Tucker3-ALS算法，它的步骤如图3所示。

Figure 3. Tucker3-ALS algorithm flowchart

图3. Tucker3-ALS算法流程图

根据 $G^{(P \times Q R)} = U^{H} X^{(I \times J K)} (Z \otimes V)$ 可知，只需要进行矩阵的逆变换就可以得到原始矩阵 $X^{(I \times J K)}$ ，即

${\hat{X}}^{(I \times J K)} = U G^{(P \times Q R)} (Z^{H} \otimes V^{H})$ (21)

对数据进行压缩之后，接着对矩阵 $G$ 的PARAFAC模型进行初始化，主要是对 $G^{(P \times Q R)}$ 、 $G^{(Q \times R P)}$ 、 $G^{(R \times P Q)}$ 进行奇异值分解。

PARAFAC模型是从三个不同的维度来展开数据阵，经过压缩后的矩阵 $G$ 可以表示为：

$G_{r} = L D_{r} (N) M^{T} + E_{r}, r = 1, \dots, R$ (22)

其中 $L (P \times F)$ 、 $M (Q \times F)$ 、 $N (R \times F)$ 。 $F$ 表示经过奇异值分解之后取得前 $F$ 个左奇异值矢量，即因子的个数。 $E_{r}$ 表示误差矩阵。在PARAFAC模型的框架下，三阶张量可以重构成二维矩阵形式。同样在COMFAC分解中，压缩阵 $G$ 可以表示为：

${\hat{G}}^{(P \times Q R)} = L {(N ⊙ M)}^{T}$ (23)

等价于：

${\hat{G}}^{(P \times Q R)} = L {(N^{*} ⊙ M^{*})}^{H}$ (24)

为了简化计算，忽略矩阵上标，将(24)代入(21)中，得到：

${\hat{X}}^{(I \times J K)} = U G^{(P \times Q R)} (Z^{H} \otimes V^{H}) \approx U L {(N^{*} ⊙ M^{*})}^{H} {(Z \otimes V)}^{H} = U L {[(Z \otimes V) (N^{*} ⊙ M^{*})]}^{H}$ (25)

根据运算性质(10)可得

${\hat{X}}^{(I \times J K)} \approx U L {[(Z \otimes V) (N^{*} ⊙ M^{*})]}^{H} = U L {(Z N^{*} ⊙ V M^{*})}^{H} = U L {(Z^{*} N ⊙ V^{*} M)}^{H}$ (26)

定义 $A \approx U L, B \approx V^{*} M, C \approx Z^{*} N$ ，最终得到

${\hat{X}}^{(I \times J K)} \approx A {(C ⊙ B)}^{T}$ (27)

最后再根据实际应用背景来计算误码率大小。

5. 应用

本文只考虑上行DS-CDMA系统，其中用户数量 $F = 6$ ，传输的数据是采用随机二进制相移键控(Binary Phase Shift Keying, BPSK)符号，扩频码长度为 $K = 4$ ，符号快照的块长度为 $J = 50$ 。在观测接收到张量 $χ \in ℂ^{I \times 50 \times 4}$ 后，首先加入不同信噪比条件下的高斯白噪声，其次进行模糊去除，接着执行COMFAC算法来对DS-CDMA系统进行盲数据检测。它们的性能通过误码率(Bit Error Rate, BER)来体现，这是通过平均1000次蒙特卡洛独立实验得到的平均值。最后，本文在相同的信噪比条件下，得到天线根数范围取4到20对应的误码率，这是由于根据PARAFAC模型的唯一性可知天线根数至少为4，否则误码率的大小脱离实际意义。在此基础上将误码率与天线根数进行数据拟合，取合适的拟合方式，拟合函数中最小的极小值点进行四舍五入得到不同信噪比条件下的最优天线的根数。如图4和图5所示是信噪比1~30条件下误码率与天线根数的拟合曲线。拟合函数如式(28)~(57)，具体见附录。

(a) (b) (c) (d) (e) (f) (g) (h) (i) (j) (k) (l) (m) (n)

Figure 4. Fit the function curve picture. (a) SNR = 1; (b) SNR = 2; (c) SNR = 3; (d) SNR = 4; (e) SNR = 5; (f) SNR = 6; (g) SNR = 7; (h) SNR = 8; (i) SNR = 9; (j) SNR = 10; (k) SNR = 11; (l) SNR = 12; (m) SNR = 13; (n) SNR = 14

图4. 拟合函数曲线图。(a) SNR = 1; (b) SNR = 2; (c) SNR = 3; (d) SNR = 4; (e) SNR = 5; (f) SNR = 6; (g) SNR = 7; (h) SNR = 8; (i) SNR = 9; (j) SNR = 10; (k) SNR = 11; (l) SNR = 12; (m) SNR = 13; (n) SNR = 14

(a) (b) (c) (d) (e) (f) (g) (h) (i) (j) (k) (l) (m) (n) (o) (p)

Figure 5. Fit the function curve picture. (a) SNR = 15; (b) SNR = 16; (c) SNR = 17; (d) SNR = 18; (e) SNR = 19; (f) SNR = 20; (g) SNR = 21; (h) SNR = 22; (i) SNR = 23; (j) SNR = 24; (k) SNR = 25; (l) SNR = 26; (m) SNR = 27; (n) SNR = 28; (o) SNR = 29; (p) SNR = 30

图5. 拟合函数曲线图。(a) SNR = 15; (b) SNR = 16; (c) SNR = 17; (d) SNR = 18; (e) SNR = 19; (f) SNR = 20; (g) SNR = 21; (h) SNR = 22; (i) SNR = 23; (j) SNR = 24; (k) SNR = 25; (l) SNR = 26; (m) SNR = 27; (n) SNR = 28; (o) SNR = 29; (p) SNR = 30

通过观察拟合曲线并且结合拟合方程来求函数的最小的极小值点，结合实际意义进行四舍五入，得到结果如表1所示。通过观察结果可知拟合方式主要有三种：Sum of Sine、Fourier、Polynomial。它们的拟合优度都达到90%以上，高信噪比的最优天线数少，并且比较稳定；低信噪比的最优天线数多，并且随着信噪比变低，它的最优天线数量不稳定；当信噪比为负数的时候，误码率趋于稳定的速度非常慢进而导致研究最优天线数失去实际意义。

Table 1. Optimal antenna statistics table

表1. 最优天线统计表

6. 结论与研究展望

6.1. 总结

本文在DS-CDMA系统中基于张量分解技术建立PARAFAC模型，固定用户数量、扩频码数、符号快照的块长度数，改变信噪比和天线数量，接着执行COMFAC算法得出误码率，将误码率与天线数量进行数据拟合，最后观察拟合函数曲线和拟合方程得出最优天线数。根据蒙特卡洛模拟和数据拟合发现拟合方式主要有三种：Sum of Sine、Fourier、Polynomial，同时拟合优度都在90%以上。此外，最优天线受信噪比的影响，高信噪比的最优天线数少，并且稳定为5根天线；低信噪比的最优天线数多，并且随着信噪比变低，它的最优天线数量不稳定；当信噪比为负数的时候，误码率趋于稳定的速度非常慢，得到的结果也失去研究意义。

6.2. 研究展望

第一，本文并没有考虑负信噪比，因为误码率趋于稳定的速度非常慢，甚至误码率的大小一直脱离实际意义，这是受算法的影响。虽然PARAFAC模型只是张量分解的一种形式，执行COMFAC算法的时候虽然极大提高了运算时间，但是算法精度和收敛速度可以进一步提高。

第二，本文只是取的离散的正整数信噪比，此外信噪比为负数的时候并没有考虑。因此，如果信噪比范围扩大，并且是连续的情况，最优天线的数量可以进一步深入研究。

第三，本文只是在DS-CDMA系统中固定用户数量、扩频码数、符号快照的块长度数来进行仿真模拟，实际生活中的问题更加复杂，误码率和用户、扩频码、天线数量等的其他函数关系可以进一步研究，甚至误码率的其他影响因素也是不容忽视的。

未来的工作聚焦于上述三点，并且模型和算法的优化更是研究的重点。

文章引用

周婧,戴泽安,张志超. 基于张量分解的DS-CDMA盲数据检测研究
Research on Blind Data Detection for DS-CDMA Based on Tensor Decomposition[J]. 建模与仿真, 2023, 12(03): 2976-2993. https://doi.org/10.12677/MOS.2023.123275

参考文献

1. Cichocki, A., Mandic, D., De Lathauwer, L., et al. (2015) Tensor Decompositions for Signal Processing Applications: From Two-Way to Multiway Component Analysis. IEEE Signal Processing Magazine, 32, 145-163. https://doi.org/10.1109/MSP.2013.2297439

2. 丁真苹. 基于张量分解的大规模MIMO半盲信道估计和空时编码技术研究[D]: [硕士学位论文]. 郑州: 郑州大学, 2021.

3. Sidiropoulos, N.D., Giannakis, G.B. and Bro, R. (2000) Blind PARAFAC Receivers for DS-CDMA Systems. IEEE Transactions on Signal Processing, 48, 810-823. https://doi.org/10.1109/78.824675

4. Sørensen, M., De Lathauwer, L. and Deneire, L. (2010) PARAFAC with Or-thogonality in One Mode and Applications in DS-CDMA Systems. 2010 IEEE International Conference on Acoustics, Speech and Signal Processing, Dallas, 14-19 March 2010, 4142-4145. https://doi.org/10.1109/ICASSP.2010.5495717

5. Fernandes, C.A.R., de Almeida, A.L.F. and da Costa, D.B. (2012) Unified Tensor Modeling for Blind Receivers in Multiuser Uplink Cooperative Systems. IEEE Signal Processing Letters, 19, 247-250. https://doi.org/10.1109/LSP.2012.2188510

6. 刘旭, 许宗泽. 正交约束PARAFAC的DS-CDMA盲接收机[J]. 应用科学学报, 2009, 27(2): 131-136.

7. Bro, R., Sidiropoulos, N.D. and Giannakis, G.B. (1999) A Fast Least Squares Algo-rithm for Separating Trilinear Mixtures. Workshop on Independent Component Analysis and Blind Separation, Aussois, 11-15 January 1999, 11-15.

8. De Lathauwer, L. and Castaing, J. (2007) Tensor-Based Techniques for the Blind Separation of DS-CDMA Signals. Signal Processing, 87, 322-336. https://doi.org/10.1016/j.sigpro.2005.12.015

9. 曾卉露. 基于复平行因子的多故障源盲分离方法研究[D]: [硕士学位论文]. 南昌: 南昌航空大学, 2021.

10. 于桂晨. 大规模张量的表示方法及应用研究[D]: [硕士学位论文]. 大连: 大连理工大学, 2020.

11. 张珽. 基于蒙特卡罗方法的通信系统误码率的仿真[J]. 无线通信技术, 2010, 19(1): 20-22+25.

12. 曾璐, 谢晓尧. 基于MATLAB扩频通信系统误码率的研究[J]. 通信技术, 2011, 44(11): 25-26+29.

13. 李开泰, 黄艾香. 张量分析及其应用[M]. 北京: 科学出版社, 2004.

14. 张贤达. 矩阵分析与应用[M]. 北京: 清华大学出版社, 2004.

15. Rao, C.R. (1970) Estimation of Heteroscedastic Variances in Lin-ear Models. Journal of the American Statistical Association, 65, 161-172. https://doi.org/10.1080/01621459.1970.10481070

16. 刘越. 基于张量分解的MIMO系统半盲信道估计算法研究[D]: [硕士学位论文]. 郑州: 郑州大学, 2016.

17. Kruskal, J.B. (1977) Three-Way Arrays: Rank and Uniqueness of Trilinear Decompositions, with Application to Arithmetic Complexity and Statistics. Linear Algebra and Its Applications, 18, 95-138. https://doi.org/10.1016/0024-3795(77)90069-6

附录

当 $S N R = 1$ 的时候，拟合函数：

$\begin{array}{l} f (x) = 2.2 \sin (0 .05336 x + 2 .301) + 0 .566 \sin (0 .2253 x + 2 .961) \\ \begin{matrix} + 0 .05511 \sin (0 .6263 x + 0 .7143) + 0 .01604 \sin (0 .9226 x + 0 .3011) \end{matrix} \end{array}$ (28)

当 $S N R = 2$ 的时候，拟合函数：

$\begin{array}{l} f (x) = 1 .339 \sin (0 .07806 x + 1 .954) + 0 .4957 \sin (0 .2684 x + 2 .444) \\ \begin{matrix} + 0 .1126 \sin (0 .6023 x + 1 .19) + 0 .03626 \sin (0 .7287 x + 2 .934) \end{matrix} \end{array}$ (29)

当 $S N R = 3$ 的时候，拟合函数：

$\begin{array}{l} f (x) = - 4 .996 \times 10^{7} + 6 .089 \times 10^{7} \cos (0 .06478 x) + 6 .092 \times 10^{7} \sin (0 .06478 x) + 3 .208 \times 10^{4} \cos (0.12956 x) \\ \begin{matrix} - 5 .475 \times 10^{7} \sin (0.12956 x) - 1 .771 \times 10^{7} \cos (0.19434 x) + 1 .768 \times 10^{7} \sin (0.19434 x) \end{matrix} \\ \begin{matrix} + 7 .811 \times 10^{6} \cos (0.25912 x) + 1 .221 \times 10^{4} \sin (0.25912 x) - 1 .054 \times 10^{6} \cos (0.3239 x) \end{matrix} \\ \begin{matrix} - 1 .059 \times 10^{6} \sin (0.3239 x) - 446 .2 \cos (0.38868 x) + 1 .326 \times 10^{5} \sin (0.38868 x) \end{matrix} \end{array}$ (30)

当 $S N R = 4$ 的时候，拟合函数：

$\begin{array}{l} f (x) = 0 .4534sin (0 .1701 x + 0.6735) + 0 .3187 \sin (0.3531 x + 1.392) + 0 .09555sin (0 .6825 x + 0 .2507) \\ \begin{matrix} + 0.0409 \sin (1.005 x - 0.4421) - \end{matrix} 0 .01361 \sin (1.005 x + 6.941) \end{array}$ (31)

当 $S N R = 5$ 的时候，拟合函数：

$f (x) = 5.764 \times 10^{5} x^{4} - 0 .003194 x^{3} + 0 .0646 x^{2} - 0 .5679 x + 1.875$ (32)

当 $S N R = 6$ 的时候，拟合函数：

$\begin{array}{l} f (x) = 0 .5778sin (0 .09163 x + 1 .849) + 0 .2456 \sin (0 .2606 x + 2 .728) + 0 .03345sin (0 .7272 x - 0.1952) \\ \begin{matrix} + 0 .01825 \sin (1.113 x - 1.875) + \end{matrix} 0 .01142 \sin (1.515 x + 2.316) \end{array}$ (33)

当 $S N R = 7$ 的时候，拟合函数：

$f (x) = - 6 .33 \times 10^{6} x^{5} + 0 .0004268 x^{4} - 0 .01117 x^{3} + 0 .1417 x^{2} - 0.8734 x + 2.129$ (34)

当 $S N R = 8$ 的时候，拟合函数：

$f (x) = - 5 .929 \times 10^{6} x^{5} + 0 .0003952 x^{4} - 0 .01023 x^{3} + 0 .1284 x^{2} - 0.7818 x + 1.87$ (35)

当 $S N R = 9$ 的时候，拟合函数：

$f (x) = - 5.849 \times 10^{6} x^{5} + 0 .0003867 x^{4} - 0 .009888 x^{3} + 0 .1219 x^{2} - 0 .7239 x + 1.674$ (36)

当 $S N R = 10$ 的时候，拟合函数：

$f (x) = - 3.964 \times 10^{6} x^{5} + 0 .0002663 x^{4} - 0 .006922 x^{3} + 0 .08681 x^{2} - 0 .5246 x + 1.233$ (37)

当 $S N R = 11$ 的时候，拟合函数：

$\begin{array}{l} f (x) = 28.31 + 12.82 \cos (0.1453 x) - 46.96 \sin (0.1453 x) - 26.39 \cos (0.2906 x) \\ \begin{matrix} \begin{matrix} \end{matrix} & - 15.53 \sin (0.2906 x) - 9.681 \cos (0.4359 x) + 9.471 \sin (0.4359 x) + 1.872 \cos (0.5812 x) \end{matrix} \\ \begin{matrix} \begin{matrix} \end{matrix} & + 3.344 \sin (0.5812 x) + 0.5181 \cos (0.7265 x) \end{matrix} - 0.1295 \sin (0.7265 x) \end{array}$ (38)

当 $S N R = 12$ 的时候，拟合函数：

$f (x) = - 2.726 \times 10^{6} x^{5} + 0.0001795 x^{4} - 0.004573 x^{3} + 0.05625 x^{2} - 0.3339 x + 0.7746$ (39)

当 $S N R = 13$ 的时候，拟合函数：

$\begin{array}{l} f (x) = 1.152 \times 10^{5} - 1.365 \times 10^{5} \cos (0.06289 x) - 1.367 \times 10^{5} \sin (0.06289 x) \\ \begin{matrix} - 132.9 \cos (0.12578 x) + 1.124 \times 10^{5} \sin (0.12578 x) + 3.082 \times 10^{4} \cos (0.18867 x) \end{matrix} \\ \begin{matrix} - 3.07 \times 10^{4} \sin (0.18867 x) - 1.01 \times 10^{4} \cos (0.25156 x) - 31.63 \sin (0.25156 x) \end{matrix} \\ \begin{matrix} + 752.6 \cos (0.31445 x) + 759.8 \sin (0.31445 x) \end{matrix} \end{array}$ (40)

当 $S N R = 14$ 的时候，拟合函数：

$\begin{array}{l} f (x) = 6.435 \times 10^{7} - 9.798 \times 10^{7} \cos (0.03191 x) - 4.414 \times 10^{7} \sin (0.03191 x) \\ \begin{matrix} + 4.093 \times 10^{7} \cos (0.06328 x) + 4.626 \times 10^{7} \sin (0.06328 x) - 6.937 \times 10^{6} \cos (0.09573 x) \end{matrix} \\ \begin{matrix} - 2.234 \times 10^{7} \sin (0.09573 x) - 6.423 \times 10^{5} \cos (0.12764 x) + 5.231 \times 10^{6} \sin (0.12764 x) \end{matrix} \\ \begin{matrix} + 2.783 \times 10^{5} \cos (0.15955 x) - 4.586 \times 10^{5} \sin (0.15955 x) \end{matrix} \end{array}$ (41)

当 $S N R = 15$ 的时候，拟合函数：

$f (x) = - 1.887 \times 10^{6} x^{5} + 0.0001242 x^{4} - 0.003164 x^{3} + 0.0389 x^{2} - 0.2308 x + 0.5362$ (42)

当 $S N R = 16$ 的时候，拟合函数：

$f (x) = - 1.642 \times 10^{- 6} x^{5} + 0 .0001063 x^{4} - 0 .002663 x^{3} + 0.03217 x^{2} - 0.1872 x + 0.4248$ (43)

当 $S N R = 17$ 的时候，拟合函数：

$\begin{array}{l} f (x) = 0.02813 \sin (0.1876 x + 0.4112) + 0.02047 \sin (0.3674 x + 1.265) \\ \begin{matrix} + 0.006501 \sin (0.7457 x - 0.3742) + 0.002799 \sin (1.13 x - 2.269) \end{matrix} \\ \begin{matrix} + 0.001412 \sin (1.559 x + 1.482) \end{matrix} \end{array}$ (44)

当 $S N R = 18$ 的时候，拟合函数：

$\begin{array}{l} f (x) = 0.03523 \sin (0.1878 x + 0.4974) + 0.02679 \sin (0.3489 x + 1.561) \\ \begin{matrix} \begin{array}{l} + 0.007204 \sin (0.7176 x + 0.03323) + 0.003803 \sin (1.128 x - 2.043) \\ + 0.002173 \sin (1.52 x + 2.349) \end{array} \end{matrix} \end{array}$ (45)

当 $S N R = 19$ 的时候，拟合函数：

$\begin{array}{l} f (x) = 0.04125 \sin (0.1882 x + 0.4663) + 0.03123 \sin (0.3094 x + 2.039) \\ \begin{matrix} + 0.004423 \sin (0.7486 x - 0.3983) + 0.003346 \sin (1.149 x - 2.535) \end{matrix} \end{array}$ (46)

当 $S N R = 20$ 的时候，拟合函数：

$\begin{array}{l} f (x) = 0.02564 \sin (0.1882 x + 0.4662) + 0.01854 \sin (0.3717 x + 1.237) \\ \begin{matrix} + 0.006017 \sin (0.7566 x - 0.5913) + 0.002827 \sin (1.164 x - 2.798) \end{matrix} \\ \begin{matrix} + 0.00291 \sin (2.754 x - 2.592) \end{matrix} \end{array}$ (47)

当 $S N R = 21$ 的时候，拟合函数：

$\begin{array}{l} f (x) = 0.02624 \sin (0.1879 x + 0.5108) + 0.01972 \sin (0.3379 x + 1.681) \\ \begin{matrix} + 0.004308 \sin (0.7513 x - 0.5372) + 0.002093 \sin (1.125 x - 2.107) \end{matrix} \\ \begin{matrix} + 0.001628 \sin (1.523 x + 1.908) \end{matrix} \end{array}$ (48)

当 $S N R = 22$ 的时候，拟合函数：

$\begin{array}{l} f (x) = 7.571 \times 10^{6} - 8.988 \times 10^{6} \cos (0.06542 x) - 9.463 \times 10^{6} \sin (0.06542 x) - 4.272 \times 10^{5} \cos (0.13084 x) \\ \begin{matrix} + 8.283 \times 10^{6} \sin (0.13084 x) + 2.88 \times 10^{6} \cos (0.19626 x) - 2.465 \times 10^{6} \sin (0.19626 x) \end{matrix} \\ \begin{matrix} - 1.176 \times 10^{6} \cos (0.26168 x) - 1.225 \times 10^{5} \sin (0.26168 x) + 1.378 \times 10^{5} \cos (0.3271 x) \end{matrix} \\ \begin{matrix} + 1.794 \times 10^{5} \sin (0.3271 x) + 3147 \cos (0.39252 x) - 1.981 \times 10^{4} \sin (0.39252 x) \end{matrix} \end{array}$ (49)

当 $S N R = 23$ 的时候，拟合函数：

$\begin{array}{l} f (x) = 1.466 \times 10^{7} - 1.775 \times 10^{7} \cos (0.0651 x) - 1.799 \times 10^{7} \sin (0.0651 x) - 2.138 \times 10^{5} \cos (0.1302 x) \\ \begin{matrix} + 1.605 \times 10^{7} \sin (0.1302 x) + 5.288 \times 10^{6} \cos (0.1953 x) - 5.08 \times 10^{6} \sin (0.1953 x) \end{matrix} \\ \begin{matrix} - 2.285 \times 10^{6} \sin (0.1953 x) - 6.167 \times 10^{4} \sin (0.2604 x) + 2.981 \times 10^{5} \cos (0.3255 x) \end{matrix} \\ \begin{matrix} + 3.191 \times 10^{5} \sin (0.3255 x) + 1599 \cos (0.3906 x) - 3.863 \times 10^{4} \sin (0.3906 x) \end{matrix} \end{array}$ (50)

当 $S N R = 24$ 的时候，拟合函数：

$\begin{array}{l} f (x) = 6.969 \times 10^{6} - 8.415 \times 10^{6} \sin (0.06514 x) - 8.572 \times 10^{6} \sin (0.06514 x) - 1.408 \times 10^{5} \cos (0.13028 x) \\ \begin{matrix} + 7.631 \times 10^{6} \sin (0.13028 x) + 2.533 \times 10^{6} \cos (0.19542 x) - 2.396 \times 10^{6} \sin (0.19542 x) \end{matrix} \\ \begin{matrix} - 1.087 \times 10^{6} \cos (0.26056 x) - 4.07 \times 10^{4} \sin (0.26056 x) + 1.398 \times 10^{5} \cos (0.3257 x) \end{matrix} \\ \begin{matrix} + 1.537 \times 10^{5} \sin (0.3257 x) + 1059 \cos (0.39084 x) - 1.838 \times 10^{4} \sin (0.39084 x) \end{matrix} \end{array}$ (51)

当 $S N R = 25$ 的时候，拟合函数：

$\begin{array}{l} f (x) = 1.937 \times 10^{6} - 2.307 \times 10^{6} \cos (0.06545 x) - 2.414 \times 10^{6} \sin (0.06545 x) - 9.6014 \times 10^{4} \cos (0.1309 x) \\ \begin{matrix} + 2.118 \times 10^{6} \sin (0.1309 x) + 7.291 \times 10^{5} \cos (0.19635 x) - 6.359 \times 10^{5} \sin (0.19635 x) \end{matrix} \\ \begin{matrix} - 3 \times 10^{5} \cos (0.2618 x) - 2.746 \times 10^{4} \sin (0.2618 x) + 3.57 \times 10^{4} \cos (0.32725 x) \end{matrix} \\ \begin{matrix} + 4.5 \times 10^{4} \sin (0.32725 x) + 701.5 \cos (0.3927 x) - 5026 \sin (0.3927 x) \end{matrix} \end{array}$ (52)

当 $S N R = 26$ 的时候，拟合函数：

$\begin{array}{l} f (x) = 6.829 \times 10^{6} - 8.184 \times 10^{6} \cos (0.0649 x) - 8.458 \times 10^{6} \sin (0.0649 x) - 2.469 \times 10^{5} \cos (0.1298 x) \\ \begin{matrix} + 7.47 \times 10^{6} \sin (0.1298 x) + 2.529 \times 10^{6} \cos (0.1947 x) - 2.289 \times 10^{6} \sin (0.1947 x) \end{matrix} \\ \begin{matrix} - 1.06 \times 10^{6} \sin (0.2596 x) - 7.096 \times 10^{4} \sin (0.2596 x) + 1.306 \times 10^{5} \cos (0.3245 x) \end{matrix} \\ \begin{matrix} + 1.548 \times 10^{5} \sin (0.3245 x) + 1828 \cos (0.3894 x) - 1.78 \times 10^{4} \sin (0.3894 x) \end{matrix} \end{array}$ (53)

当 $S N R = 27$ 的时候，拟合函数：

$\begin{array}{l} f (x) = 4.4 \times 10^{6} - 5.23 \times 10^{6} \cos (0.065 x) - 5.489 \times 10^{6} \sin (0.065 x) - 2.332 \times 10^{5} \cos (0.13 x) \\ \begin{matrix} + 4.809 \times 10^{6} \sin (0.13 x) + 1.662 \times 10^{6} \cos (0.195 x) - 1.436 \times 10^{6} \sin (0.195 x) \end{matrix} \\ \begin{matrix} - 6.805 \times 10^{5} \cos (0.26 x) - 6.679 \times 10^{4} \sin (0.26 x) + 8.016 \times 10^{4} \cos (0.325 x) \end{matrix} \\ \begin{matrix} + 1.028 \times 10^{5} \sin (0.325 x) + 1711 \cos (0.39 x) - 1.139 \times 10^{4} \sin (0.39 x) \end{matrix} \end{array}$ (54)

当 $S N R = 28$ 的时候，拟合函数：

$\begin{array}{l} f (x) = 9.601 \times 10^{6} - 1.146 \times 10^{7} \cos (0.06537 x) - 1.193 \times 10^{7} \sin (0.06537 x) - 425600 \cos (0.13074 x) \\ \begin{matrix} + 1.05 \times 10^{7} \sin (0.13074 x) + 3.595 \times 10^{6} \cos (0.19611 x) - 3.181 \times 10^{6} \sin (0.19611 x) \end{matrix} \\ \begin{matrix} - 1.491 \times 10^{6} \cos (0.26148 x) - 1.22 \times 10^{5} \sin (0.26148 x) + 1.802 \times 10^{5} \cos (0.32685 x) \end{matrix} \\ \begin{matrix} + 2.216 \times 10^{5} \sin (0.32685 x) \end{matrix} + 3131 \cos (0.39222 x) - 2.508 \times 10^{4} \sin (0.39222 x) \end{array}$ (55)

当 $S N R = 29$ 的时候，拟合函数：

$\begin{array}{l} f (x) = 5.7 \times 10^{6} - 6.854 \times 10^{6} \cos (0.06525 x) - 7.038 \times 10^{6} \sin (0.06525 x) - 1.662 \times 10^{5} \cos (0.1305 x) \\ \begin{matrix} + 6.237 \times 10^{6} \sin (0.1305 x) + 2.093 \times 10^{6} \cos (0.19575 x) - 1.931 \times 10^{6} \sin (0.19575 x) \end{matrix} \\ \begin{matrix} - 8.856 \times 10^{5} \cos (0.261 x) - 4.767 \times 10^{4} \sin (0.261 x) + 1.113 \times 10^{5} \cos (0.32625 x) \end{matrix} \\ \begin{matrix} + 1.275 \times 10^{5} \sin (0.32625 x) + 1224 \cos (0.3915 x) - 1.491 \times 10^{4} \sin (0.3915 x) \end{matrix} \end{array}$ (56)

当 $S N R = 30$ 的时候，拟合函数：

$\begin{array}{l} f (x) = 3.044 \times 10^{6} - 3.592 \times 10^{6} \cos (0.0654 x) - 3.825 \times 10^{6} \sin (0.0654 x) - 2.09 \times 10^{5} \cos (0.1308 x) \\ \begin{matrix} + 3.327 \times 10^{6} \sin (0.1308 x) + 1.174 \times 10^{6} \cos (0.1962 x) - 9.705 \times 10^{5} \sin (0.1962 x) \end{matrix} \\ \begin{matrix} - 4.71 \times 10^{5} \cos (0.2616 x) - 5.994 \times 10^{4} \sin (0.2616 x) + 5.316 \times 10^{4} \cos (0.327 x) \end{matrix} \\ \begin{matrix} + 7.351 \times 10^{4} \sin (0.327 x) + 1539 \cos (0.3924 x) - 7885 \sin (0.3924 x) \end{matrix} \end{array}$ (57)

期刊菜单