本研究提出了一种可重复的生物信息学分析流程，整合转录组学、机器学习和免疫浸润分析，以鉴定慢性阻塞性肺疾病（COPD）中潜在的诊断生物标志物和调控网络。
需要JoVE订阅才能观看此内容。 请登录或开始免费试用
研究文章
* These authors contributed equally
本研究提出了一种可重复的生物信息学分析流程，整合转录组学、机器学习和免疫浸润分析，以鉴定慢性阻塞性肺疾病（COPD）中潜在的诊断生物标志物和调控网络。
慢性阻塞性肺疾病（COPD）是一种以持续性气流受限和慢性炎症为特征的进行性呼吸系统疾病，但N4-乙酰胞苷（ac4C）RNA修饰在其发病机制中的作用仍基本未被探索。本研究旨在从已发表的数据库中系统筛选与ac4C相关的基因（ac4C-RGs），并分析其在COPD中的调控网络，从而识别潜在的生物标志物，用于后续机制研究，而不预设任何特定基因与ac4C修饰之间存在直接调控关系。通过转录组数据鉴定差异表达基因（DEGs），并采用加权基因共表达网络分析（WGCNA）挖掘关键共表达模块。对DEGs、显著模块和ac4C-RGs进行交叉分析。使用LASSO回归、XGBoost和随机森林算法筛选关键基因，随后基于逻辑回归构建诊断模型。通过受试者工作特征（ROC）曲线分析、曲线下面积（AUC）及其95%置信区间、校准曲线评估和决策曲线分析（DCA）评价模型性能。共鉴定出160个重叠基因，其中六个核心基因（PTRF、PRKCDBP、UPP1、TOR3A、FAM168B和B4GALT2）被三种机器学习算法一致选出。该诊断模型表现出良好的区分能力，在训练集、内部测试集和外部验证集中的AUC分别为0.766、0.759和0.723。调控网络分析提示可能存在竞争性内源RNA（ceRNA）轴及转录因子相互作用，而免疫浸润谱分析显示关键基因与多种免疫细胞亚群之间存在显著相关性。药物-基因相互作用分析和分子对接结果表明，氟尿嘧啶、卡培他滨和5-苄基环尿苷可能与UPP1具有良好的预测结合亲和力。综上所述，PTRF、PRKCDBP、UPP1、TOR3A、FAM168B和B4GALT2被鉴定为COPD中潜在的ac4C相关生物标志物，可能参与免疫和代谢调控，为未来的功能研究和治疗探索提供了基础。
慢性阻塞性肺疾病（COPD）是一种慢性且具有异质性的呼吸系统疾病，其特征是由于肺泡和气道结构异常导致进行性气流受限1,2。蛋白酶与抗蛋白酶失衡、氧化应激、慢性炎症以及细胞衰老构成了COPD的核心病理生理机制，导致肺组织的结构破坏和功能损害3,4。此外，COPD受多种危险因素影响，包括长期吸烟、环境污染、职业暴露、呼吸道感染以及遗传易感性5,6。COPD已对全球经济造成重大负担，预计在2020年至2050年间每年将占全球GDP的0.111%7。尽管目前的治疗策略（如支气管扩张剂、吸入性糖皮质激素、肺康复和长期氧疗）可缓解症状，但阻止疾病进展仍具挑战性。由于临床表现存在显著异质性，患者的预后差异较大8。因此，迫切需要新的诊断生物标志物和预后指标，以提升COPD的管理水平并改善患者生存率。
RNA修饰是指RNA分子发生的化学改变，可改变RNA的结构和功能，从而调控基因表达9,10。常见的RNA修饰包括N6-甲基腺苷（m6A）、假尿苷（Ψ）、5-甲基胞嘧啶（m5C）和ac4C11。ac4C修饰在维持mRNA稳定性及促进mRNA翻译方面发挥关键作用12,13。NAT10是目前唯一已知的可催化ac4C修饰的真核生物酶，其活性对ac4C修饰的形成至关重要14。研究表明，氧化应激、细胞衰老、炎症与ac4C修饰之间存在显著相关性。例如，炎症性肠病（IBD）患者结肠组织中的CD4+ T淋巴细胞表现出NAT10水平显著升高15。NAT10通过增强趋化因子CCL2和CXCL1的ac4C乙酰化，促进巨噬细胞和中性粒细胞的浸润，从而加剧炎症损伤16。细胞对氧化应激的响应可能涉及ac4C修饰，因为在氧化应激条件下ac4C水平显著上升。此外，NAT10通过ac4C修饰稳定TGFB1 mRNA，促进PM2.5诱导的肺纤维化，进而触发上皮-间质转化17。然而，ac4C修饰在慢性阻塞性肺疾病（COPD）中的作用仍 largely 未被探索，凸显了该领域进一步研究的必要性。
本研究利用 GEO 数据库鉴定慢性阻塞性肺疾病（COPD）患者与对照组之间的差异表达基因（DEGs）。随后，将来自多组学数据整合的已发表 ac4C 相关基因（ac4C-RGs）列表18 与差异表达基因进行整合，以识别重叠的候选基因。鉴于 ac4C 修饰已知可影响炎症和氧化应激——这两个过程在 COPD 中至关重要——我们假设与 ac4C 调控网络相关的基因可能在 COPD 中发生异常调控。然而，这些重叠基因并不被视为 NAT10 的直接底物或由 ac4C 直接调控的基因，而是被认为与 ac4C 相关网络有关的候选基因。通过多种机器学习算法筛选关键基因，随后构建并验证诊断模型。接着，确定了与 COPD 发病机制相关的关键通路，并预测了潜在的靶向药物。最后，通过 RT-qPCR 实验验证关键基因的表达水平，为 COPD 的病理生理机制提供了初步见解，并为治疗探索指出了潜在方向。
访问受限。请登录或开始试用以查看此内容。
机构审查委员会声明
本研究遵循《赫尔辛基宣言》进行。研究方案经深圳市罗湖区中医院伦理委员会批准（批件号：2024-LHQZYYYXLL-KY-039），所有参与者在入组前均签署了书面知情同意书。本方案所用研究工具与材料的详细信息见于 材料表.
数据来源与处理
与慢性阻塞性肺疾病（COPD）相关的基因表达数据集来自基因表达综合数据库（Gene Expression Omnibus, GEO）。GSE54837 数据集被用作转录组数据集，GSE112811 数据集则作为验证集（表 1）。ac4C 相关基因（ac4C-RGs）通过文献检索获得18。使用 R 软件包 limma 鉴定 COPD 组与对照组之间的差异表达基因（DEGs）。当 |log2FC| > 0 且 p < 0.05 时，差异表达基因被认为具有统计学显著性。通过绘制火山图以可视化基因表达变化的整体分布情况。
WGCNA 的构建
使用R对GSE54837数据集进行WGCNA分析，以鉴定与慢性阻塞性肺疾病（COPD）相关的模块。在网络构建之前，采用hclust函数结合平均连接法和欧氏距离度量进行层次聚类分析，识别并剔除离群样本。最优软阈值幂次（β = 10）被选定以实现无标度拓扑拟合指数 R2 ≥ 0.85，以平衡无标度拓扑结构与平均连接度。构建邻接矩阵并将其转换为拓扑重叠矩阵（TOM）。采用动态树切割算法（deepSplit = 2，minClusterSize = 50）识别基因模块。对具有特征基因相关性的模块 > 0.75 的模块随后使用 mergeCloseModules 函数进行合并。然后使用 Pearson 相关系数将模块特征基因与临床特征（COPD 状态、年龄、性别和吸烟状态）进行相关性分析，以鉴定与 COPD 相关的模块用于后续分析。
重叠基因的筛选、富集分析及蛋白质相互作用网络分析
使用 R 软件包 ggvenn 生成维恩图，以鉴定差异表达基因（DEGs）、MEsalmon 模块基因和 ac4C 相关基因（ac4C-RGs）之间的重叠基因。利用 R 软件包 clusterProfiler，基于基因本体（Gene Ontology, GO）和京都基因与基因组百科全书（Kyoto Encyclopedia of Genes and Genomes, KEGG）数据库对重叠基因进行功能富集分析。从 STRING 数据库（https://string-db.org/）获取蛋白质-蛋白质相互作用（PPI）信息，用于分析重叠基因在蛋白质水平上的相互作用关系。使用 Cytoscape 软件对所得的 PPI 网络进行可视化。
通过机器学习鉴定关键基因
应用了三种机器学习方法：最小绝对收缩与选择算子（LASSO）回归、极端梯度提升（XGBoost）以及随机森林（RF）。LASSO回归使用glmnet包实现，并采用10折交叉验证来确定最优惩罚参数λ。type.measure参数设置为“deviance”，family参数设置为“binomial”。通过λmin准则选择最优λ值，该准则可使交叉验证偏差最小化，最终筛选出17个基因。XGBoost使用xgboost包进行，设定以下超参数：nrounds = 100，max_depth = 6，eta = 0.3，subsample = 0.8，colsample_bytree = 0.8，eval_metric = “logloss”。特征重要性按增益（gain）指标排序，选取前30个基因。随机森林使用randomForest包实现，ntree = 200。特征重要性按Gini指数平均下降值排序，同样选取前30个基因。将三种机器学习方法筛选出的基因取交集，以确定后续分析中的关键基因。
构建并评估用于风险预测的逻辑回归模型
GSE54837 数据集被随机分为训练集（70%）和测试集（30%）。使用 MASS 软件包中的 glm 函数，在训练集上构建逻辑回归模型，以关键基因的表达水平作为输入特征。采用 pROC 软件包生成的 ROC 曲线评估模型性能。通过 2,000 次自助重复抽样计算 AUC 的 95% 置信区间。使用 rms 软件包生成的校准曲线（基于 1,000 次自助重采样）评估模型校准度。利用 dca 软件包进行决策曲线分析（DCA），以评估在不同阈值概率范围内的净临床获益。使用 rms 软件包中的 nomogram 函数构建列线图，以实现个体化风险预测。
回归方程为：
logit(P) = 0.5823 + 0.6010 × UPP1 - 0.6563 × PTRF + 0.3853 × B4GALT2 - 0.3972 × FAM168B + 0.1848 × PRKCDBP - 0.4787 × TOR3A. (1)
此处，P 表示 COPD 的预测概率，每个系数代表相应基因表达值对 COPD 对数优势的贡献。
表达分析、GeneMANIA 网络和分子调控网络
使用Wilcoxon秩和检验比较GSE54837数据集中慢性阻塞性肺疾病（COPD）组与对照组之间的基因表达水平。采用ggplot2软件包生成箱线图以可视化表达水平的分布情况，图中叠加显示中位数、四分位距（IQR）及各个数据点。利用GeneMANIA构建基因网络并预测功能相互作用。搜索时采用默认参数：物种 = Homo sapiens，最大相关基因数 = 20。将所得网络下载并进行可视化展示，边的颜色表示相互作用类型。构建竞争性内源RNA（ceRNA）网络以研究转录后调控机制。使用两个独立数据库预测靶向六个关键基因的miRNA：DIANA-microT（得分 ≥ 0.8）和miRanda（得分 ≥ 140，能量 ≤ −20 kcal/mol）。取两个数据库共同预测出的miRNA交集用于构建miRNA-mRNA配对关系。随后，利用StarBase数据库预测靶向这些miRNA的lncRNA。使用Cytoscape构建并可视化lncRNA-miRNA-mRNA调控网络。采用ChIP-X富集分析版本3（ChEA3）预测转录调控关系。对于每个预测出转录因子（TF）的关键基因，筛选出富集得分最高的前10个转录因子。在Cytoscape中构建TF-靶基因调控网络。
基因集富集分析与免疫细胞浸润评估
使用 clusterProfiler 软件包进行基因集富集分析（GSEA），以探究每个关键基因的生物学功能。对于每个关键基因，根据中位值将样本分为高表达组和低表达组。使用 limma 对两组间进行差异表达分析，并根据带符号的 log₂ 倍数变化对所得基因列表进行排序。采用 gseGO 函数对 GO 生物过程术语、gseKEGG 函数对 KEGG 通路进行 GSEA 分析，参数设置如下：minGSSize = 10，maxGSSize = 500，pvalueCutoff = 0.05，nPerm = 1,000。使用 GSVA 软件包中实现的单样本基因集富集分析（ssGSEA）估算了 28 种免疫细胞类型的相对丰度。包含 28 种免疫细胞类型标志基因的 curated 基因集特征矩阵来自既往文献19。对每个样本，应用 gsva 函数，参数设置为 method = "ssgsea"，ssgsea.norm = TRUE，kcdf = "Gaussian"。使用 cor.test 函数计算 ssGSEA 富集评分与六个关键基因表达水平之间的 Spearman 相关系数。p 值采用 Benjamini-Hochberg 方法进行多重检验校正。利用 pheatmap 软件包将相关性矩阵以热图形式可视化。
药物预测、分子对接和疾病关联分析
利用 DrugBank 数据库鉴定了靶向关键基因的潜在治疗化合物。在 Cytoscape 中构建了“关键基因-靶向药物”相互作用网络，以可视化预测的药物-基因相互作用。使用 CB-Dock2 平台进行分子对接，以评估结合亲和力。人源 UPP1 的三维蛋白结构从蛋白质数据库中获取（PDB ID: 7B8T）。药物分子结构（SMILES 格式）从 PubChem 获取。采用 AutoDock Vina 引擎进行对接，结果按结合自由能（ΔG，单位为 kcal/mol）排序。使用 PyMOL 可视化对接复合物。利用比较毒物基因组学数据库（CTD）分析关键基因与环境暴露相关的人类疾病之间的关联。对每个基因单独查询，并提取关联性最强的前十种疾病，使用雷达图进行可视化。
RT-qPCR 实验方案
从中国深圳罗湖区中医院的八名慢性阻塞性肺疾病（COPD）患者和八名健康对照者中采集外周静脉血样本。COPD 的诊断依据全球慢性阻塞性肺疾病倡议（GOLD）标准，定义为支气管扩张剂使用后 FEV1/FVC < 0.7。 < 0.70。对照组由年龄和性别相匹配、无呼吸系统疾病史且肺功能检查正常的健康志愿者组成（FEV1% 预计值 ≥ 80% 且 FEV1/FVC ≥ 0.70）。患者的基线信息如表所示 表 2从COPD血液样本中使用血液RNA提取试剂盒提取总RNA。对于cDNA合成，取500 ng总RNA，按照试剂盒提供的说明书，使用含基因组DNA去除步骤的cDNA合成试剂盒进行逆转录。将所得cDNA稀释至150 ng/μL。μL.
采用基于SYBR Green的qPCR预混液在实时PCR系统上进行RT-qPCR。每个10 μL反应体系包含5 μL 2x SYBR Green预混液、0.5 μL正向引物（10 μM）、0.5 μL反向引物（10 μM）、1 μL稀释后的cDNA（15 ng/μL）以及3 μL无核酸酶水。扩增程序为：95 °C预变性5 min，随后进行40个循环（95 °C变性10 s，60 °C退火延伸30 s），最后从60 °C升至95 °C进行熔解曲线分析，以验证扩增特异性。所有反应均设置技术重复三次。以β-actin作为内参基因。各目的基因的引物扩增效率通过标准曲线稀释系列验证，效率范围为90%至110%。基因表达水平经β-actin归一化后，采用2-ΔΔCt法计算相对表达量。COPD组与对照组之间的统计学差异分析采用Mann-Whitney U检验。
统计分析
使用 Cytoscape 创建网络可视化图，使用 R 软件进行统计分析。除非另有说明，对于非正态分布的数据采用 Mann-Whitney U 检验，对于正态分布的数据采用 Student t 检验来比较两组之间差异。以 p < 0.05 为具有统计学显著性。
访问受限。请登录或开始试用以查看此内容。
重叠基因的鉴定、富集分析及蛋白质相互作用网络构建
从 GSE54837 数据集中共鉴定出 3,371 个差异表达基因（DEGs），其中包括 1,675 个上调基因和 1,696 个下调基因。前 10 个显著上调和下调的基因如图所示 图1A对 GSE54837 数据进行了层次聚类分析（补充图 1A），并应用10的软阈值化幂以确保无标度网络拓扑结构（图1B）。采用动态树切割法构建基因共表达模块，最小模块大小设为50个基因，每个模块被赋予一种独特的颜色（补充图1B）。具有特征基因相关性的模块 > 0.75 随后被合并（补充图1C, 图1C），从而得到14个不同的模块。根据模块特征基因与临床特征之间的皮尔逊相关性分析，MEsalmon模块（包含5,226个基因）与慢性阻塞性...
访问受限。请登录或开始试用以查看此内容。
一种保守的RNA修饰ac4C主要存在于信使RNA（mRNA）和转运RNA（tRNA）中，可增强mRNA的稳定性及翻译效率20。NAT10是目前唯一已知介导ac4C修饰的RNA乙酰转移酶21。研究表明，NAT10在慢性阻塞性肺疾病（COPD）患者肺上皮细胞中表达上调。NAT10的敲低会破坏线粒体功能和转录组响应22。基于整合的多组学分析，本研究鉴定了六个与COPD密切相关的关键基因，并构建了一个诊断模型，该模型显示出中等诊断效能，具有进一步研究的潜在价值。进一步分析揭示了这些基因在转录调控、ceRNA网络以及免疫微环境中的关键作用。此外，还预测了潜在的靶向药物，为COPD的发病机制提供了新的见解，并有助于推动个体化精准治疗策略的发展。通过加权基因共表达网络分析（WGCNA）和差异表达分析，我们筛选出在COPD中差异表达且同时存在于既往发表的ac4C相关基因列表中的基因。需要强调的是，这些基因的筛选依据是其与ac4C调控网络的关联性，而非已证实与NAT10或ac4C乙酰化之间存在机制性联系，最终共获...
访问受限。请登录或开始试用以查看此内容。
作者声明无利益冲突。研究中所有受试者均签署了知情同意书。
我们感谢深圳市罗湖区中医院为本研究提供了必要的临床设施和行政支持。最后，我们衷心感谢所有参与本研究的患者和健康志愿者；他们的贡献对本工作至关重要。本研究得到了深圳市医疗卫生“三名”工程（编号：SZZYSM202401018）、罗湖区重点专科基金（编号：LX202402021）和罗湖区重点专科基金（编号：LX202302064）的资助。
访问受限。请登录或开始试用以查看此内容。
| 姓名 | 公司 | 目录编号 | 评论 |
|---|---|---|---|
| β-肌动蛋白引物 | Tsingke | N/A | 正向引物：5’-CATGTACGTTGCTATCCAGGC-3’ 反向引物：5’-CTCCTTAATGTCACGCACGAT-3’ |
| B4GALT2 引物 | Tsingke | N/A | 正向引物：5’-GGGCAGACTGCTGATCGAG-3’ 反向引物：5’-CCGGTGTCTAAAGGGGATGAT-3’ |
| CB-Dock2 | LabShare | Online | 分子对接 |
| clusterProfiler | Bioconductor | v4.14.6 | 富集分析 |
| Cytoscape | Cytoscape Consortium | v3.8.3 | 网络可视化 |
| DrugBank | University of Alberta | Online | 药物预测 |
| FAM168B 引物 | Tsingke | N/A | 正向引物：5’-TCTGGGGTTCCCTATGCAAAT-3’ 反向引物：5’-GTAGGATTCGCTCCAGGATACA-3’ |
| glmnet | CRAN | v4.1 | LASSO 回归 |
| GSVA | Bioconductor | v1.52.3 | ssGSEA 分析 |
| Hifair III 第一链 cDNA 合成超混合液 | YEASEN | 11141ES | cDNA 合成 |
| Hieff RTPCR SYBR Green 主混合液 | YEASEN | 11201ES | qPCR 扩增 |
| limma | Bioconductor | v3.54.0 | 差异表达分析 |
| LightCycler 480 II 系统 | Roche | LightCycler 480 II | 实时 PCR |
| PTRF 引物 | Tsingke | N/A | 正向引物：5’-GGGCCGTAGACCAGATCCA-3’ 反向引物：5’-CTTGCTCACCGTATTGCTCGT-3’ |
| PRKCDBP 引物 | Tsingke | N/A | 正向引物：5’-CACGTTCTGCTCTTCAAGGAG-3’ 反向引物：5’-TGTACCTTCTGCAATCCGGTG-3’ |
| R 软件 | R Foundation | v4.4.2 | 统计计算 |
| randomForest | CRAN | v4.7 | 随机森林 |
| RNA isolater MolPure 血液 RNA 提取试剂盒 | YEASEN | 19241ES50 | RNA 提取 |
| STRING 数据库 | EMBL | Online | 蛋白质相互作用网络 |
| TOR3A 引物 | Tsingke | N/A | 正向引物：5’-CCCTTGCTCTGTCGTTCCAC-3’ 反向引物：5’-CCCGTCCCGATACAGGTTC-3’ |
| UPP1 引物 | Tsingke | N/A | 正向引物：5’-CTGTCAGTCATGGTATGGGCA-3’ 反向引物：5’-GAGCACCGGGCATAGTACA-3’ |
| WGCNA | CRAN | v1.72 | 共表达网络 |
| xgboost | CRAN | v1.7 | XGBoost 算法 |
访问受限。请登录或开始试用以查看此内容。