颜佳佳论文
多智能体随机网络结构的实时精确估算
基于扰动理论的特征向量估算方法
设原矩阵为 $A$,扰动后矩阵为 $A+\zeta C$(扰动矩阵 $\zeta C$,$\zeta$是小参数),令其第 $i$ 个特征值、特征向量分别为 $\lambda_i,x_i$ 和 $\tilde\lambda_i,\tilde x_i$。
特征向量的一阶扰动公式:
$$ \Delta x_i =\tilde x_i - x_i \;\approx\; \zeta \sum_{k\neq i} \frac{x_k^T\,C\,x_i}{\lambda_i - \lambda_k}\;x_k, $$- 输出:对应第 $i$ 个特征向量修正量 $\Delta x_i$。
特征值的一阶扰动公式:
$$ \Delta\lambda_i = \tilde\lambda_i - \lambda_i \;\approx\;\zeta\,x_i^T\,C\,x_i $$ **关键假设:**当扰动较小( $\zeta\ll1$) 且各模态近似正交均匀时,常作进一步近似 $$ x_k^T\,C\,x_i \;\approx\; x_i^T\,C\,x_i \; $$ 正交: $\{x_k\}$ 本身是正交基,这是任何对称矩阵特征向量天然具有的属性。均匀:我们把 $C$ 看作“不偏向任何特定模态”的随机小扰动——换句话说,投影到任何两个方向 $(x_i,x_k)$ 上的耦合强度 $x_k^T,C,x_i\quad\text{和}\quad x_i^T,C,x_i$ 在数值量级上应当差不多,因此可以互相近似。
因此,将所有的 $x_k^T C x_i$ 替换为 $x_i^T C x_i$:
$$ \Delta x_i \approx \zeta \sum_{k\neq i} \frac{x_i^T C x_i}{\lambda_i - \lambda_k} x_k = \zeta (x_i^T C x_i) \sum_{k\neq i} \frac{1}{\lambda_i - \lambda_k} x_k = \sum_{k\neq i} \frac{\Delta \lambda_i}{\lambda_i - \lambda_k} x_k \tag{*} $$ $$ \Delta x_i \approx\sum_{k\neq i} \frac{\Delta \lambda_i}{\lambda_i - \lambda_k} x_k \tag{*} $$问题:
-
当前时刻的邻接矩阵
$$ A^{(1)}\in\mathbb R^{n\times n},\qquad A^{(1)},x_i^{(1)}=\lambda_i^{(1)},x_i^{(1)},\quad |x_i^{(1)}|=1. $$ -
下一时刻的邻接矩阵
$$ A^{(2)}\in\mathbb R^{n\times n}, $$ 已知它的第 $i$ 个特征值 $\lambda_i^{(2)}$(卡尔曼滤波得来). 求当前时刻的特征向量 $x_i^{(2)}$。
下一时刻第 $i$ 个特征向量的预测为
$$ \boxed{ x_i^{(2)} \;=\; x_i^{(1)}+\Delta x_i \;\approx\; x_i^{(1)} +\sum_{k\neq i} \frac{\lambda_i^{(2)}-\lambda_i^{(1)}} {\lambda_i^{(1)}-\lambda_k^{(1)}}\; x_k^{(1)}. } $$ 通过该估算方法可以依次求出下一时刻的所有特征向量。多智能体随机网络特征值滤波建模
1. 状态转移模型
系统的特征值向量 $\lambda_k$(状态向量)随时间演化的动态方程为:
$$ \lambda_k = \lambda_{k-1} + w_{k-1} $$- 参数说明:
- $\lambda_k \in \mathbb{R}^{r \times 1}$:$k$ 时刻的特征值向量,$r$ 为特征值个数。
- $w_{k-1} \sim \mathcal{N}(0, {Q})$:过程噪声,均值为零,协方差矩阵为对角阵 $\mathbf{Q} \in \mathbb{R}^{r \times r}$(因特征值独立)。
- 简化假设:
- 状态转移矩阵 $\mathbf{A}$ 和控制输入矩阵 $\mathbf{B}$ 为单位阵或零(无外部控制输入),故模型简化为随机游走形式。
2. 测量模型
$v_k$
观测到的特征值向量 $z_k$ 为:
$$ z_k = \lambda_k + v_k $$- 参数说明:
- $z_k \in \mathbb{R}^{r \times 1}$:观测向量,维度与状态向量相同。
- $v_k \sim \mathcal{N}(0, \mathbf{R})$:测量噪声,协方差 $\mathbf{R} \in \mathbb{R}^{r \times r}$ 为对角阵(噪声独立)。
- 简化假设:
- 测量矩阵 $\mathbf{H}$ 为单位阵,即观测直接反映状态。
3. 噪声协方差矩阵的设定
- $\mathbf{Q}$ 和 $\mathbf{R}$ 为对角矩阵,对角元素由特征值方差确定: $$ \mathbf{Q} = \text{diag}(2\sigma_1^2, 2\sigma_2^2, \dots, 2\sigma_r^2), \quad \mathbf{R} = \text{diag}(\sigma_1^2, \sigma_2^2, \dots, \sigma_r^2) $$ 其中 $\sigma_i^2$ 为第 $i$ 个特征值的初始方差(由引理3-1推导)。
网络结构优化
直接SNMF分解(无优化)
- 输入矩阵:原始动态网络邻接矩阵 $A$(可能稠密或高秩)
- 处理流程:
- 直接对称非负矩阵分解:$A \approx UU^T$
- 通过迭代调整$U$和旋转矩阵$Q$逼近目标
- 存在问题:
- 高秩矩阵需要保留更多特征值($\kappa$较大)
- 非稀疏矩阵计算效率低
先优化再SNMF(论文方法)
- 优化阶段(ADMM):
- 目标函数:$\min_{A_{\text{opt}}} (1-\alpha)|A_{\text{opt}}|* + \alpha|A{\text{opt}}|_1$
- 输出优化矩阵$A_{\text{opt}}$
- SNMF阶段:
- 输入变为优化后的$A_{\text{opt}}$
- 保持相同分解流程但效率更高
网络优化中的邻接矩阵重构问题建模与优化
可行解集合定义(公式4-4)
$$ \Omega = \left\{ A \middle| A^T = A,\, A \odot P = A_{\text{pre}} \odot P,\, A \odot A_{\max}' = 0 \right\} $$- $A^T = A$:确保邻接矩阵对称
- $A \odot P = A_{\text{pre}} \odot P$:掩码矩阵$P$ ,原来已有的连接不变,只让优化原本没有连接的地方。($\odot$为Hadamard积)
- $A \odot A_{\max}' = 0$:功率约束矩阵$\ A_{\max}'$ 禁止在原本无连接的节点间新增边
限制矩阵$A_{\max}'$的定义(公式4-5)
$$ A'_{\max, ij} = \begin{cases} 0, & \text{若 } A_{\max, ij} \ne 0 \\ 1, & \text{若 } A_{\max, ij} = 0 \end{cases} $$- $A_{\max, ij}$表示在最大发射功率下哪些节点对之间能连通(非零表示可连通,零表示即便满功率也连不通)
- $A_{\max}'$在“连不通”的位置上是1,其他位置是0。通过$A_{\max}'$标记禁止修改的零元素位置
- 对于所有满足 $A'{\max,ij}=1$ 的位置(即物理不可连通的节点对),必须有 $A{ij}=0$,即始终保持断开
原始优化目标(公式4-6)
$$ \min_{A} \, (1-\alpha)\, \text{rank}(A) + \alpha \|A\|_0 $$ $\|A\|_0$ 表示矩阵 $A$ 中非零元素的个数- 目标:平衡低秩性($\text{rank}(A)$)与稀疏性($|A|_0$)
- 问题:非凸、不可导,难以直接优化
凸松弛后的目标(公式4-7)
$$ \min_{A} \, (1-\alpha)\, \|A\|_* + \alpha \|A\|_1 $$- 核范数$|A|_*$:奇异值之和,替代$\text{rank}(A)$
- L1范数$|A|_1$:元素绝对值和,替代$|A|_0$
- 性质:凸优化问题,存在全局最优解
求解方法
- 传统方法: 可转化为**半定规划(SDP)**问题,使用内点法等求解器。但缺点是计算效率低,尤其当矩阵规模大(如多智能体网络节点数 $n$ 很大)时不可行。
- 改进方法: 采用ADMM(交替方向乘子法)结合投影和对偶上升的方法,适用于动态网络(矩阵频繁变化的情况)。
ADMM核心算法
变量定义与作用
- 输入变量:
- $A_{pre}$:初始邻接矩阵(优化前的网络拓扑)。
- $P$:对称的0-1矩阵,用于标记 $A_{pre}$ 中非零元素的位置(保持已有边不变)。
- $A'{max}$:功率最大时的邻接矩阵的补集($A'{maxij} = 1$ 表示 $A_{maxij} = 0$,即不允许新增边)。
- $\alpha$:权衡稀疏性($L_1$ 范数)和低秩性(核范数)的系数。
- iters:ADMM迭代次数。

算法步骤详解
(S.1) 更新原始变量 $A$(对应ADMM的$x$步)
- 代码行4-17:通过内层循环(投影和对偶上升)更新 $A$。
- 行4-11:
- 通过内层循环(行8-11)迭代更新 $R$,本质是梯度投影法:
- $temp_R^{k+1} = M - X^k \odot A_{\text{pre}}$(计算残差)。
- $X^{k+1} = X^k + \beta(A_{\text{pre}} \odot temp_R^{k+1})$(梯度上升步,$\beta$ 为步长)。
- 本质:通过迭代强制 $A$ 在 $P$ 标记的位置与 $A_{pre}$ 一致。
- 通过内层循环(行8-11)迭代更新 $R$,本质是梯度投影法:
- 行13-17:将 $A$ 投影到 $A \odot A'_{\text{max}} = 0$ 的集合。
- 类似地,通过内层循环(行14-17)更新 $Y$:
- $temp_A^{k+1} = D_2 - Y^k \odot A'_{\text{max}}$(残差计算)。
- $Y^{k+1} = Y^k + \gamma(A'_{\text{max}} \odot temp_A^{k+1})$(对偶变量更新)。
- 类似地,通过内层循环(行14-17)更新 $Y$:
- 行4-11:
(S.2) 更新辅助变量 $Z_1, Z_2$(对应ADMM的$z$步)
通过阈值操作分离目标函数的两部分:
- 行18-19:分别对核范数和 $L_1$ 范数进行阈值操作:
- $Z_1^{t+1} = T_r(A^{t+1} + U_1^t)$:
$T_r(\cdot)$ 是奇异值阈值算子(核范数投影),对$A + U1$ 做SVD分解,保留前 $r$ 个奇异值。作用:把自己变成低秩矩阵=》强制 $A$ 低秩。 - $Z_2^{t+1} = S_{\alpha}(A^{t+1} + U_2^t)$:
$S_{\alpha}(\cdot)$ 是软阈值算子($L_1$ 范数投影),将小于 $\alpha$ 的元素置零。把自己变成稀疏矩阵=》促进 $A$ 的稀疏性。
- $Z_1^{t+1} = T_r(A^{t+1} + U_1^t)$:
(S.3) 更新 拉格朗日乘子$U_1, U_2$(对应ADMM的对偶上升)
- 行20-21:通过残差 $(A - Z)$ 调整拉格朗日乘子 $U_1, U_2$:
- $U_1^{t+1} = U_1^t + A^{t+1} - Z_1^{t+1}$(核范数约束的乘子更新)。
- $U_2^{t+1} = U_2^t + A^{t+1} - Z_2^{t+1}$($L_1$ 范数约束的乘子更新)。
- 作用:惩罚 $A$ 与辅助变量 $Z1, Z2$ 的偏差(迫使$A$更贴近$Z$),推动收敛。
网络结构控制
- 核心目标:将优化后的低秩稀疏矩阵 $A$ 转化为实际网络参数(如功率、带宽),并维持动态网络的连通性和稳定性。
- 具体实现:
- 通过PID控制调整发射/接收功率,使实际链路带宽匹配矩阵 $A$ 的优化值。
- 结合CSMA/CA协议处理多节点竞争,确保稀疏网络下的高效通信。
优化模型(4.2节)
↓ 生成目标带宽矩阵A
香农公式 → 计算目标Pr → 自由空间公式 → 计算目标Pt
↓
PID控制发射机(AGC电压) → 实际Pt ≈ 目标Pt
↓
PID控制接收机(AAGC/DAGC) → 实际Pr ≈ 目标Pr
↓
实际带宽 ≈ Aij (闭环反馈)
-
发射机:
- 功能:将数据转换为无线信号并通过天线发射。
- 关键参数:发射功率($P_t$)、天线增益($G_t$)、工作频率(决定波长$\lambda$)。
- 控制目标:通过调整AGC电压,动态调节发射功率,以匹配优化后的带宽需求(矩阵$A$中的$A_{ij}$)。
-
接收机:
- 功能:接收无线信号并转换为可处理的数据。
- 关键参数:接收功率($P_r$)、噪声($N_0$)、天线增益($G_r$)。
- 控制目标:通过AAGC/DAGC增益调整,确保接收信号强度适合解调,维持链路稳定性。
具体步骤
步骤1:生成目标带宽矩阵 $A$(4.2节优化模型)
-
数学建模:
-
通过凸松弛优化问题(公式4-7)得到低秩稀疏矩阵 $A$: $$ \min_A (1-\alpha) |A|_* + \alpha |A|_1 \quad \text{s.t.} \quad A \in \Omega $$
-
约束集 $\Omega$ 确保矩阵对称性、保留原有链路($A \odot P = A_{\text{pre}} \odot P$)、禁止不可达链路($A \odot A'_{\max} = 0$)。
-
-
物理意义:
- 非零元素 $A_{ij}$ 直接表示 目标信道带宽 $C_{ij}$(单位:bps),即: $$ A_{ij} = C_{ij} = W \log_2\left(1 + \frac{P_r}{N_0 W}\right) \quad \text{(香农公式4-10)} $$
步骤2:从带宽 $A_{ij}$ 反推功率参数
-
接收功率 $P_r$ 计算:
-
根据香农公式解耦: $$ P_r = (2^{A_{ij}/W} - 1) N_0 W $$
-
输入:噪声 $N_0$、带宽 $W$、目标带宽 $A_{ij}$。
-
-
发射功率 $P_t$ 计算:
-
通过自由空间公式(4-11): $$ P_t = \frac{P_r L (4\pi d)^2}{G_t G_r \lambda^2} $$
-
输入:距离 $d$、天线增益 $G_t, G_r$、波长 $\lambda$、损耗 $L$。
-
-
逻辑分支:
- 若 $A_{ij} \neq A_{\text{pre}ij}$(需调整链路):
- 计算 $P_r$ 和 $P_t$;
- 若 $A_{ij} = 0$(无连接):
- 直接设 $P_r = P_t = 0$。
- 若 $A_{ij} \neq A_{\text{pre}ij}$(需调整链路):
步骤3:发射机功率调整(图4-2a)
- 定义目标:$P_t$(来自步骤2)。
- 测量实际:通过传感器获取当前发射功率 $P_{t,\text{actual}}$。
- 计算偏差:$e(t) = P_t - P_{t,\text{actual}}$。
- PID调节:通过AGC电压改变发射功率,逼近 $P_t$。
步骤4:接收机功率调整(图4-2b)
- 定义目标:$P_r$(来自步骤2)。
- 测量实际:检测空口信号功率 $P_{r,\text{actual}}$。
- 计算偏差:$e(t) = P_r - P_{r,\text{actual}}$。
- PID调节:
- 调整AAGC(模拟增益)和DAGC(数字增益),持续监测直至 $|e(t)| < \epsilon$。
基于谱聚类的无人机网络充电
(1) 谱聚类分组Spectral_Clustering(表5.1)
目标:将无人机和充电站划分为 $K$ 个簇,使充电站位于簇中心。
步骤:
- 输入:带权邻接矩阵 $A$(权值=无人机间距离)、节点数 $N$、充电站数 $K$。
- 拉普拉斯矩阵:
$$L = D - A, \quad D_{ii} = \sum_j A_{ij}$$ - 归一化:
$$L_{norm} = D^{-\frac{1}{2}}LD^{\frac{1}{2}}$$ - 谱分解:求 $L_{norm}$ 前 $K$ 小特征值对应的特征向量矩阵 $V \in \mathbb{R}^{N \times K}$。
- 聚类:对 $V$ 的行向量进行 k-means 聚类,得到标签 $\text{labels}$。
输出:每个无人机/充电站的簇编号 $\text{labels}$。
(2) 无人机选择充电站(表5-2)
目标:电量低的无人机前往对应簇中心的充电站。
步骤:
- 周期性运行(间隔 $\Delta t$):
- 通过
Push_Sum
协议获取所有无人机位置Positions
。 - 计算距离矩阵 $A$。
- 通过
- 动态聚类:调用
Spectral_Clustering(A)
更新簇标签。 - 充电触发:若电量 $E < P_{th}$,向簇中心请求坐标 $\text{CS_point}$ 并前往。
关键公式:
$$A_{ij} = \| \text{Position}_i - \text{Position}_j \|_2$$(3) 充电站跟踪算法(表5-3)
目标:充电站动态调整位置至簇中心。
步骤:
- 周期性运行(间隔 $\Delta t$):
- 同无人机算法获取 $A$ 和
labels
。
- 同无人机算法获取 $A$ 和
- 定位簇中心:
- 充电站根据编号匹配簇标签。
- 计算簇内无人机位置均值:
$$\text{CS_point} = \frac{1}{|C_k|} \sum_{i \in C_k} \text{Position}_i$$
其中 $C_k$ 为第 $k$ 簇的节点集合。
- 移动至新中心并广播位置。
(4) 算法改进
替换通信协议:用第3章的卡尔曼滤波 替代 Push_Sum
,获取特征值、特征向量重构全局矩阵 $A$,减少消息传递。
基于T-GAT的无人机群流量预测
TCN
流量矩阵 $X \in \mathbb{R}^{N \times T}$,其中:
- $N$:无人机节点数量(例如10架无人机)。
- $T$:时间步数量。
- 每个元素 $X_{i,t}$ 表示第 $i$ 个节点在时间 $t$ 的总流量(如发送/接收的数据包数量或带宽占用)。
流量矩阵的形状
假设有3架无人机,记录5个时间步的流量数据,矩阵如下:
$$ X = \begin{bmatrix} 100 & 150 & 200 & 180 & 220 \\[6pt] 50 & 75 & 100 & 90 & 110 \\[6pt] 80 & 120 & 160 & 140 & 170 \end{bmatrix} $$- 行 ($N=3$):每行代表一架无人机的历史流量序列(例如第1行表示无人机1的流量变化:100 → 150 → 200 → 180 → 220)。
- 列 ($T=5$):每列代表所有无人机在同一时间步的流量状态(例如第1列表示在时间 $t_1$ 时,三架无人机的流量分别为:[100, 50, 80])。
TCN处理流量矩阵:
- 卷积操作 TCN 的每个卷积核会滑动扫描所有通道(即所有无人机)的时序数据。 例如,一个大小为 3 的卷积核会同时分析每架无人机连续 3 个时间步的流量(例如从 $t_1$ 到 $t_3$),以提取局部时序模式。
- 输出时序特征
经过多层扩张卷积和残差连接后,TCN 会输出一个高阶特征矩阵 $H_T^l$,其形状与输入类似(例如
(1, 3, 5)
),但每个时间步的值已包含了:- 趋势信息:流量上升或下降的长期规律。
TCN的卷积核仅在单个通道内滑动,计算时仅依赖该节点自身的历史时间步。节点间的交互是通过后续的**图注意力网络(GAT)**实现的。
与 GAT 的衔接
- TCN 输出的特征矩阵 $H_T^l$ 会传递给 GAT 进行进一步处理。
- 时间步对齐:通常取最后一个时间步的特征(例如
H_T^l[:, :, -1]
)作为当前节点特征。 - 空间聚合:GAT 根据邻接矩阵计算无人机间的注意力权重,例如考虑“无人机2的当前流量可能受到无人机1过去3分钟流量变化的影响”。
LSTM+GAT训练过程说明(RWP网络节点移动预测)
先LSTM后GAT侧重点在时序特征的提取;先GAT后LSTM侧重点在空间特征的提取。
1. 数据构造
输入数据:
- 节点轨迹数据:
- 每个节点在1000个时间单位内的二维坐标 $(x, y)$,形状为 $[N, 1000, 2]$($N$个节点,1000个时间步,2维特征)。
- 动态邻接矩阵序列:
- 每个时间步的节点连接关系(基于距离阈值或其他规则生成),得到1000个邻接矩阵 $[A_1, A_2, \dots, A_{1000}]$,每个 $A_t$ 的形状为 $[2, \text{num_edges}_t]$(稀疏表示)。
滑动窗口处理:
- 窗口大小:12(用前12个时间步预测第13个时间步)。
- 滑动步长:1(每次滑动1个时间步,生成更多训练样本)。
- 生成样本数量:
- 总时间步1000,窗口大小12 → 可生成 $1000 - 12 = 988$ 个样本。
样本格式:
- 输入序列 $X^{(i)}$:形状 $[N, 12, 2]$($N$个节点,12个时间步,2维坐标)。
- 目标输出 $Y^{(i)}$:形状 $[N, 2]$(第13个时间步所有节点的坐标)。
- 动态邻接矩阵:每个样本对应12个邻接矩阵 $[A^{(i)}_1, A^{(i)}2, \dots, A^{(i)}{12}]$(每个 $A^{(i)}_t$ 形状 $[2, \text{num_edges}_t]$)。
2. 训练过程
模型结构:
-
LSTM层:
- 输入:$[N, 12, 2]$($N$个节点的12步历史轨迹)。
- 输出:每个节点的时序特征 $[N, 12, H]$($H$为LSTM隐藏层维度)。
- 关键点:LSTM独立处理每个节点的时序,节点间无交互。
-
GAT层:
- 输入:取LSTM最后一个时间步的输出 $[N, H]$(即每个节点的最终时序特征)。
- 动态图输入:使用第12个时间步的邻接矩阵 $A^{(i)}_{12}$(形状 $[2, \text{num_edges}]$)。
- 输出:通过图注意力聚合邻居信息,得到空间增强的特征 $[N, H']$($H'$为GAT输出维度)。
-
预测层:
- 全连接层将 $[N, H']$ 映射到 $[N, 2]$,预测下一时刻的坐标。
训练步骤:
- 前向传播:
- 输入 $[N, 12, 2]$ → LSTM → $[N, 12, H]$ → 取最后时间步 $[N, H]$ → GAT → $[N, H']$ → 预测 $[N, 2]$。
- 损失计算:
- 均方误差(MSE)损失:比较预测坐标 $[N, 2]$ 和真实坐标 $[N, 2]$。
- 反向传播:
- 梯度从预测层回传到GAT和LSTM,更新所有参数。
3. 数据维度变化总结
步骤 | 数据形状 | 说明 |
---|---|---|
原始输入 | $[N, 1000, 2]$ | $N$个节点,1000个时间步的$(x,y)$坐标。 |
滑动窗口样本 | $[N, 12, 2]$ | 每个样本包含12个历史时间步。 |
LSTM输入 | $[N, 12, 2]$ | 输入LSTM的节点独立时序数据。 |
LSTM输出 | $[N, 12, H]$ | $H$为LSTM隐藏层维度。 |
GAT输入(最后时间步) | $[N, H]$ | 提取每个节点的最终时序特征。 |
GAT输出 | $[N, H']$ | $H'$为GAT输出维度,含邻居聚合信息。 |
预测输出 | $[N, 2]$ | 下一时刻的$(x,y)$坐标预测。 |
4. 关键注意事项
-
动态图的处理:
- 每个滑动窗口样本需匹配对应时间步的邻接矩阵(如第 $i$ 到 $i+11$ 步的 $[A^{(i)}1, \dots, A^{(i)}{12}]$),但GAT仅使用最后一步 $A^{(i)}_{12}$。
- 若图结构变化缓慢,可简化为所有窗口共享 $A^{(i)}_{12}$。
-
数据划分:
- 按时间划分训练/验证集(如前800个窗口训练,后188个验证),避免未来信息泄露。