0 引言

       成人型弥漫性胶质瘤作为中枢神经系统中最常见的原发性恶性肿瘤，是导致原发性脑肿瘤死亡的主要原因^[1]。由于其具有复杂的时空异质性，尽管当前存在多种治疗手段，包括手术切除、放疗、化疗以及前沿的溶瘤病毒免疫诱导疗法等^{[2, 3]}，患者的总体预后情况仍不容乐观。2021年世界卫生组织（World Health Organization, WHO）中枢神经系统肿瘤分类第五版中，依旧强调了异柠檬酸脱氢酶（isocitrate dehydrogenase, IDH）在胶质瘤分类与预后评估中的核心地位。具体而言，相较于IDH野生型胶质瘤，IDH突变型胶质瘤的患者往往展现出更为积极的预后趋势^[4]。相较于“金标准”病理学，MRI在量化成人型弥漫性胶质瘤全局异质性中具有独特优势，以及无创下预测IDH分子状态对于指导治疗并改善患者预后具有重要意义^[5]。

       影像组学在肿瘤领域中的应用取得了快速进展^[6]，即可以从影像图像中提取定量的手工特征来揭示肿瘤的内部信息^{[7, 8]}。研究表明，多模态MRI影像组学在预测胶质瘤病理学特征上的贡献，已得到了广泛认可^[9]。然而，既往的研究大多侧重于病变区的全局特征提取^[10]，即强调整个病变组织的关联性而忽略了“真实情况”下肿瘤微环境异质性的空间局限性^[11]。因此，可以考虑通过将复数的局部特征汇总以更好地表达全局异质性。熵是研究最多的手工特征之一，其主要表征体素信号分布的变化。通过熵值的计算，可以对特定区域的异质性进行量化分析，从而为临床诊断提供重要的参考依据^[12]。已有研究证实，熵在食管癌^[13]、乳腺癌^[14]和非小细胞肺癌^[15]中表现出优异的诊断性能，但在胶质瘤领域的研究尚显不足。

       在大多数研究使用全局特征表征肿瘤异质性的基础上，使用整体肿瘤作为感兴趣区将难以避免地进一步降低研究结果的可靠程度。近年来，生境分析的概念被引入医学领域，它对肿瘤及其毗邻区域内的每一个体素进行细致的异质性描述，并通过聚类方法汇集相似体素以解释生物多样性问题，为描述肿瘤异质性提供了一种有效手段^[16]。目前，生境成像分析已应用于成人型弥漫性胶质瘤的单中心小样本基因亚型预测及生存预后方面^{[17, 18]}，一个多中心来源的样本将有利于构建模型的泛化。

       因此，在本项回顾性多中心研究中，我们基于多模态磁共振成像的栖息地局部熵值，旨在表征成人型弥漫性胶质瘤的异质性，然后开发和验证一个用于预测IDH分子状态的综合模型。

1 材料与方法

1.1 一般资料

       在这项研究中，我们回顾性收集和分析了来自内蒙古医科大学附属医院、加州大学弗朗西斯科分校（University of Pennsylvania, UPENN）以及癌症基因图谱（The cancer genome atlas, TCGA）的数据，共计533例受试者。其中，包括了在2017年7月至2024年2月于内蒙古医科大学附属医院就诊的60例受试者，本研究遵守《赫尔辛基宣言》，经内蒙古医科大学附属医院伦理委员会批准，免除受试者知情同意，批准文号：KY2023064。UCSF和TCGA 数据为公共数据库资源，391例的UCSF数据收集日期为2015年至2021年，82例的TCGA 数据收集日期为2008年至2013年；二者的完整数据已分别获得UCSF和TCGA的伦理委员会批准且均可在癌症影像档案（The Cancer Imaging Archive, TCIA）获得，免除受试者知情同意。纳入标准：（1）遵照WHO CNS5经病理诊断为成人型弥漫性胶质瘤；（2）有完善的术前MRI资料且图像质量良好；（3）在扫描前未接受过治疗（包括手术、放疗、化疗和联合放化疗）。排除标准：（1）扫描时间与手术时间间隔大于30天；（2）复发患者。

1.2 MRI扫描

       内蒙古医科大学附属医院受试者均于术前在3.0 T磁共振扫描仪（MAGNETOM Skyra; Siemens Healthcare, Erlangen, Germany）上采用32通道头颈联合线圈进行扫描。扫描序列包括轴位反转恢复T1WI、轴位快速自旋回波T2WI、轴位扩散加权成像（diffusion-weighted imaging, DWI）、T2液体衰减反转恢复（fluid attenuated inversion recovery, FLAIR）和对比增强T1加权成像（contrast enhanced T1-weighted imaging, T1WI-CE）[于静脉注射钆喷酸葡胺注射液（Gadovist, Bayer Healthcare）0.1 mmol/kg后进行]。相关采集协议见表1。UCSF及TCGA公共数据库相应的医学成像数据在TCIA（https://www.cancerimagingarchive.net/）中下载。

点击查看表格

表1  MRI扫描序列参数
Tab. 1  Different sequence scanning parameters of MRI

1.3 图像预处理和全病变区分割

       首先，使用SimpleITK对除扩散以外的所有图像进行N4偏置场校正，以减少信号强度不均匀性。然后，基于生成对抗网络预训练的超分辨率重建模型^[19]提升图像分辨率。适当倍率的超分辨率重建被认为有利于后续诊断模型的构建，因此在这项研究中我们将原始分辨率提升了2倍。随后，我们使用RIAS医学影像预处理工具（v0.2.2, https://github.com/lisherlock/RMIT）对图像依次进行了各向同性重采样（1×1×1）和z-score标准化。我们在MATLAB平台中使用SPM12的默认参数（B-spline）进行不同模态之间的配准。

       在图像分割之前，首先使用公开的深度学习算法进行脑提取（https://www.github.com/ecalabr/brain_mask）。然后，将之前预处理后的图像使用改进的nnU-Net深度学习分割框架^[20]进行全病变区域（肿瘤实体加瘤周水肿）的自动分割。最后，由一名低年资（具有2年工作经验）神经放射科医师使用3D-Slicer（https://www.slicer.org）软件进行了手动校正，最终由另一名高年资医师（具有30年工作经验）许可。

1.4 图像离散化和熵值转换

       我们首先对训练集进行了全局图像离散化^[21]，即根据整个队列的实际最大最小值进行离散化直方图，然后通过等间距Bin宽将Bin数调整为32。考虑到全病变区域的大小，使用3×3×3的矩阵进行全病变区域局部熵值的获取。熵值在影像组学中是一个衡量信息混乱程度的量化指标，它反映了图像中纹理的复杂程度或像素分布的模式。具体来说，熵值越大，表示图像纹理的不均匀程度和复杂度越高。可视化熵图见图1、2。最后，通过训练集衍射离散化直方图信息到测试集中，以避免信息泄漏的可能。

点击查看大图

图1  不同MRI序列经熵值转化后的结果可视化图。图为一例IDH突变型少突胶质细胞瘤的受试者，1A~1F为不同MRI序列图像（从左至右依次为T1、T2、T1-CE、T2-FLAIR、DWI、ADC），1G~1L为不同MRI序列图像所对应的熵图（例如1A对应的熵图为1G，1B对应的熵图为1H）。CE：对比增强；FLAIR：液体衰减反转恢复；DWI：扩散加权成像；ADC：表观扩散系数。
Fig. 1  The results visualization after entropy transformation of different MRI sequences. There is an IDH-mutated oligodendroglioma patient in the picture. 1A-1F represent images from different MRI sequences (from left to right, they are T1, T2, T1-CE, T2-FALIR, DWI and ADC), and 1G-1L correspond to the entropy maps of those MRI sequences (for example, the entropy map corresponding to 1A is 1G, and the entropy map corresponding to 1B is 1H). CE: contrast enhancement; FLAIR: fluid attenuation inversion recovery; DWI: represents diffusion-weighted imaging; ADC: apparent diffusion coefficient.

点击查看大图

图2  不同MRI序列经熵值转化后的结果可视化图。图为一例IDH野生型胶质母细胞瘤的受试者，2A~2F为不同MRI序列图像（从左至右依次为T1、T2、T1-CE、T2-FLAIR、DWI、ADC），2G~2L为不同MRI序列图像所对应的熵图（例如2A对应的熵图为2G，2B对应的熵图为2H）。CE：对比增强；FLAIR：液体衰减反转恢复；DWI：扩散加权成像；ADC：表观扩散系数。
Fig. 2  The results visualization after entropy transformation of different MRI sequences. There is an IDH wild-type glioblastoma patient in the picture. 2A-2F represent images from different MRI sequences (from left to right, they are T1, T2, T1-CE, T2-FALIR, DWI and ADC), and 2G-2L correspond to the entropy maps of those MRI sequences (for example, the entropy map corresponding to 2A is 2G, and the entropy map corresponding to 2B is 2H). CE: contrast enhancement; FLAIR: fluid attenuation inversion recovery; DWI: diffusion-weighted imaging; ADC: apparent diffusion coefficient.

1.5 生境异质性划分和特征提取

       采用K-means无监督聚类对T2-FLAIR和T1-CE图像进行全病变区的二维分割，聚类中心的取值范围为2~5。然后，使用Calinski-Harabasz、Silhouette Coefficient和Davies-Bouldin Index对分割好的结果进行评价及筛选。对于特征提取，在将分割好的生境亚区域映射到其他模态上后，获取全部模态不同生境亚区的16个一阶特征，包括10分位数、90分位数、四分位距、范围、最大值、最小值、均值、中位数、偏度、峰度、方差、均方根误差、平均绝对偏差、鲁棒平均绝对偏差、能量和总能量，这些特征的获取并不需要再次的图像离散化。此外，还获取了全病变区的以上全局特征，以用于机器学习建模后的比较。

1.6 模型构建

       FeAture Explorer（FAE, 0.6.0,https://github.com/salan668/FAE）被用于模型建立。在预测模型的开发过程中考虑了多条流水线组合，包括1种数据平衡方法（SMOTE）、3种特征归一化方法（MinMax、Z-Score和Mean）、2种降维方法（PCA和PCC）、4种特征选择方法（ANOVA、KW、RFE和Relief）和10种机器模型（AB、AE、DT、GP、LDA、LR、LRLASSO、NB和SVM），共计 240 条基本流水线。本研究以UCSF公共数据库为训练集，并通过10折交叉验证进行内部验证，其余2个数据库均作为独立测试集（本院数据为测试集1，TCGA数据库为测试集2）。模型特征数量的选择和超参数的调整通过交叉验证完成。

1.7 统计学方法

       计量数据以均数±标准差表示。组间年龄的比较采用单因素ANOVA，分类变量比较采用χ²检验、Fisher's检验或Mann-Whitney U检验。所有统计分析都是双侧的，P<0.05表示差异具有统计学意义。以受试者工作特征（receiver operating characteristic, ROC）曲线及其曲线下面积（area under the curve, AUC）、准确度、敏感度和特异度等指标评价模型诊断效能，然后使用DeLong检验比较模型间的差异性。模型与实际结果之间的偏差通过校准曲线可视化。采用决策曲线分析最终预测模型在不同阈值概率下的净效益，以增加在临床实践中实际应用的可能性。

2 结果

2.1 临床资料

       在不同数据集之间年龄差异具有统计学意义（P=0.010），性别差异无统计学意义（P=0.209）。各数据集的临床资料统计描述详见表2。

点击查看表格

表2  不同数据集临床资料的比较
Tab. 2  Comparision of clinical data in different datasets

2.2 生境异质性划分结果评估

       聚类中心为2或3时，对应的Calinski-Harabasz分别为95 080和100 379，Silhouette系数分别为0.477和0.422，Davies-Bouldin指数分别为0.741和0.810。聚类中心为4或5时，对应的Calinski-Harabasz分别为87 522和76 425，Silhouette系数分别为0.364和0.326，Davies-Bouldin指数分别为0.917和0.983。Calinski-Harabasz值越大、Silhouette系数越接近1以及Davies-Bouldin指数越接近0往往意味着聚类的结果更理想。在权衡4种不同聚类结果的3个评价指标后，因聚类2和3的结果更相近且优，聚类4和5的结果较差，所以在后续研究中将仅分析聚类2与3。

2.3 影像组学结果

2.3.1 特征筛选结果

       本研究采用多模态序列构建影像组学模型，即整合了6个序列的影像学特征。针对不同栖息地区域（全病变区、聚类2和聚类3），我们分别构建了最佳流水线（表3）。全病变区的多模态影像组学模型由6个特征构成，其中系数的绝对值最大的特征是T1-CE的90分位数（系数为2.459），即该特征对模型贡献度最大，其他特征的系数范围为-0.594~0.275。聚类2区域的多模态影像组学模型中由7个特征构成，其中系数的绝对值最大的特征是聚类区域1的T1-CE的均值（系数为21.114），其他特征的系数范围为-8.953~4.635。聚类3区域的多模态影像组学模型由8个特征构成，其中系数的绝对值最大的特征是聚类区域1的T1-CE的均值（系数为8.592），其他特征的系数范围为-5.534~7.492。特征的具体分布见图3。

点击查看大图

图3  特征相关性热图。IDH为异柠檬酸脱氢酶；rMAD：稳健的平均绝对偏差；FLAIR：液体衰减反转恢复；CE：对比增强；DWI：扩散加权成像；ADC：表观扩散系数。
Fig. 3  Heatmaps of feature correlation. IDH: isocitrate dehydrogenase; rMAD: robust average absolute deviation; FLAIR: fluid attenuation inversion recovery; CE: contrast enhancement; DWI: diffusion-weighted imaging; ADC: apparent diffusion coefficient.

点击查看表格

表3  用于多模态模型构建的特征选择
Tab. 3  Selected features for model construction

2.3.2 模型诊断性能评估及临床可行性

       成人型弥漫性胶质瘤不同栖息地区域的多模态影像组学模型在预测IDH分子状态方面均表现出优异的诊断性能（图4）。在全病变区的多模态影像组学模型中，训练集AUC为0.942（95% CI：0.919~0.966），验证集为0.933，测试集1为0.861，测试集2为0.847。在聚类2区域的多模态影像组学模型中，训练集AUC为0.974（95% CI：0.959~0.989），验证集为0.963，测试集1为0.739，测试集2为0.754。在聚类3区域的多模态影像组学模型中，训练集AUC为0.963（95% CI：0.945~0.980），验证集为0.950，测试集1为0.864，测试集2为0.769。当聚类数为2时，两个独立测试集的敏感度较高（分别为95.2%和80.0%）；当聚类数为3时，两个独立测试集的特异度较高（分别为72.2%和89.2%）（表4）。根据DeLong检验，聚类2的测试集1与聚类3的测试集1之间差异存在统计学意义（P=0.046），全病变区的测试集2与聚类2的测试集2之间差异存在统计学意义（P=0.023），其余集合间差异无统计学意义（P>0.05）。三个模型的校准曲线和决策分析曲线表明了模型高且相似的预测一致性和临床适用性（图5）。

点击查看大图

图5  不同栖息地区域多模态影像组学模型的决策曲线分析（5A）和校准曲线（5B）。
Fig. 5  Decision curve analysis (5A) and calibration curve (5B) of the multi-modal radiomics model in different habitat areas.

点击查看大图

图4  不同栖息地区域多模态影像组学模型受试者工作特征曲线。AUC：曲线下面积。
Fig. 4  Receiver operating characteristic curves of the multi-modal radiomics model in different habitat areas. AUC: area under the curve.

点击查看表格

表4  多模态模型对于IDH分子状态的预测性能
Tab. 4  Prediction performance of multi-modal model in predicting IDH molecular status

3 讨论

       本研究结合局部熵值的获取和K-means无监督聚类算法，进行了成人型弥漫性胶质瘤异质性的可视化和量化，然后建立和验证了用于预测成人型弥漫性胶质瘤的IDH分子状态的多模态机器学习模型。结果全病变区、聚类2和聚类3的三种模型均显示出了良好的诊断性能（训练集AUC值范围为0.942~0.974，测试集AUC值范围为0.739~0.864）。本研究将局部特征与生境分析联合应用，实现了对肿瘤空间异质性评估的精准描述，这为了解成人型弥漫性胶质瘤的生物学行为提供了新的研究视角。

3.1 图像预处理：离散化与熵值转换

       图像预处理策略的选择是影响研究结果可靠性和可重复性的关键因素^[22]。不同的预处理方法可能导致特征提取结果的显著差异，从而对预测模型的性能产生重要影响。在这项研究中，我们选择通过调节Bin宽以确定Bin数为固定值的离散化策略，Bin数离散化通过将强度值转换为固定数量的离散区间来确保影像特征的可比性和一致性。本研究选择以32作为固定Bin数来调节Bin宽，是因为32作为常用Bin数在大量研究中展现出了良好的可重复性^[22]。XU等^[23]指出，与基于Bin宽度的离散化方法相比，Bin数离散化方法在特征稳定性方面具有明显优势，特别是在多中心研究中能够更好地控制数据变异。在完成图像离散化处理后，我们进一步引入了熵值转换这一关键步骤。熵作为异质性评估的代表性特征，可表征体素信号分布的变化。近年来，熵参数在肿瘤影像分析中的价值得到了广泛验证^{[24, 25]}。在一项荟萃分析^[26]中，熵在预测原发性胰腺导管腺癌方面表现良好，被认证为不良临床预后的指标。尽管熵能捕捉到常规影像可能遗漏掉的微观特征，但其仍然是一个抽象的概念。我们通过构建熵图，将熵值可视化，从而更直观地反映肿瘤内部的异质性分布。在熵图的定量分析中，我们发现不同序列间的高熵体积分布（即高熵所对应的像素体积）有所不同，但整体区间大同小异，这些高熵区域可能对应着肿瘤内部的特定生物学过程，如细胞增殖活跃区、坏死区或血管生成区等。然而，这些影像特征与潜在生物标志物（如基因表达谱、蛋白质组学特征）之间的确切关联机制仍不明确。未来的研究需要结合多组学数据，从分子层面深入阐释熵值变化背后的生物学机制，这将为建立更精准的影像-生物学关联模型提供重要理论基础。

3.2 栖息地成像在成人型弥漫性胶质瘤中的应用

       生境分析的基础在于肿瘤内部的生物学特性与其周围微环境之间存在着密切的相互作用关系。传统的影像分析方法往往将肿瘤视为一个整体^[27]，这可能无法充分捕获肿瘤内部复杂的病理学特征和空间异质性。相比之下，栖息地成像方法通过识别和划分肿瘤内部的异质性区域，允许我们更细致地探索不同生境区域对胶质瘤IDH分子状态的影响。ZHANG等^[28]采用K-means聚类算法，基于T1WI、T2WI及T1-CE三个常规MRI序列，成功将肿瘤划分为2~4个子区域。该方法在预测IDH突变状态方面表现出优异的性能，在训练集和外部验证集中分别获得了0.918~0.994和0.905的AUC值，显著优于传统整体分析方法^[29]（训练集和验证集分别获得了0.885和0.864的AUC值）。VERMA等^[30]进一步扩展了该方法的应用范围，他们从肿瘤的增强区域、肿瘤周围T2-FLAIR高信号区及坏死区的亚区域中提取影像组学特征，成功构建了稳健的无进展生存期预测模型。在本研究中，我们采用T1-CE和T2-FLAIR对肿瘤亚区域进行分割，其中T1-CE序列能够准确反映肿瘤的血供情况和新血管生成状态，这对于理解肿瘤的血管生成机制和代谢特征具有重要价值；而T2-FLAIR序列则能够敏感地检测肿瘤细胞密度变化和周围水肿程度，为评估肿瘤侵袭性提供重要信息。通过整合这两个序列的信息，我们能够获得更全面、更精确的肿瘤微环境特征。

3.3 多模态模型的构建与比较

       本研究中采用十种不同的分类器进行机器学习模型构建，其中逻辑回归（logistic regression, LR）模型及支持向量机（support vector machines, SVM）模型在预测IDH分子状态时表现最为突出，这两种模型在既往的影像组学研究中也展现出优异的诊断性能^[31]。具体而言，LR模型在处理具有线性关系的数据信息时表现出明显优势，而SVM通过寻找最优超平面，使得两个类别的数据点之间间隔最大化，这一特性使得SVM在处理高维、非线性特征空间数据时具有独特优势。在特征选择方面，先前的研究多以提取全局特征为主，KASAP等^[32]使用对比增强T1-CE图像开发全局特征模型，训练集平均AUC值为0.846，测试集平均AUC值为0.842。但其样本量相对较小，且为单中心研究，使得模型的泛化能力受限。现有研究表明一阶特征与肿瘤特征的关联性较高阶特征更为显著，其中能量参数通常显示出最强的放射性组织关联^[33]。本研究发现“均值”“能量”及“峰度”这三个一阶特征对于IDH预测的综合贡献度较高，这表示这些特征可能隐藏着胶质瘤某些生物异质性信息，这一发现与WEI等^[34]研究结果相吻合。他们认为，与IDH突变型胶质瘤相比，IDH野生型胶质瘤通常表现出更强的异质性和更高的细胞密度，其对应的“均值”特征就越低，提示肿瘤存在较高程度的坏死；而“能量”特征更高，这可能反映了肿瘤内部代谢活动的增加；“峰度”作为一个衡量数据分布极端性的指标，IDH野生型胶质瘤由于代谢途径的改变，往往伴随着更显著的血管增生、坏死以及更高的细胞增殖活性，这些病理特征在影像上可能表现为“峰度”参数的异常升高。

3.4 本研究的局限性

       本研究具有几个局限性。首先，本研究虽为多中心研究，但涉及扫描设备及参数存在的差异可能导致部分机器学习模型出现欠拟合现象，从而影响了模型的临床适用性。其次，本研究对数据图像信息只进行了熵值转换，是因为特征转换计算过程涉及大量运算，耗时较长，难以满足临床实际应用的时效性要求。但相较于高阶特征，熵参数的解释性更强。最后，考虑到生境分割的可解释性和泛化难易，本研究采用了基于K-means的无监督聚类算法和2~5个聚类中心，而高维度图像信息的降维和结果与临床的关联性将是未来生境研究的一大难点。

4 结论

       综上所述，本研究通过常规磁共振成像的栖息地局部熵值特征表征了成人型弥漫性胶质瘤的异质性。最终，构成的三种模型（全病变区、聚类2及聚类3）在预测成人型弥漫性胶质瘤IDH分子状态方面均表现出了良好的诊断性能和临床适用性。局部特征和生境分析的联合应用将为非侵入性评估其他多种病理异常区域提供新的思路和方法。