帕金森病（Parkinson’s disease, PD）是常见的神经退行性疾病，主要表现为运动迟缓、震颤和肌肉僵硬等运动症状^[1]. 患者的临床表现具有异质性，主要在症状的严重程度、病程进展速度及对治疗的反应上^[2]. 根据临床表现的运动症状，帕金森病通常被分为3个亚型^[3]. TD亚型患者的病程进展较慢，震颤症状较突出，且对药物治疗的反应较好^[4]；PIGD亚型患者的病情恶化速度较快，常伴有姿势控制障碍和严重的非运动症状^[5]； IT亚型患者表现出混合型的特征^[6]. 现有的研究大多集中于亚型的临床表现和治疗反应，对于不同亚型背后的神经机制的研究尚未深入. 对帕金森病亚型神经机制的研究对于深入理解不同亚型的临床表现及其潜在机制具有重要意义^[7].

近年来研究表明，运动功能障碍并非帕金森病的唯一特征，深层次的脑功能异常可能在其中扮演重要角色^[8]. 脑电图（electroencephalogram, EEG）是非侵入性、时间分辨率高的技术，可以捕捉大脑的活动特征^[9]. 基于脑电构建脑网络，对不同脑区的功能连接进行分析，是脑相关疾病的常用手段^[10]. Aljalal等^[11]研究发现，帕金森病患者在低频段活动增强，在高频段活动低于健康个体. Yang等^[12]发现，PIGD患者在低频段连接增强，在额叶和顶叶区域表现显著. Orcioli-silva等^[13]发现，在特定任务下（如障碍规避）PIGD亚型的额叶皮层激活水平更高，这与其运动的自动调节机制有关. 通过EEG信号进一步分析功能连接网络^[14]，可以探究不同亚型的神经机制. 现有的EEG研究多关注帕金森患者与健康个体之间的差异，对静息态帕金森病亚型间脑功能差异的研究尚不充分^[15].

由于帕金森病亚型在症状和脑电特征上的显著异质性，传统的参数模型在分析这些信号时存在局限性，尤其是在应对帕金森患者脑网络的非平稳性和高变异性时^[16]. 非参数分析方法不依赖模型假设，能够更灵活地捕捉动态的神经信号特征，适合应用于帕金森脑电的复杂变化^[17]. 此外，本研究引入与多种神经功能密切相关的个体化α频段分析^[18-19]，探索帕金森病不同亚型的脑网络特征，以期揭示亚型在神经功能上的差异性，分析这些特征与临床表现的潜在联系.

本实验共招募38名帕金森患者，均来自天津市环湖医院. 经临床评估后，根据统一帕金森病评分量表（unified parkinson’s disease rating scale, UPDRS）中的震颤评分与强直评分的比值^[20]，对患者进行亚型分组，这一分类方法在国内外已被广泛应用^[5]. 比值大于1.5者纳入TD亚型，共10名患者，平均年龄为（64.2±4.54）岁；比值小于或等于1.2的人纳入PIGD亚型，共17名患者，平均年龄为（65.12±3.97）岁；比值为1.0~1.5的人^[21]纳入IT亚型，共11名患者，平均年龄为（63.08±3.99）岁. 本实验已由河北工业大学生物医学伦理委员会审核批准通过，所有受试者均签署了知情同意书. 实验受试者的年龄、性别、受教育程度经统计学检验，不具有统计学差异（p > 0.05）.

使用Neuracle公司NeuSen W系列32通道无线脑电采集系统采集脑电数据，所有受试者处于清醒和放松的状态，所有实验均在14:00—17:00于安静的脑电图室内进行. 根据国际标准导联10-20系统放置电极. 数据采样频率为1 000 Hz. 采集患者闭眼状态下5 min的脑电数据.

为了全面评估帕金森病患者的临床表现及其随疾病进展的变化，分析多种临床量表，以详细地监测患者的运动症状、非运动症状、认知功能及生活质量的多个方面. 帕金森评定量表使用UPDRS. 该量表是评估帕金森患者临床症状最全面的工具，包括以下4个部分：精神和行为的评估—第一部分（mentation, behavior, and mood, UPDRS-Ⅰ）、日常生活活动的影响—第二部分（activities of daily living, UPDRS-Ⅱ）、运动评估—第三部分（motor examination, UPDRS-Ⅲ）以及并发症—第四部分（complications of therapy, UPDRS-Ⅳ）. 认知功能评估：使用简易智力状态检查量表（mini-mental state examination, MMSE）和蒙特利尔认知评估量表（Montreal cognitive assessment, MoCA）2种工具，评估患者的整体认知状态. 这2种量表旨在评估患者的整体认知状态，以便排除认知功能对其他评估结果的潜在影响.

通过EEGLAB对脑电信号进行预处理，对原始脑电数据进行带通滤波（频段为0.5~50 Hz）去除工频干扰. 利用独立成分分析（independent component analysis，ICA）算法，去除信号中明显的伪迹信号，如心电、肌电、眼电成分. 降低采样率至256 Hz并进行重参考. 根据神经振荡活动的分布频段，将脑电信号按θ（4~7 Hz）、IAF频段、β（13~30 Hz）对脑电信号进行分频处理，提取脑电特征.

为了分析亚型间表现出的异质性问题，采用个体化的方法分析不同亚型的α频段：通过提取每个受试者的个体化峰值α频率和峰值α功率，实现个体化的频段分析. 功率谱密度^[22]（power spectral density, PSD）的计算基于快速傅里叶变换（fast Fourier transform, FFT）或多窗口谱分析，具体公式为

(1)$ {{\mathrm{PS D}}}(f)=\frac{1}{T} \int_{-T / 2}^{T / 2}\left|x(t) {{\mathrm{exp}}} \;{\left( -{\mathrm{i}} 2 {\text{π}} f t\right)}\right|^2 {\mathrm{d}} t. $

为了提取个体化α频段，选取枕区电极（O1、O2、Oz、PO3、PO4）的功率谱峰值位置，计算受试者的IAF^[23]. 峰值频率PAF的计算公式为

(2)$ {\mathrm{P A F}}={\displaystyle\sum_{i=1}^n f_i P_i}\left/{\displaystyle\sum_{i=1}^n P_i}\right. . $

式中：${f_i}$为频率，${P_i}$为该频率下的功率.

为了解决模型驱动脑网络研究的局限性问题，采用基于非参数脑网络的方法，构建静息态脑网络. 对信号进行连续小波变换，基于变换后的信号构建交叉谱矩阵，使用Wilson分解来估计传递函数和噪声协方差矩阵. 小波变换采用复Morlet小波变换估计信号的时频分布，定义如下：

(3)$ \psi (t) = \frac{1}{{\sqrt {{\text{π}} {f_{\mathrm{b}}}} }}{{\mathrm{exp}}\left(-{{{t^2}}}/{{{f_{\mathrm{b}}}}}\right)}{{\mathrm{exp}}\;{\left( {\mathrm{i}}2{\text{π}} {f_{\mathrm{c}}}t\right)}}. $

式中：${f_{\mathrm{b}}}$为带宽参数，${f_{\mathrm{c}}}$为中心频率.

对信号${x_n}(t)$的连续小波变换定义为

(4)$ {W_{{x_n}}}(t,s) = \frac{1}{{\sqrt s }} \int_{ - \infty }^\infty {x_n}(\eta ){\psi ^ * }\left(\frac{{\eta - t}}{s}\right){\mathrm{d}}\eta . $

对于每个频率尺度${{s}}$，构建交叉谱矩阵$ {{\mathrm{\bf{WS}}}}(f)$，其中${{\rm{WS}}}{_{ij}}(f)$为频率${{f}}$下通道${{i}}$和${{j}}$的小波变换系数的乘积：

(5)$ {\mathrm{{WS}}}{_{ij}}(f) = {W_{{x_i}}}(t,f)W_{{x_j}}^ * (t,f). $

谱矩阵分解选择Wilson分解，用于从交叉谱矩阵${\mathrm{\bf{WS}}}$中估计传递函数和噪声协方差矩阵. 基于得到的传递函数矩阵H和噪声协方差矩阵Σ，计算格兰杰因果关系，公式如下：

(6)$\begin{split} &F_{i \rightarrow j}(t, f)=\\&\ln \left(\frac{{\mathrm{W S}}_{j j}(t, f)}{{\mathrm{W S}}_{j j}(t, f)-\left({\boldsymbol{\varSigma}}_{i i}(t)-\dfrac{{\boldsymbol{\varSigma}}_{i j}(t)^2}{{\boldsymbol{\varSigma}}_{j j}(t)}\right)\left|H_{i j}(t, f)\right|^2}\right).\end{split}$

式中：${F_{i,j}}(f)$表示在频率${{f}}$下信号${{j}}$对信号${{i}}$的因果影响. 通过计算各通道间的非参数脑网络因果关系，得到脑电信号的功能连接网络.

当处理复杂脑网络时，选择适当的阈值至关重要，它决定了网络的稀疏性和保留的连接强度，这些会影响网络的功能特性. 选择效率成本优化^[24]的方法进行阈值选择. 效率成本优化是用于网络分析的评价准则，它能够有效地平衡网络的连通性和稀疏性. 阈值为0~1.0，步长为0.01，对于每一个阈值，各自计算脑网络的特征指标，并计算出阈值对应的质量指数. 质量指数越大，网络结构越优，选择合适的质量指数J，确定最优阈值.

质量指数的定义为

(7)$ J(\rho ) = \frac{{{E_{\mathrm{g}}}+{E_{\mathrm{l}}}}}{\rho }. $

式中：${E_{\mathrm{g}}}$为该网络的全局效率，描述整个网络中信息传递的效率；${E_{\mathrm{l}}}$为网络的局部效率，描述网络中单个节点及其邻居间信息传递的效率；$\rho $为网络的连接密度，等于网络中实际连接数与可能连接数的比例.

基于最优阈值筛选后的脑网络，可以利用图论方法计算关键网络的特征. 这些网络特征有助于揭示脑网络的拓扑属性和功能变化，从而深入理解帕金森病不同亚型的脑网络特性.

杠杆中心性^[25]是脑网络分析中的重要指标，用于识别网络中的关键节点. 该度量评估节点与其邻居节点之间的相对连接程度，从而衡量节点在网络中的影响力，评估网络的信息处理效率. 杠杆中心性的定义如下：

(8)$ L_{{i}}=\frac{1}{k_{i}} \displaystyle\sum_{j \in N(i)} \frac{k_{i}-k_{j}}{k_{i}+k_{j}} . $

式中：${k_i}$为节点${{i}}$的节点度数（即连接到节点${{i}}$的边的数量）；$N(i)$为节点${{i}}$的相邻节点的集合，通过比较节点${{i}}$与相邻节点的节点度差异来计算，从而衡量节点${{i}}$在网络中的影响力.

效率密度^[26]是衡量网络整体信息传递能力的一个指标，它通过结合全局效率和网络密度来描述网络信息的整合效率. 全局效率的计算公式为

(9)$ E_{{\mathrm{g }}}=\frac{1}{N(N-1)} \displaystyle\sum_{i \neq j} {l_{i j}}^{-1}. $

式中：${{N}}$为网络中的节点数，${l_{ij}}$为节点${{i}}$到节点${{j}}$的最短路径距离.

网络密度用于反映实际存在的连接边数与网络可能的最大边数之间的比例，计算公式为

(10)$ d(G) = \frac{{2L_{\mathrm{n}}}}{{N(N - 1)}}. $

式中：${{L_{\mathrm{n}}}}$为网络中的实际连接边数.

效率密度的计算公式为

(11)$ E(G) = d(G) \cdot {E_{{\mathrm{g}}}}. $

本研究中的数据统计分析主要是对脑网络特征值和临床量表评分（包括UPDRS、MMSE、MoCA等）进行计算，基于Shapiro-Wilk检验计算数据是否符合正态分布. 若结果符合正态分布，则使用单因素方差分析（analysis of variance, ANOVA）；若结果不符合正态分布，则基于Kruskal-Wallis秩和检验进行多组比较. 每个亚型的临床量表数据分别与对应亚型的不同脑网络特征数据进行Pearson相关性分析. 将显著性水平设定为0.05.

分别对不同亚型的年龄、受教育年限、MMSE和MoCA量表进行统计学检验，如表1所示. 其中，Q为问卷得分. 结果显示，3组间无显著性差异（见表1），表明3种亚型在整体认知功能上无显著区别，排除认知差异性带来的影响. UPDRS-II表现出显著性差异（见表1），这表明不同亚型在日常活动中的障碍程度更高.

由于帕金森病不同亚型在电生理特征上可能存在差异性，个体化研究常见于α频段峰值频率^[19]，分析各亚型在α频段的个体化峰值频率有助于揭示潜在的亚型间神经机制.

对PIGD、TD和IT亚型所有患者的α频段峰值频率进行统计分析，比较不同亚型下所有患者的平均值及其标准差，探讨不同亚型的神经活动模式.

如图1所示为3种帕金森病亚型在α频段的功率谱和个体化峰值频率. 其中，f为频率，P为功率谱密度. 具体来看，PIGD亚型的个体化峰值α频率平均为(9.54±0.55) Hz，TD组的平均频率为(8.99±0.61) Hz，IT亚型的平均频率为(9.64±0.77) Hz. 此外，PIGD与TD亚型的个体化α频率差异显著（p = 0.04），表明2种亚型在神经活动模式上存在一定的差异.

为了确保脑网络分析结果的稳定性和准确性，在不同阈值设置下评估网络特征的变化，选择最佳的阈值进行后续分析. 由于脑网络的稀疏性和连通性会受到阈值设定的影响，合适的阈值能够平衡网络连接的密度和信息传递的效率，从而更准确地反映各亚型的脑网络特征.

采用效率成本优化的方法，计算不同阈值下的质量指数. 基于式（7）计算不同阈值对应的质量指数，用以量化网络效率和连接密度的平衡状态. 通过绘制质量指数J与阈值T的关系曲线（见图2）可知，当质量指数为0.37时达到最大值.

分析2种关键的图论特征：效率密度E和杠杆中心性L（见图3、4）. 图4中，*表示p<0.05. 这些特征不仅提供网络拓扑属性的重要信息，而且反映脑网络在帕金森病不同亚型患者中的功能变化.

E是综合指标，反映了网络中节点的局部聚集能力与全局信息传递效率. E由代表相同节点和边数的随机网络的全局效率和网络密度计算.

如图3所示，在θ、IAF和β频段下分析PIGD、TD和IT 3种帕金森病亚型的E. 研究发现，在β频段中，PIGD亚型与TD亚型间的E表现出统计学上的显著差异（p = 0.04），说明在该频段下，不同亚型的脑网络功能的整体传输容量存在显著的差别. 相比之下，IT亚型在θ、IAF和β频段的E变化较小，显示其脑网络特征在这些频段内相对稳定.

杠杆中心性是通过评估不同脑区与其邻居的连接情况来衡量节点在网络中的相对重要性. 这一指标对于识别可能对网络行为产生不成比例影响的关键脑区至关重要. 在帕金森病患者脑网络中，具有异常L的区域可能在病理信号的传播中扮演关键角色，或者可能代表疾病进展中的代偿活动区域.

在本研究中，不同亚型之间在多个频段和脑区的L上显示出显著差异. 如图4（a）~（e）所示，在θ频段，PIGD亚型在额叶的L显著高于IT亚型（p < 0.05），而在顶叶区域，IT亚型显著低于PIGD亚型和TD亚型（p < 0.05）. 在IAF频段，TD亚型在中央叶的L显著高于PIGD亚型（p < 0.05），而在同一频段内，PIGD亚型和TD亚型在枕叶的L显著高于IT亚型（p < 0.05）.

分析脑网络特征（L和E）与帕金森病不同亚型患者的临床量表之间的相关性.

结果显示，PIGD亚型在图5（a）即IAF频段下的E与UPDRS-Ⅰ评分呈显著负相关（r = −0.47，p = 0.04）. 此外，在θ频段下，如图5（b）所示，PIGD亚型的颞叶L与UPDRS-Ⅱ评分呈显著正相关（r = 0.61，p = 0.008）. 如图5（c）所示，在枕叶区域，PIGD亚型在θ频段下的L与UPDRS-Ⅱ评分呈负相关（r= −0.57，p = 0.01）. 如图5（d）所示，对于TD亚型，在β频段下，顶叶的L与UPDRS-Ⅱ评分呈显著正相关（r = 0.66，p = 0.03）.

本研究深入探讨了帕金森病不同亚型在脑电活动和脑网络特征方面的差异，以及这些特征与临床症状之间的相关性. 通过非参数脑网络分析，研究发现不同帕金森病亚型在多频段下的脑电活动存在一定的差异，特别是在L和E特征上，这些特征与临床症状表现出相关性.

分析E可知，不同的帕金森病亚型在不同频段下的脑网络活动有一定的差异. PIGD亚型的E在β频段下显著低于TD亚型，显示出PIGD亚型在高频段的脑网络信息整合传输功能受到更大的影响. PIGD亚型临床表现更加严重^[5]，这可能与其在高频段的脑电整体传输能力低于其他亚型有关，该结果反映出其在姿势控制相关神经网络中的功能受损^[27]. 与此相对，TD亚型在β频段下的E保持较高的水平，提示其脑电活动在高频段较活跃，尤其是与震颤症状相关. 这一发现与震颤相关临床量表的正相关性结果吻合，且此表现与Orcioli-silva等^[12]发现TD亚型患者在β频段活动的显著变化结果一致，表明TD亚型β频段脑电活动的增强可能与震颤控制机制的受损密切相关.

不同的帕金森病亚型在各脑区的L分析显示出明显差异，在脑网络功能连通性上表现出不同，如图4（a）~（e）所示. 在θ频段下，PIGD亚型在额叶区域的L显著高于IT亚型，这表明额叶在PIGD亚型中可能发挥更重要的作用，PIGD亚型临床表现出的姿势不稳与步态困难与额叶的神经调控密切相关. 此外，较高的L提示该区域在运动功能代偿中扮演较核心的角色. 相比之下，TD亚型在IAF频段中央叶的L显著高于PIGD亚型，表明TD亚型患者的中央运动皮层处理震颤信号时的网络连接活跃性更强. IT亚型的分析结果显示，L在不同频段和脑网络特征中的变化较小，这一现象可能表明IT亚型患者的脑网络在各频段下保持一定的一致性.

至于个体化脑电部分，研究初期针对帕金森病患者的黄金研究频段-β频段进行亚型间的频段分析，个体化β频段结果与常规β及常规α在不同亚型的峰值频率分析表现均无明显差异（p > 0.05）. 不同亚型在IAF频段下中央叶和枕叶的L存在显著差异，结合IAF频段峰值频率的结果：TD亚型在该频段下的表现显著低于其他2种亚型，不同亚型在IAF频段的多个特征表现出显著差异，这些发现进一步表明，个体化脑电特性在帕金森亚型研究中具有潜在的应用价值. 如图4（e）所示的枕叶部分是脑电图中α频段的主要出现区域. 近年来，针对IAF和枕叶的研究逐渐增多，IAF作为内源性频率，与认知功能、注意力调节密切相关，已在一些电生理刺激中作为刺激频率选择的关键参数展现出应用价值^[28]. 枕叶是偏头痛这类脑部病变的治疗中可用的电生理刺激靶点^[29]，本研究中发现帕金森病亚型中IAF的差异性，提示枕区的IAF研究可能成为未来亚型精确治疗的重要切入点，优化刺激参数选择，改善帕金森病患者的症状.

通过进一步的相关性分析，研究发现脑网络特征与多种临床量表之间存在密切的联系，提示特定脑区的网络功能变化可能对亚型化的症状表现及疾病进展起到重要作用. 图5（a）中， PIGD亚型下IAF频段的E与UPDRS-I评分呈显著负相关，表明脑网络整体传输效率的变化与生活行为有关. 在图5（b）的θ频段，PIGD亚型的颞叶L与日常生活运动评分呈显著正相关，指示颞叶网络节点的增强与日常功能改善相关，颞叶本身与语言理解、听觉处理以及情绪调节等日常生活功能有关^[30]，这一发现为通过认知训练或非药物治疗手段优化PIGD亚型的颞叶功能提供了可能方向. 如图5（d）所示，对于TD亚型，β频段下顶叶的L与运动症状评分显著正相关，这可能反映了顶叶在震颤控制中的关键作用. Vaccari等^[31]的研究发现，顶叶与运动控制有关，顶叶受损主要表现出感觉障碍，这与本研究发现帕金森亚型间枕叶部分脑网络特征的显著性差异变化结果具有一致性，验证了TD亚型运动震颤症状的严重程度与顶叶功能的变化密切相关.

通过分析帕金森病患者的脑网络特征与临床表现之间的关联，揭示了不同亚型在脑电活动和脑网络特征上的差异性，为理解帕金森病的神经生物学基础提供了新的视角. 研究发现，不同亚型在额叶、顶叶和枕叶区域的脑网络特征变化显著，提示这些关键脑区的功能变化可能在亚型特异性症状表现中起到了重要作用. 这些结果为未来针对帕金森病亚型的辅助康复治疗提供了基础支持.

未来的工作应着重于进一步验证这些脑网络特征在帕金森病亚型中的稳定性和可靠性，结合亚型随病程进展的动态变化进行深入研究，帮助优化 不同亚型患者的康复干预策略，推动帕金森病亚型的研究变得更加精细化，并向临床应用转化的道路迈进. 同时，进一步扩大样本量，引入其他类型的数据（如核磁共振数据），以丰富数据的维度. 这不仅可以从多个角度探讨帕金森病亚型间的差异，为理解不同亚型的神经机制提供新的线索，而且可以提高统计分析的效能，增强结果的稳定性与可重复性，为更全面地覆盖帕金森病不同亚型患者的临床异质性提供可能.

本研究针对帕金森病亚型的高度异质性问题，使用非参数脑网络分析方法，应对传统分析方法在处理脑电信号非平稳性时的局限性. 通过这种方法，研究能够更全面地分析不同亚型患者在多个频段下的脑电活动和脑网络功能特征，特别是在L和E这类图论指标上，展现出较好的区分能力. 这一分析不仅揭示了脑网络特征与运动症状之间的复杂关系，而且结合个体化频段分析，捕捉到了不同频段下更具个体差异的神经活动模式，为理解亚型间的神经机制提供了更深刻的见解. 研究结果显示，脑网络特征在区分不同亚型的帕金森病患者时起到了重要作用，为未来个体化治疗和辅助康复策略的制定提供了基础参考. 本研究的现有样本已满足分析要求，但适当增加并引入多模态数据可以提高结果的稳健性，为更深入的帕金森病亚型研究提供更强的支持.

综上所述，本研究表明非参数脑网络分析在帕金森病亚型研究中具有显著的优势，结合个体化频段方法，针对帕金森病亚型间的异质性问题，提供了更灵活和精准的分析手段. 利用该方法，能够捕捉到不同亚型在脑电活动和脑网络功能上的显著特征差异，为理解帕金森病亚型的运行机制提供了新的视角. 此外，本研究的结果为后续客观识别帕金森病亚型提供了一定的支持，也为其他神经精神类疾病的亚型研究提供了借鉴和参考，展示了该方法广泛应用的潜力.