哮喘是一种常见慢性异质性气道疾病，以咳嗽、咳痰、喘息等呼吸道症状为特征^[1]。气道炎症作为引发哮喘症状的关键环节，可导致气道高反应、黏液高分泌及重塑等过程^[2]。气道上皮作为抵御过敏原、病原体和其他有害物质进入机体的首道防线，在启动哮喘气道炎症过程中至为关键。

线粒体作为一种高度动态的双层膜细胞器，不仅负责细胞能量供应，还参与细胞基因组稳定、活性氧产生、离子稳态、应激反应及细胞命运决策等过程，通过不断融合分裂、生物发生及自噬等方式，实现自身质量控制。一般而言，细胞为了维持线粒体网络及内环境稳定，往往通过自噬方式选择性包裹和降解受损线粒体，完成线粒体内容物降解过程^[3]。研究表明，哮喘气道上皮细胞中线粒体功能存在障碍^[4-6]，上皮细胞线粒体自噬已引起学者们的广泛关注。已发现蟑螂过敏原处理的气道上皮线粒体自噬增加^[7]，脂溶性木脂素、芝麻素通过减少线粒体自噬和凋亡减轻哮喘气道炎症^[8]，环境污染物及过敏原诱发的线粒体功能障碍与自噬导致的生物能量学调节紊乱，可能是哮喘表型变化的重要机制^[9]。然而，线粒体自噬在哮喘中的调控机制及重要意义目前仍缺乏深入认识，该方向的深入挖掘可为哮喘治疗新靶点的研发提供有效依据。

线粒体自噬受一系列高度保守的调节因子调控，这些调节因子被称为线粒体自噬相关基因 （mitophagy⁃related gene，MRG）。本研究利用机器学习算法构建了基于 22 个 MRG 的哮喘预测模型，以受试者工作特征（receiver operating charac⁃ teristic，ROC）曲线、校准曲线、列线图、决策曲线分析和外部数据集验证等方法评估其效能；通过一致性聚类分析、根据 MRG 不同特征鉴别出两种哮喘亚型，比较每种亚型年龄、性别、生物学功能及特殊通路的差异；选用加权基因共表达网络分析法（weighted correlation network analysis，WGC⁃ NA）确定关键基因，利用药物基因组学数据库 （connectivity map，cMap）识别针对不同亚型可能的潜在小分子药物。

1 材料和方法

1.1 材料

1.1.1 数据的获得与处理

本研究选取基因表达综合数据库（gene expres⁃ sion omnibus data base，GEO）中GSE76226、GSE63142 及GSE67472 数据集的健康受试者（70例）与哮喘患者（289例）气道上皮刷检样本转录组数据，对其进行去批次处理、主成分分析评估去除情况；选取 GSE179156数据集验证预测模型性能。

线粒体自噬相关基因集“REACTOME_MI⁃ TOPHAGY.V7.5.1.Gmt”来自Reactome数据库。

1.1.2 动物

6~8周龄健康无特定病原体（specific pathogen free，SPF）级BALB/c雌鼠（约20 g）、用于原代气道上皮细胞提取所需的 SPF 级 1~3 d C57BL/6 乳鼠均购自南京医科大学医药实验动物中心，置于南京医科大学实验动物中心饲养，屏障中心适应1周后开始实验。本研究获南京医科大学动物福利伦理委员会批准（动物伦理编号：1906008）。

1.2 方法

1.2.1 基因集变异分析

采用基因集变异分析（gene set variation analy⁃ sis，GSVA）研究哮喘患者与健康受试者之间、不同哮喘亚型之间MRG整体表达差异。

1.2.2 差异表达MRG鉴定

运用“limma”软件包挖掘样本间差异表达基因 （differentially expressed gene，DEG），阈值为P <0.05； 取DEG与MRG交集作为哮喘中差异表达的MRG。

1.2.3 预测模型建立及验证

随机森林模型对 22 个 MRG 进行筛选，多因素 Logistic回归建立预测模型；构建ROC曲线，计算预测模型ROC曲线下面积（area under curve，AUC）；列线图模型，以校准曲线和决策曲线分析（decision curve analysis，DCA）评估该模型的预测能力；同时应用外部数据集进一步验证该模型的有效性。

1.2.4 共识聚类分析

根据 22 个 MRG 表达谱，对 289 例哮喘患者样本进行无监督聚类分析；根据一致聚类结果，确定最佳聚类数。

1.2.5 富集分析与加权基因共表达网络分析

采用基因本体论（gene ontology，GO）及京都基因与基因组百科全书（kyoto encyclopedia of genes and genomes，KEGG）通路分析法对不同哮喘亚型 DEG 进行富集，使用WGCNA识别每个亚型的核心基因。

1.2.6 潜在小分子药物鉴定

选取各亚型核心基因与DEG交集基因，将差异倍数最大的前 150 个上调与下调基因输入至 cMap 数据库中，药物评分范围为-100~100，评分越低表明该药物在相应哮喘亚型治疗方面的潜力越大。

1.2.7 急性哮喘小鼠模型构建

SPF级6~8周龄BALB/c小鼠适应性喂养后，随机分为哮喘组和对照组。哮喘组第0、7、14天腹腔注射2 mg/kg 鸡卵清白蛋白（ovalbumin，OVA）+80 mg/kg 氢氧化铝进行致敏，第20、21、22天鼻内滴注4 mg/kg OVA 进行激发；对照组 PBS 代替，处理方式同上。各组小鼠末次激发后次日取小鼠肺组织进行后续相关实验，所有小鼠麻醉后迅速进行操作。

1.2.8 免疫组化

SPF 级 6~8 周龄 BALB/c 小鼠肺组织石蜡标本切片，经脱蜡、水化、高压热修复后，3% H₂O₂灭活内源性过氧化物酶，室温孵育15 min，PBS冲洗3次 （5 min/次）；10%山羊血清封闭，室温孵育1 h，甩干，滴加线粒体外膜转位酶 5（translocase of outer mito⁃ chondrial membrane5，TOMM5）多克隆抗体（25607⁃1 ⁃AP，上海Proteintech公司），4℃过夜；室温下PBS冲洗3次（5 min/次），滴加生物素标记的山羊抗兔IgG，室温孵育 1 h；PBS 冲洗 3 次（5 min/次），显微镜下 DAB着色，蒸馏水充分洗涤；依次以苏木素复染40 s、蒸馏水冲洗、梯度乙醇脱水、二甲苯透明化、中性树胶固定封片；最后显微镜下观察，出现棕黄/褐色视为阳性染色。

1.2.9 气液交界（air⁃liquid interface，ALI）模型

C57BL/6乳鼠处死后酒精消毒，剪T字形切口、镊子取出肺组织，消化后提取原代肺球细胞。收集培养约2周的肺球，离心后加入消化酶解离，再次离心弃上清液、加入 PneumaCult^TM⁃Ex Plus Medium 培养基重悬，调整密度接种至 Transwell 小室上室，下室加入适量培养基，置于细胞培养箱（37℃、5% CO₂）培养；细胞长满上室后，仅于下室加入适量PneumaCult^TM⁃ ALI Medium 培养基、每 2 d 换液，第 21 天于下室予 10 ng/mL 白介素（interleukin，IL）⁃13 刺激，动态观察小鼠气道上皮形成及形态变化。

1.3 统计学方法

统计分析由 R 软件（版本 4.3.1）及 GraphPad Prism 9 软件处理，数据使用均数±标准差（$\stackrel{-}{x}\pm s$）表示，不同分组间的比较采用t检验。P <0.05为差异有统计学意义。

2 结果

2.1 哮喘气道中线粒体自噬失调

GSE76226、GSE63142及GSE67472气道上皮刷检样本数据集分析结果示，相比健康受试者，哮喘患者 MRG 整体表达显著上调（ 图1A、B）； GSE179156数据集进一步验证了该结果（图1C）；火山图显示，MRG 中共富集出 8 个 DEG，其中 1 个下调[自噬接头受体蛋白 1（sequestosome1，SQSTM1/ p62）]，7 个上调[TOMM5、FUN14 结构域蛋白 1 （FUN14 domain containing1，FUNDC1）、线粒体外膜转位酶 22（translocase of outer mitochondrial mem⁃ brane22，TOMM22）、泛素⁃核糖体结合蛋白 52（ubiq⁃ uitin A ⁃ 52 residue ribosomal protein fusion product， UBA52）、线粒体融合蛋白2（mitofusin⁃2，MFN2）、泛素B（ubiquitin B，UBB）、磷酸甘油酸变位酶5（phos⁃ phoglycerate mutase5，PGAM5），图1D]。

2.2 TOMM5免疫组化验证

通过分析比较哮喘患者与健康个体间的MRG，选取统计学差异最为显著的DEG即TOMM5进行后续验证。

小鼠肺组织免疫组化结果显示，气道上皮细胞表达 TOMM5，且 OVA 激发的哮喘小鼠气道上皮 TOMM5表达显著高于对照组（图2A），与分析结果一致。

ALI 模型观察小鼠气道上皮形成及形态变化，发现IL⁃13可显著上调小鼠原代气道上皮细胞中的 TOMM5蛋白水平（图2B）。

2.3 预测模型的建立和评估

基于基因表达数据与疾病状态或表型间的关联程度及基因间的相关性、冗余性、模型的预测性能等，运用多种机器学习算法（包括随机森林、逐步回归和多变量Logit模型）对22个MRG进行分析，根据多种因素综合考虑，筛选出与疾病最相关的7个特征基因 TOMM5、FUNDC1、TOMM22、SQSTM1、 PGAM5、MFN2、核糖体蛋白 S27a（ribosomal protein S27a，RPS27A）构建哮喘预测模型（图3A、B）；ROC 评估预测性能，结果显示AUC 值0.795，提示该模型对哮喘发病的预测性能良好（图3C）；随后基于上述 7个基因的表达模式，构建评分模型以计算哮喘确诊概率（图3D）；构建校正曲线和 DCA 曲线以评估和验证模型性能，结果显示，校正曲线中模拟曲线和实际曲线与理想线较为匹配（图3E），DCA曲线显示使用所有特征基因构建的模型净效益最高（图3F），体现出模型的可靠性。以上结果均提示，该模型在哮喘预测方面具有潜力，为进一步研究和应用提供了有力支持。

2.4 基于线粒体自噬特征的哮喘分型

基于 22 个 MRG 表达谱，对哮喘患者转录组学数据进行共识聚类分析，确定最佳聚类数 k=2（图4A），将哮喘患者分为 2 个亚型；同时采用热图 （图4B）形式比较两种亚型在线粒体自噬基因表达、临床特征（年龄和性别）以及诊断模型性能上的差异。

GSVA 分析表明两种亚型 MRG 整体表达上存在差异（图4C），亚型 1 中 RPS27A、UBB、MFN2、ATG5、ATG12、MAP1LC3B、CSNK2A1、TOMM20、 UBA52、TOMM22、TOMM5 及 FUNDC1 表达显著高于亚型2（图4D），应用预测模型对亚型1的预测准确性更高（图4E）；年龄方面两亚型间差异无统计学意义（图4F），但需指出的是亚型1中女性患者数量略多于男性患者（图4G）。

2.5 两种亚型GO及KEGG富集分析

对筛选出的DEG同时进行GO及KEGG富集分析。GO 分析结果示，两亚型间肽酶活性调控的功能特征差异显著，尽管存在一定共性，但在特定的生物过程中差异明显；亚型1显著富集于与体液免疫反应相关功能（图5A），亚型2则富集于与蛋白质三聚化、亚硝基化、金属离子转运以及对缺氧反应的相关功能（图5B）。

KEGG 分析结果示，两种亚型均与 IL ⁃17（由 Th17细胞和3型固有淋巴细胞产生的IL⁃17可刺激中性粒细胞性气道炎症）信号通路相关，但每种亚型均有各自独特的通路；亚型1特异富集病毒蛋白与细胞因子及受体相互作用、金黄色葡萄球菌感染、趋化因子信号等通路，亚型2特异富集与卡波西肉瘤相关的疱疹病毒感染、结核病及精氨酸生物合成相关通路（图5C）。

2.6 小分子药物预测

使用WGCNA识别与哮喘不同亚型相关的关键基因模块（图6A~C），以基于模块特征基因的连接性度量值作为筛选标准，确定与各亚型关系最密切的上调和下调基因；基于上述差异基因，cMap数据库中筛选出各亚型潜在小分子药物，结合cMap分析结果及相关文献，确定哮喘亚型1与哮喘亚型2最具希望的预测药物分别为 XMD8⁃92 和 Verrucarin⁃A（图6D，表1）。

3 讨论

哮喘由表观遗传调控和环境因素暴露相互作用驱动，影响着全球超过3亿的成人及儿童。作为一种常见呼吸系统异质性疾病，其表型多样、潜在发病机制并不明确^[10]。目前基于辅助性T细胞2 型 （T helper type2，Th2）状态，哮喘被分为Th2⁃高型及Th2⁃低型。尽管其分型不同，但临床症状如气短、喘息和咳嗽等却非常相似。因此，利用传统临床与生理参数，从生物学角度界定哮喘不同亚型以辅助临床诊疗，面临着重重困难^[11]。

线粒体功能障碍常见于各种急慢性肺部疾病^[12]。研究发现，哮喘患者气道上皮中可出现线粒体肿胀甚至碎裂，加剧氧化损伤、诱发上皮凋亡^[13]；抑制线粒体自噬可减少蟑螂过敏原诱导的气道上皮细胞因子释放^[7]；芝麻素可通过调节线粒体自噬、减轻哮喘气道炎症水平^[8]。系列研究结果均提示，线粒体自噬可能在哮喘的发生发展中扮演着重要角色。本研究通过系统分析健康受试者和哮喘患者气道上皮转录组数据，探究哮喘与线粒体自噬之间的关系，为哮喘的诊断和治疗提供理论依据。

研究结果显示哮喘患者MRG整体表达上调，其中7个上调显著、1个下调显著，证实哮喘患者气道上皮确实存在线粒体自噬状态改变。免疫组化进一步验证发现，OVA 诱导的哮喘小鼠气道上皮与 IL⁃13处理的原代气道上皮细胞的TOMM5表达均显著上调、与转录组学数据一致，表明该转录组数据生物信息分析可靠。

继而使用机器学习算法，从22个 MRG中筛选出 7个疾病最相关基因以构建哮喘预测模型，通过验证发现预测效能良好，同时为探寻疾病治疗新靶点提供了宝贵线索。其中，TOMM5与TOMM22是线粒体外膜转运酶复合物组分，可促进核编码蛋白质进入线粒体^[14-15]，在调节线粒体稳态中起重要作用。 FUNDC1 作为线粒体外膜蛋白，通过直接与 LC3 相互作用介导受体依赖性线粒体自噬，通过募集线粒体动力相关蛋白1（dynamin related protein 1，DRP1） 驱动线粒体分裂以应对缺氧应激^[16-17]。SQSTM1/ p62被广泛应用作为细胞自噬流量标志物^[18]，有证据显示慢性哮喘小鼠存在肺组织细胞自噬不足、 p62上调现象，敲降哮喘气道平滑肌中p62可显著抑制细胞增殖和迁移^[19]。PGAM5通常位于线粒体内膜、是一种丝氨酸/苏氨酸磷酸酶，应激条件下可稳定 PTEN诱导假定激酶1（PTEN⁃induced putative kinase1，PINK1）、募集 E3 泛素连接酶 PARKIN，降解受损线粒体^[20-21]，在肺纤维化中被认为是线粒体稳态的重要调节因子^[22]。MFN2具有诱导线粒体自噬等其他非融合性功能，而线粒体融合蛋白1（mitofusin⁃1， MFN1）的泛素化则促进线粒体降解、阻止受损线粒体融合^[23-24]。RPS27A构成了核糖体40S小亚基，在核糖体的生物发生中十分关键、参与调节p53的活性，进而可能影响p53依赖的细胞应激反应和细胞周期调控^[25]。

根据MRG表达谱进一步将哮喘分为两种亚型，分析了亚型间基因表达、临床特征和通路富集方面的差异。两亚型间MRG表达模式不同，亚型1的7个诊断基因中 5 个（RPS27A、MFN2、TOMM22、TOMM5、 FUNDC1）表达显著高于亚型2，提示两种亚型气道上皮的线粒体自噬状态存在显著差异。

在特定生物学过程和通路上，两种亚型亦不同。亚型1和病毒蛋白及细胞因子与细胞因子受体的相互作用、金黄色葡萄球菌感染以及趋化因子信号通路有关。作为哮喘发病机制中的重要因素，病毒及细菌感染被认为参与了哮喘的发生发展。人鼻病毒作为哮喘儿童常见的病毒类型，可导致儿童和成人哮喘的发展及恶化^[26]；感染鼻病毒的上皮细胞则对呼吸道细菌如金黄色葡萄球菌的易感性高，已证实金黄色葡萄球菌定植与哮喘患病率之间存在显著关联^[27-28]。趋化因子作为哮喘气道上皮细胞信号活化的重要因素，在哮喘发生与进展中作用关键。受损上皮通过释放 C⁃C 基序趋化因子 11 募集嗜酸性粒细胞，释放大量促炎介质促发并加重哮喘病程^[29]，其与配体间的相互作用在哮喘发病与治疗中不容忽视。因此，研究结果为微生物与哮喘间的联系及趋化因子在哮喘中的意义提供了进一步支持。

最后通过 WGCNA 分析以筛选小分子靶向药物，发现XMD8⁃92可能作为哮喘亚型1的靶向药物，是细胞外信号调节激酶 5（extracellular signal⁃regu⁃ lated kinase5，ERK5）的抑制剂。ERK5是参与细胞存活、增殖、分化、抗凋亡信号传导以及上皮⁃间充质转化等过程的重要分子^[30]，对哮喘气道重塑起重要作用^[31]，靶向抑制其活性或可使患者获益。Verrucarin⁃ A作为哮喘亚型2可能的靶向药物，是一种蛋白质合成抑制剂，可通过干扰肽基转移酶活性而抑制蛋白质的生物合成，或可具有防止过敏性鼻炎患者哮喘发展的潜能^[32]。

综上，本研究为MRG用于哮喘的预测、分型、个体化治疗提供了新的依据。通过分析哮喘中MRG 的表达模式变化，建立了具有跨不同数据集且高预测性能的哮喘预测模型；分别在哮喘动物及细胞模型层面证实 TOMM5 基因上调，揭示线粒体自噬在哮喘发生发展中的重要意义；通过分析MRG表达谱特征进行不同哮喘分型，揭示了哮喘的高度异质性。虽然本研究样本量相对较小，但所获结果仍具有重要价值，为深入理解哮喘的病理生理机制及研发新的治疗策略提供了新的线索。

参考文献