md_files/科研/高飞论文.md

17 KiB
Raw Blame History

高飞论文

证明特征值序列为平稳的时间序列

问题设定

  • 研究对象
    \{\lambda_1(A)_t\}_{t\in\mathbb Z} 是随时间变化的随机对称矩阵 A_t 的最大特征值序列(如动态网络的邻接矩阵)。
  • 目标
    证明 \{\lambda_1(A)_t\}二阶(弱)平稳的时间序列,即
    1. $E[\lambda_1(A)_t]=\mu_1$(与 t 无关);
    2. $\operatorname{Var}[\lambda_1(A)_t]=\sigma_1^2<\infty$(与 t 无关);
    3. \operatorname{Cov}(\lambda_1(A)_t,\lambda_1(A)_{t-k})=\gamma(k) 只依赖滞后 $k$。

关键假设

  • 矩阵统计特性(引理 1

    • A_tN\times N 实对称随机矩阵;元素 \{a_{ij}\}_{i\le j} 相互独立且有界:$|a_{ij}|\le K$。

    • 非对角元素:$E[a_{ij}]=\mu>0,\ \operatorname{Var}(a_{ij})=\sigma^2$;对角元素:$E[a_{ii}]=v$。

    • N 足够大时

      
      E[\lambda_1(A_t)]\approx(N-1)\mu+v+\tfrac{\sigma^2}{\mu}\equiv\mu_1,\qquad
      \operatorname{Var}[\lambda_1(A_t)]\approx2\sigma^2\equiv\sigma_1^2 .
      

说明:

  • $\sigma^2$
    这是随机矩阵 A_t 的非对角线元素 a_{ij} (i \neq j) 的方差,即

    
    \text{Var}(a_{ij}) = \sigma^2.
    

    根据引理1的假设所有非对角线元素独立同分布均值为 $\mu$,方差为 $\sigma^2$。

  • $\sigma_1^2$
    这是最大特征值 \lambda_1(A_t) 的方差,即

    
    \text{Var}[\lambda_1(A_t)] \equiv \sigma_1^2.
    

    N 足够大时,\sigma_1^2 近似为 $2\sigma^2$。

  • 时间序列模型
    去中心化序列

    
      \tilde z_t:=\lambda_1(A)_t-\mu_1
    

    假设其服从 AR(1)

    
      \tilde z_t=\rho\,\tilde z_{t-1}+\varepsilon_t,\qquad
      \varepsilon_t\stackrel{\text{i.i.d.}}{\sim}\text{WN}(0,\sigma_\varepsilon^{2}),\ \ |\rho|<1,
    

    \varepsilon_t 与历史 \{\tilde z_{s}\}_{s<t} 独立。

证明主特征值序列平稳

(1) 均值恒定性的推导

  • 去中心化后 $E[\tilde z_t]=0$。因此
    
      E[\lambda_1(A)_t]=E[\tilde z_t]+\mu_1=\mu_1,
    
    t 无关,满足第一条。

(2) 方差恒定

AR(1)模型定义为:


z_t = \rho z_{t-1} + \varepsilon_t

\begin{aligned}
z_t &= \varepsilon_t + \rho z_{t-1} \\
&= \varepsilon_t + \rho (\varepsilon_{t-1} + \rho z_{t-2}) \\
&= \varepsilon_t + \rho \varepsilon_{t-1} + \rho^2 \varepsilon_{t-2} + \cdots \\
&= \sum_{j=0}^\infty \rho^j \varepsilon_{t-j}
\end{aligned}

\text{Var}(z_t) = \text{Var}\left( \sum_{j=0}^\infty \rho^j \varepsilon_{t-j} \right)= \sum_{j=0}^\infty \rho^{2j} \text{Var}(\varepsilon_{t-j})

由于\text{Var}(\varepsilon_{t-j}) = \sigma_\varepsilon^2 对所有 j 成立,


= \sigma_\varepsilon^2 \sum_{j=0}^\infty \rho^{2j}=\frac{\sigma_\varepsilon^2}{1-\rho^2}
  • $|\rho| < 1$ 是保证级数收敛和方差有限的充要条件。

根据引理1$\text{Var}[\lambda_1(A_t)] \approx 2\sigma^2 = \sigma_1^2$。为使模型与理论一致,可设:


\sigma_\varepsilon^2 = (1 - \rho^2) \cdot 2\sigma^2

此时:


\text{Var}[\tilde{z}_t] = 2\sigma^2 = \sigma_1^2

(3) 协方差仅依赖滞后 $k$

对 $k\ge0$


  \gamma(k):=\operatorname{Cov}(\tilde z_t,\tilde z_{t-k})
            =\rho^{k}\sigma_{\tilde z}^{2},

仅含 k 而与 t 无关;于是


  \operatorname{Cov}(\lambda_1(A)_t,\lambda_1(A)_{t-k})=\gamma(k),

满足第三条。

(4) 平稳性的核心条件

  1. |ρ| < 1 是关键条件
    • 直观上:\rho 越小,当前特征值对过去的依赖越弱;
    • \rho=\pm1 会让方差发散,不可能稳态。
  2. 噪声独立性\varepsilon_t 为白噪声,确保新信息与历史无关。

证明剩余特征值平稳(大模型说不可取):

1. 收缩操作Deflation的严格定义

A_t 的谱分解为:


A_t = \sum_{i=1}^N \lambda_i u_i u_i^\top,

其中 $\lambda_1 \geq \lambda_2 \geq \dots \geq \lambda_N$,且 \{u_i\} 是标准正交基。

  • 第一次收缩
    定义剩余矩阵 $A_{t,2} = A_t - \lambda_1 u_1 u_1^\top$,其性质为:

    • 特征值:$\lambda_2, \lambda_3, \dots, \lambda_N$(即移除 \lambda_1 后剩余特征值不变)。
    • 特征向量:u_2, \dots, u_N 保持不变(因 u_1 与其他特征向量正交)。
  • k 次收缩
    递归定义:

    
    A_{t,k+1} = A_{t,k} - \lambda_k u_k u_k^\top,
    

    剩余矩阵 A_{t,k+1} 的特征值为 $\lambda_{k+1}, \dots, \lambda_N$。

每次收缩移除当前主成分,剩余矩阵的特征值是原始矩阵中未被移除的部分。


2. 剩余特征值的统计特性

目标:证明 \{\lambda_k(A_t)\}_{t \in \mathbb{Z}}k \geq 2 也是弱平稳的。

(1) 均值恒定性
  • 剩余矩阵的期望
    由线性性:

    
    E[A_{t,k+1}] = E[A_t] - \sum_{i=1}^k E[\lambda_i u_i u_i^\top].
    

    A_t 的元素分布时不变,且 \lambda_iu_i 的期望稳定(由主特征值的平稳性保证),则 E[A_{t,k+1}]t 无关。

  • 特征值期望
    对剩余矩阵 $A_{t,k+1}$,其主特征值 \lambda_{k+1}(A_t) 的期望近似为:

    
    E[\lambda_{k+1}(A_t)] \approx (N-k-1)\mu + v + \frac{\sigma^2}{\mu} \equiv \mu_{k+1},
    

    其中 (N-k-1)\mu 是剩余非对角元素的贡献(假设每次收缩后非对角元素统计特性不变)。

(2) 方差恒定性
  • 剩余矩阵的方差
    收缩操作通过正交投影移除 $\lambda_k u_k u_k^\top$,因此剩余矩阵 A_{t,k+1} 的元素方差仍为 $\sigma^2$(对角元素可能需调整)。
    由引理1的推广

    
    \text{Var}[\lambda_{k+1}(A_t)] \approx 2\sigma^2 \equiv \sigma_{k+1}^2.
    
  • 动态模型
    假设去中心化序列 \tilde{z}_{k+1,t} = \lambda_{k+1}(A_t) - \mu_{k+1} 服从AR(1)

    
    \tilde{z}_{k+1,t} = \rho_{k+1} \tilde{z}_{k+1,t-1} + \varepsilon_{k+1,t}, \quad |\rho_{k+1}| < 1,
    

    稳态方差为:

    
    \sigma_{\tilde{z}_{k+1}}^2 = \frac{\sigma_{\varepsilon_{k+1}}^2}{1-\rho_{k+1}^2} = \sigma_{k+1}^2.
    
(3) 协方差仅依赖滞后 $m$
  • 协方差函数:
    
    \gamma_{k+1}(m) = \text{Cov}(\tilde{z}_{k+1,t}, \tilde{z}_{k+1,t-m}) = \rho_{k+1}^{|m|} \sigma_{\tilde{z}_{k+1}}^2.
    
    仅依赖 $m$,与 t 无关。

3. 递推证明的完整性

  1. 归纳基础
    k=1 时(主特征值),平稳性已证。

  2. 归纳假设
    假设 \lambda_k(A_t) 的平稳性成立,即:

    • $E[\lambda_k(A_t)] = \mu_k$(常数),
    • $\text{Var}[\lambda_k(A_t)] = \sigma_k^2$(有限),
    • $\text{Cov}(\lambda_k(A_t), \lambda_k(A_{t-m})) = \gamma_k(m)$。
  3. 归纳步骤

    • 通过收缩操作,\lambda_{k+1}(A_t) 成为 A_{t,k+1} 的主特征值。
    • A_{t,k+1} 满足与 A_t 相同的统计假设(独立性、有界性、时不变性),则 \lambda_{k+1}(A_t) 的平稳性可类比主特征值的证明。

JB-test

JB-testJarque-Bera test 是一种用于检验样本数据是否服从正态分布的统计假设检验方法。这个检验特别适用于判断数据的偏度skewness和峰度kurtosis是否符合正态分布的特性。

正态分布具有以下特性:

  • 偏度Skewness 为 $0$,表示数据的分布是对称的。
  • 峰度Kurtosis 为 $3$,表示数据的峰度是"中等"的。

JB-test的统计量

Jarque-Bera统计量的计算公式为


JB = \frac{n}{6} \left( S^2 + \frac{(K - 3)^2}{4} \right)

其中:

  • n 是样本的大小。
  • S 是样本的偏度skewness衡量分布的对称性。
  • K 是样本的峰度kurtosis衡量分布的尖峭程度。

JB-test的分布和检验步骤

  • 零假设($H_0$:数据服从正态分布。
  • 备择假设($H_1$:数据不服从正态分布。

在进行检验时,首先计算 JB 统计量,然后将其与卡方分布进行比较:

  • JB 统计量的分布近似于自由度为 2 的卡方分布(当样本量较大时)。
  • 如果 JB 统计量的值大于临界值(根据设定的显著性水平,比如 $0.05$),则拒绝零假设,即认为数据不符合正态分布。
  • 如果 JB 统计量的值小于临界值,则无法拒绝零假设,即认为数据服从正态分布。

结论

  • 如果 JB 统计量接近 $0$,说明数据的偏度和峰度与正态分布的期望非常接近,数据可能符合正态分布。
  • 如果 JB 统计量远离 $0$,则说明数据的偏度或峰度与正态分布的特征差异较大,数据不符合正态分布。

特征值精度预估

1. 噪声随机变量与协方差

符号 含义
w_i i过程噪声样本
v_j j观测噪声样本
Q 过程噪声的真实方差(协方差矩阵退化)
R 观测噪声的真实方差(协方差矩阵退化)

说明

  • 在矩阵形式的 Kalman Filter 中,通常写作

    
      w_k\sim\mathcal N(0,Q),\quad v_k\sim\mathcal N(0,R).
    
  • 这里为做统计检验,把 w_i, v_j 当作样本,Q,R 就是它们在标量情况下的方差。


2. 样本统计量

符号 含义
N_w,\;N_v 过程噪声样本数和观测噪声样本数
\bar w 过程噪声样本均值
\bar v 观测噪声样本均值
s_w^2 过程噪声的样本方差估计
s_v^2 观测噪声的样本方差估计

定义


\bar w = \frac1{N_w}\sum_{i=1}^{N_w}w_i,\quad
s_w^2 = \frac1{N_w-1}\sum_{i=1}^{N_w}(w_i-\bar w)^2,

\bar v = \frac1{N_v}\sum_{j=1}^{N_v}v_j,\quad
s_v^2 = \frac1{N_v-1}\sum_{j=1}^{N_v}(v_j-\bar v)^2.

3. 方差比的 F 分布区间估计

  1. 构造 F 统计量

    
    F = \frac{(s_w^2/Q)}{(s_v^2/R)}
        = \frac{s_w^2}{s_v^2}\,\frac{R}{Q}
        \sim F(N_w-1,\,N_v-1).
    
  2. 置信区间(置信度 $1-\alpha$
    查得

    
    F_{L}=F_{\alpha/2}(N_w-1,N_v-1),\quad
      F_{U}=F_{1-\alpha/2}(N_w-1,N_v-1),
    

    
    \begin{align*}
    P\Big\{F_{\rm L}\le F\le F_{\rm U}\Big\}=1-\alpha \quad\Longrightarrow \quad P\Big\{F_{\rm L}\,\le\frac{s_w^2}{s_v^2}\,\frac{R}{Q}\le F_{\rm U}\,\Big\}=1-\alpha.
    \end{align*}
    
  3. 解出 \frac{R}{Q} 的区间

    
    P\Bigl\{\,F_{L}\,\frac{s_v^2}{s_w^2}\le \frac{R}{Q}\le F_{U}\,\frac{s_v^2}{s_w^2}\Bigr\}=1-\alpha.
    

    
      \theta_{\min}=\sqrt{\,F_{L}\,\frac{s_v^2}{s_w^2}\,},\quad
      \theta_{\max}=\sqrt{\,F_{U}\,\frac{s_v^2}{s_w^2}\,}.
    

4. 卡尔曼增益与误差上界

在标量情况下即状态和观测均为1维卡尔曼增益公式可简化为


K = \frac{P_k H^T}{HP_k H^T + R} = \frac{HP_k}{H^2 P_k + R}

针对我们研究对象,特征值滤波公式的系数都属于实数域。$P_{k-1}$是由上次迭代产生,因此可以$FP_{k-1}F^T$看作定值,则$P_k$的方差等于$Q$的方差,即:


\text{var}(P_k) = \text{var}(Q)

c = H, $m = 1/H$(满足 $cm = 1$),则:


K = \frac{cP_k}{c^2 P_k + R} = \frac{1}{c + m(R/P_k)} \quad R/P_k\in[\theta_{\min}^2,\theta_{\max}^2].

则极值为


K_{\max}=\frac{1}{c + m\,\theta_{\min}^2},\quad
K_{\min}=\frac{1}{c + m\,\theta_{\max}^2}.

通过历史数据计算预测误差的均值:


E(x_k' - x_k) \approx \frac{1}{M} \sum_{m=1}^{M} (x_k^{l(m)} - x_k^{(m)})\\

定义误差上界


\xi
=\bigl(K_{\max}-K_{\min}\bigr)\;E\bigl(x_k'-x_k\bigr)
=\Bigl(\tfrac1{c+m\,\theta_{\min}^2}-\tfrac1{c+m\,\theta_{\max}^2}\Bigr)
\,E(x_k'-x_k).

若令 $c,m=1$,可写成


\xi
=\frac{(\theta_{\max}-\theta_{\min})\,E(x_k'-x_k)}
      {(c^2+\theta_{\min})(c^2+\theta_{\max})}.

量化噪声方差估计的不确定性,进而评估卡尔曼滤波器增益的可能波动,并据此给出滤波误差的上界.

基于时空特征的节点位置预测

在本模型中,整个预测流程分为两大模块:

  • GCN 模块:主要用于从当前网络拓扑中提取每个节点的空间表示**。这里的输入主要包括:

    • 邻接矩阵 $A$:反映网络中节点之间的连通关系,维度为 $N \times N$,其中 N 表示节点数。(可通过第二章网络重构的方式获取)
    • 特征矩阵 $H^{(0)}$:一般是原始节点的属性信息,如历史位置数据,其维度为 $N \times d$,其中 d 是初始特征维度。
  • LSTM 模块:用于捕捉节点随时间变化的动态信息,对每个节点的历史运动轨迹进行序列建模,并预测未来时刻的坐标。
    其输入通常是经过 GCN 模块处理后,每个节点在一段时间内获得的时空融合特征序列,维度一般为 $N \times T \times d'$,其中 T 表示时间步数,d' 是经过 GCN 后的特征维度。

GCN 模块

输入

  • 邻接矩阵 $A$:维度 $N \times N$。在实际操作中,通常先加上自环形成

    
    \hat{A} = A + I.
    
  • 特征矩阵 $H^{(0)}$:维度 $N \times d$,每一行对应一个节点的初始特征(例如历史采样的位置信息或其他描述)。

图卷积操作

常用的图卷积计算公式为:


H^{(l+1)} = \sigma \Bigl(\tilde{D}^{-1/2}\tilde{A}\tilde{D}^{-1/2} H^{(l)} W^{(l)} \Bigr)

其中:

  • \tilde{A} = A + I 为加上自环后的邻接矩阵,
  • \tilde{D}\tilde{A} 的度矩阵,定义为 $\tilde{D}{ii} = \sum{j}\tilde{A}_{ij}$
  • H^{(l)} 表示第 l 层的节点特征,初始时 H^{(0)} 就是输入特征矩阵;
  • W^{(l)} 是第 l 层的权重矩阵,其维度通常为 $d_l \times d_{l+1}$(例如从 d 到 $d'$
  • \sigma(\cdot) 是非线性激活函数,例如 ReLU 或 tanh。

经过一层或多层图卷积后,可以得到最终的节点表示矩阵 $H^{(L)}$(或记为 $X$),维度为 $N \times d'$。
其中:

  • 每一行 x_i \in \mathbb{R}^{d'} 表示节点 i 的空间特征,这些特征综合反映了其在网络拓扑中的位置及与邻居的关系。

输出

  • GCN 输出:形状为 $N \times d'$;若将模型用于时序建模,则对于每个时间步,都可以得到这样一个节点特征表示。
  • 这里 d'>d 。1.高维嵌入不仅保留了绝对位置信息还包括了网络拓扑信息。2.兼容下游LSTM任务需求。

LSTM 模块

输入数据构造

在时序预测中,对于每个节点,我们通常有一段历史数据序列。假设我们采集了最近 T 个时刻的数据,然后采用“滑动窗口”的方式,预测 $T+1$、 T+2...

  • 对于每个时刻 $t$,节点 i 的空间特征 x_i^{(t)} \in \mathbb{R}^{d'} 由 GCN 得到;

  • 将这些特征按照时间顺序排列,得到一个序列:

    
    X_i = \bigl[ x_i^{(t-T+1)},\, x_i^{(t-T+2)},\, \dots,\, x_i^{(t)} \bigr] \quad \in \mathbb{R}^{T \times d'}.
    

对于整个网络来说,可以将数据看作一个三维张量,维度为 $(N, T, d')$。

LSTM 内部运作

LSTM 通过内部门控机制(遗忘门 $f_t$、输入门 i_t 和输出门 $o_t$)来更新其记忆状态 C_t 和隐藏状态 $h_t$。公式如下

  • 遗忘门

    
    f_t = \sigma(W_f [h_{t-1},\, x_t] + b_f)
    
  • 输入门和候选记忆

    
    i_t = \sigma(W_i [h_{t-1},\, x_t] + b_i) \quad,\quad \tilde{C}_t = \tanh(W_C [h_{t-1},\, x_t] + b_C)
    
  • 记忆更新

    
    C_t = f_t \odot C_{t-1} + i_t \odot \tilde{C}_t
    
  • 输出门和隐藏状态

    
    o_t = \sigma(W_o [h_{t-1},\, x_t] + b_o), \quad h_t = o_t \odot \tanh(C_t)
    

其中,x_t 在这里对应每个节点在时刻 t 的 GCN 输出特征;
[h_{t-1},\, x_t] 为连接后的向量;

LSTM 的隐藏状态 $h_i \in \mathbb{R}^{d'' \times 1}$(其中 d'' 为 LSTM 的隐藏单元数)捕捉了时间上的依赖信息。

输出与预测

最后,经过 LSTM 处理后,我们在最后一个时间步获得最终的隐藏状态 h_t 或使用整个序列的输出接着通过一个全连接层FC层将隐藏状态映射到最终的预测输出。

  • 全连接层转换公式
    
    \hat{y}_i = W_{\text{fc}} \cdot h_t + b_{\text{fc}}
    

其中,假设预测的是二维坐标(例如 xy 坐标),$W_{\text{fc}} \in \mathbb{R}^{2 \times d''}$,输出 \hat{y}_i \in \mathbb{R}^2 表示节点 i 在未来某个时刻(或下一时刻)的预测坐标。

image-20250411142712730

若整个网络有 N 个节点,则最终预测结果的输出维度为 $N \times 2$(或 $N \times T' \times 2$,如果预测多个未来时刻)。

疑问

该论文可能有点问题每个节点只能预测自身未来位置无法获取全局位置信息。如果先LSTM后GCN可能可以