0 引言

       世界卫生组织最新的一项报告显示，世界上每8人中就有1人罹患精神障碍^[1]。而在脑疾病的诊断中，大部分依靠主观量表，还未找到一种有效的生物标志物用以揭示其结构与功能机制。近年来，MRI技术的高速发展为我们认识各类脑疾病打开了新的视角。MRI根据成像方式分为功能像与结构像^[2]，传统脑结构研究主要以基于体素的形态学测量（voxel-based morphometry, VBM）为主，常利用灰质体积来表征脑结构特征^[3]或脑疾病病灶^[4]。VBM的应用虽然十分广泛，但其自动归一化算法在处理脑影像数据时可能会改变异常沟回的拓扑结构，导致提取的形态学特征发生错位^[5]。其次，VBM的方法侧重于脑体积计算，无法直接测量复杂的皮层拓扑结构^[6]。因此，近些年的研究更多地开始使用基于表面的形态学分析（surface-based morphometry, SBM）对大脑结构进行分析。SBM采用的皮质表面网格球形配准提高了大脑配准的准确性，该方法使定位更加精准^[7]的同时，将许多大脑内部拓扑结构提升到表面，挖掘出更多深层特征^{[8, 9]}。

       然而既往许多研究仅将SBM的单一指标均值作为特征，不能很好地反映脑区特性^{[10, 11, 12]}。影像组学是一种从医学影像中高通量地提取大量可复现信息的技术^[13]。它可以获得感兴趣区（region of interest, ROI）的成像特征，对其进行全面的可重复评估^[14]，定量并客观地描述异质性^[15]。基于特征的经典影像组学方法^[16]已经被广泛地运用于生物医学工程领域^[17]。如在脑疾病领域，影像组学已被应用于脑胶质瘤的研究，用于预测脑胶质瘤的生存情况^[18]，并基于脑胶质瘤数据建立了生物学功能和脑影像组学特征间的对应关系^[19]。在精神疾病的研究中，影像组学特征也被用于表征脑区间的相互关系^[20]，发掘特定功能所涉及的具体区域^[21]。上述研究说明影像组学作为一种行之有效的特征提取方法，对捕获异质性及各类疾病的预测都具有良好的效果。考虑到SBM和影像组学的优势，两者的结合必将挖掘更深层次的脑部影像特征，其难点在于如何实现三维皮层表面到二维图像的映射。

       综上，本研究将SBM方法与影像组学相结合，以不同的脑区作为不同的ROI，创新性地将表面形态学指数代替图像灰度值作为特征的提取对象，提出了一种将脑图谱从球面映射到图像矩阵的完整方法。此外，为了对特征进行合理选择和分析，利用2D影像组学特征的复测信度评估特征有效性。最后，结合多动症数据集进行多动症症状指标的预测与相关分析，来进一步验证这种脑影像组学特征提取方法的临床价值。以上工作旨在为脑影像的研究提供新的视角与方法，并为临床应用提供参考。

1 材料与方法

1.1 研究对象

       本文所有数据来自于Kennedy Krieger脑成像研究中心的kirby-21数据集（https://rdrr.io/github/adigherman/kirby21.mricloud/man/Kirby-21-MRICloud-Data.html）^[22]与ADHD-200项目（http://preprocessed-connectomes-project.org/adhd200）中北京大学中心的数据。本研究遵守《赫尔辛基宣言》，所有被试均签署了知情同意书，并获得相关伦理机构的批准。

       kirby-21数据集包含了21名健康志愿者（无神经系统疾病史）的数据，年龄22~61岁，他们在两周的时间间隔内完成了两次MRI扫描。ADHD-200多动症数据集中来自北京大学数据中心的研究对象为222例，其中95例为多动症患者，127例为健康被试，年龄7~27岁。

1.2 图像采集

       kirby-21数据集中每个被试两次T1WI图像均使用荷兰Philips Achieva 3.0 T MR扫描仪及配套的8通道相控阵（SENSitivity Encoding, SENSE）头线圈采集。参数如下：分辨率1.0 mm×1.0 mm×1.2 mm，采集时间6 min，TR 6.7 ms，TE 3.1 ms，TI 842 ms，翻转角8°，视野240 mm×256 mm，层数170，SENSE加速度因子2。

       ADHD-200数据集北京大学数据中心的T1WI图像使用德国Siemens 3.0 T MAGNETOM Trio Tim扫描仪及8通道头线圈采集。参数如下：分辨率为1.3 mm×1.0 mm×1.3 mm，采集时间8 min，TR 2530 ms、TE 3.37 ms、TI 1100 ms，翻转角9°，视野250 mm×250 mm，层数128，矩阵256×256。

1.3 MRI数据预处理

       使用基于Linux系统的Freesurfer 5.3对T1WI数据做预处理与配准分割，提取皮质厚度（cortical thickness, CT）、灰质体积（gray matter volume, GMV）和基于表面的皮层平均曲率（cortical mean curvature, CC）、皮层表面积（cortical surface area, CSA）共计四个形态指数。

       主要包括以下几个技术步骤：局部自适应分割、CT和中央表面估计、拓扑校正、球面映射与配准。其中对CT的处理采用基于投影厚度（projection-based thickness, PBT）的方法^[23]，它允许适当处理部分体积信息与沟回模糊，以及没有明确重建的脑沟不对称^[24]。在初始表面重建完成后使用球谐函数进行拓扑校正^[25]和表面细化，最后将中心表面利用球面映射匹配到Freesurfer的"FsAverage"模板上完成配准^[26]。

1.4 最优插值球面映射

       为将三维脑表面特征作为2D影像组学的图像输入，本研究提供一种基于球面局部投影的最优化插值方法。方法主要包括脑区提取、三维曲面到二维平面的降维映射、映射坐标嵌入图像矩阵三个步骤（图1）。

1.4.1 提取脑区

       基于Desikan-Killiany（DK）脑图谱，大脑被划分为68个有效脑区（左右脑各34个）。利用FreeSurfer已经重建于球体表面的DK脑图谱^[27]，对每一个脑区选择自己的中心点，从脑图和球面点集中得到点的三维坐标、所属脑区与近邻关系。遍历共计326 448个有效值点，删除所有非同区邻点关系并整理出新的元胞关系表与最大邻点数。

1.4.2 平面投影

       对每个脑区计算其形心并规定中心为近形心点（距离形心欧氏距离最短），过该点作曲面的切平面。对区内所有点有通过自身与球心的空间直线，该线与切平面的交点即为该点在切平面上的投影点。

       为弥补两点式方程的缺陷（当原坐标位于任何一坐标轴平面，即三维坐标中出现任何0值使得分式无意义），选择性地改变自变量方程的形式，去除投影的零值边界。

       完成空间平面投影后需将三维平面映射到二维平面——这涉及到一个二维直角坐标系的选取：对每个脑区的投影平面规定投影中心（零点）到其最近邻点的向量方向为X轴正方向，选取该向量所确定的平面束中法向量为Y轴的平面作为分界面，任意选择一个除中心和最邻近点外的第三投影点，其所在分界侧面作为Y轴正方向。

       对脑区内所有投影点计算与零点的欧氏距离和与X轴正向量的夹角，通过极坐标系转换得到投影点二维映射的坐标值，利用分界面作为方向的布尔选择器确定y值正负。最后为将二维坐标嵌入矩阵，选择合适的尺度对二维坐标值四舍五入进行预规整，并为映射点规整后的重合记录重合表。

1.4.3 投影坐标嵌入图像矩阵

       预规整虽然极大地提高了嵌入图像矩阵的效率，但由于数据本身并不整齐会造成坐标值由于靠近同一整数而发生挤压，在投影图像任意部分产生空洞点。且球面局部投影的效果总是随着投影面片面积的增加而下降，这导致投影结果中心稠密而边缘疏散，在边缘产生大量疏散空洞（图1）。

       为剥离重合值点、富集边缘疏散点，先计算重合点中二维投影坐标最接近该整数点的作为不动点，其余为剥离点。将重合点的8邻点划分为16个区，根据剥离点所属象限与斜率得到所属扇区，每个邻点（落值点）包括两个扇区。依此选择得到的为最优落值点，按序检查落值点是否为空，为空则剥离，不为空则进行次优匹配直至所有重合值均从整数点上剥离，如图2所示。最后根据原始近邻关系对边缘疏散点作提拉。

点击查看大图

图2  十六区最优匹配。2A：十六区最优匹配的区域划分方式、邻点包括的扇区与匹配方向。2B：本研究最优化插值方法的流程图示意。
Fig. 2  Sixteen-region optimal matching algorithm. 2A: The region division method, the sectors included in the neighbor points and the matching direction of the sixteen-region optimal matching. 2B: The flowchart of the optimal interpolation method in this study.

点击查看大图

图1  最优插值球面映射过程。
Fig. 1  Steps for optimal interpolation spherical mapping.

1.5 特征提取

       使用开源Python软件包Pyradiomics从以CT、GMV、CC与CSA为灰度值的图像矩阵中提取2D影像组学特征。为提高特征提取的精度和数目，需要去除掩膜空洞并对图像做小波变换、梯度变换、平方变换、平方根变换、对数变换与指数变换，基于每一种形式提取一阶特征、灰度共生矩阵特征（gray level co-occurrence matrix, GLCM）、灰度游程矩阵特征（gray level run length matrix, GLRLM）、相邻灰度差矩阵特征（neighbouring gray tone difference matrix, NGTDM）、灰度相关矩阵特征（gray level dependence matrix, GLDM）、灰度大小区域矩阵特征（gray level size zone matrix, GLSZM）。

1.6 特征评估

       基于kirby-21数据集，以DK脑图谱划分的单独脑区为单位，计算21例被试先后两次测试的组内相关系数（intra-class correlation coefficient, ICC）值来评价特征的复测信度，以此筛去不可信特征。ICC借助多因素方差分析中的行变量效应与列变量效应评估一致性^[28]，对本研究需考虑评定者系统误差选择合适度量的Model以确定信度系数的计算方法。形式上，对于给定的测量重复k次实验，ICC计算见式（1）。

       其中MSR是行变量均方（被试者间平方和），MSw是组内变异均方（对象内平方和）。ICC系数范围为[-1, 1]，ICC越接近1，衡量的结果越可靠，小于0.10（或为负）为无一致性。在0.10~1.00的区间内，可靠性得分又可分为：差（ICC≤0.25）、低（0.25<ICC≤0.40）、一般（0.40<ICC≤0.60）、良好（0.60<ICC≤0.75）与优秀（0.75<ICC≤1.00）^[29]。

       根据不同表面形态指数、不同类型影像组学特征、不同脑区的划分，对比复测信度可靠性得分在各级别的占比，并进行统计分析，以此实现对不同组别脑影像组学特征的评估。在不同形态指数下影像组学特征的复测信度分析中，对不同形态指数在复测信度不同可靠性级别的比率，采用one-way ANOVA检验基于不同形态指数的提取是否具有显著差异，事后检验则采用HSD（honestly significant difference）检验来得到具有显著性差异的组别。在不同影像组学特征类型间复测信度差异的分析中，对于同种形态指数下，不同类型的影像组学特征，采用one-way ANOVA（事后分析采用HSD）检验它们在10种滤波方式的ICC均值，得到具有显著性差异的组别。类似地，在不同脑区影像组学特征复测信度差异的分析中，对不同脑区在4种形态指数下各自的ICC均值序列采用one-way ANOVA（事后分析采用HSD）分析得到具有显著性差异的区域。以上统计学方法均以P<0.05为差异具有统计学意义。

1.7 多动症症状指标预测

       为进一步验证本研究所设计的方法在临床方面的应用价值，基于ADHD-200多动症数据，采用评估后的脑影像组学特征进行多动症症状指标的预测。使用复测信度对影像组学特征进行初步筛选后，再分别利用相关系数法、支持向量机的递归特征消除（support vector machine-recursive feature elimination, SVM-RFE）进行特征选择，对两种特征选择方式均选用前15%作为输入。数据集划分上使用留一法（leave-one-out cross-validation, LOOCV）进行交叉验证,最后采用线性回归实现对注意力缺陷指数、过动指数两种多动症症状指标的预测（步骤如图3），线性回归模型对于给定的含有n个元素的数据集，假定因变量Y与自变量X的向量之间存在线性关联，并且这种线性关系中还包括误差变量ε，以解释实际观测与模型预测之间的潜在偏差，见式（2）。

       以向量和矩阵的形式表示上式，可以由Y=Xβ+ε通过最小二乘法得到关于β的解[式（3）]。

       以脑区为单位得到线性回归的预测结果数值序列后，计算和真实结果序列的Pearson相关系数并采用Z检验来检验显著性（以P<0.05为具有统计学意义），进行后续统计分析。

点击查看大图

图3  基于脑区影像组学特征的多动症症状指标预测。
Fig. 3  Prediction of attention deficit hyperactivity disorder (ADHD) symptom indicators based on brain region radiomics features.

1.8 映射方法的对比与准确性

       对比了基于球面局部投影的映射方法与基于主成分分析（principal components analysis, PCA）的降维方法在DK脑图谱脑区映射中的效果，并将嵌入矩阵后重合点的占比作为评价嵌入矩阵难易度的指标，对脑区顶点数量与重合点占比进行相关分析。

       在验证本研究映射方法准确性方面，采用对所有被试每个脑区映射前后CT均值进行相关分析的方式（计算皮尔逊相关系数，并采用Z检验来检验显著性，统计显著性水平设置为P<0.05）来检验表面形态特征的改变程度。

2 结果

2.1 脑影像组学特征的复测信度差异

2.1.1 不同形态指数间的复测信度差异

       由于脑区之间的形态、功能区域和面积差别明显，因此在同一级别的对比中方差均较大，而对于ICC所评估的特征复测信度，基于CT、GMV、CC、CSA的提取的ICC均值分别为0.76、0.83、0.79与0.83，其中基于GMV和CSA提取的影像组学特征相较其他形态指数具有更高的优秀率，与CT与CC组差异具有统计学意义（P<0.05），如表1所示。

点击查看表格

表1  全脑68区关于4个形态指数的ICC多级比率
Tab. 1  ICC multi-level ratios of four morphological indices in 68 brain regions

2.1.2 不同影像组学特征类型间的复测信度差异

       综合全脑指标，应用了多种变换的影像组学特征对于不同形态指数表现出不同的特异性：常用的一阶特征仅在皮层厚度上表现出优秀的可复测性，对比其他类型影像组学特征在所有滤波方式下的平均优秀率，其优秀特征比其他5类特征平均多出10.85%；而对于其他形态指数，一阶特征的优秀率处于平均水平，在GMV、CC、CSA三种形态指数上分别为45.55%、29.58%、39.93%（该三个区域的平均水平为46.64%、29.22%、41.43%，表1）。GLCM特征的性能与一阶特征相仿，仅在基于CT的提取中较为优秀，与其他4类影像组学特征（GLRLR、GLSZM、NGTDM、GLDM）差异具有统计学意义（P<0.05）。GLRLM特征则与它们相反——在除CT以外的形态指数上表现出较好的可复测性，综合对比其他类型影像组学特征在所有滤波方式下的平均优秀率，在GMV、CC、CSA三种形态指数上分别比其他类型特征优秀率平均高出9.91%、8.63%与7.75%。即便对于CT，GLRLM优秀率也比平均水平高2.94%（平均水平为34.64%，表1）。从总体而言，除了基于CT的提取，各种类型影像组学特征的复测信度在另三种形态指数提取上表现的差异不具有统计学意义（P>0.05）。

       从图像变换方式来看，原始图像和梯度变换图像所产生的特征对影像组学特征可复测性的贡献较高，分别比平均水平高出24.44%和14.33%（综合四种形态学指标，影像组学特征的总体优秀率为37.98%）。而对原始图像作指数变换所提取的影像组学特征则普遍不具有可复测性。各形态指数下不同类型影像组学特征的ICC优秀率如图4、表2所示。

点击查看大图

图4  不同类型影像组学特征在ICC评估中的表现。图中x轴分别代表原始图像经过的10类滤波器。ICC优秀的评估标准为ICC>0.75。4A：基于CT的影像组学特征提取；4B：基于GMV的影像组学特征提取；4C：基于CC的影像组学特征提取；4D：于CSA的影像组学特征提取。ICC：组内相关系数；CT：皮层厚度；GMV：灰质体积；CC：平均曲率；CSA：皮层表面积；GLCM：灰度共生矩阵特征；GLRLM：灰度游程矩阵特征；GLSZM：灰度大小区域矩阵特征；NGTDM：相邻灰度差矩阵特征；GLDM：灰度相关矩阵特征。
Fig. 4  Performance of different types of radiomics features in ICC assessment. In the figure, the x-axis represents the 10 types of filters that the original image passes through. The excellent ICC evaluation standard is ICC>0.75. 4A: Radiomics feature extraction based on CT; 4B: Radiomics feature extraction based on GMV; 4C: Radiomics feature extraction based on CC; 4D: Radiomics feature extraction based on CSA. ICC: intra-class correlation coefficient; CT: cortical thickness; GMV: gray matter volume; CC: cortical mean curvature; CSA: cortical surface area; GLCM: gray level co-occurrence matrix; GLRLM: gray level run length matrix; GLSZM: gray level size zone matrix; NGTDM: neighbouring gray tone difference matrix; GLDM: gray level dependence matrix.

点击查看表格

表2  不同类型影像组学特征在不同形态指数下的ICC优秀率
Tab. 2  ICC excellent rate of different types of radiomics features under different morphological indices

2.1.3 不同脑区间的复测信度差异

       对单独脑区而言，右脑提取的影像组学特征与左脑上提取的影像组学特征在可复测性上的表现比较平衡。对于全部四种形态指数，它们在基于左脑第7区（l_entorhinal，左脑内嗅皮层）、左脑第34区（l_temporalpole，左脑颞极）、右脑第7区（r_entorhinal，右脑内嗅皮层）、右脑第33区（r_frontalpole，右脑额极）、右脑第34区（r_temporalpole，右脑颞极）提取的特征所具有复测性相较于其他区显著降低（P<0.05）（图5）。

点击查看大图

图5  不同脑区的影像组学特征在ICC评估中的表现。其中第1区与第5区为空区，标记为灰色（无效区）与蓝色（0值区），因此有效脑区总数为68个。图中有效脑区的色度值由浅到深反映该区ICC均值水平。ICC：组内相关系数；CT：皮层厚度；GMV：灰质体积；CC：平均曲率；CSA：皮层表面积。
Fig.5  Performance of radiomics features of different brain regions in ICC assessment. Area 1 and area 5 are empty areas, marked as gray (invalid area) and blue (0 value area), so the total number of valid brain areas is 68. The chromaticity value of the brain area in the figure depends on the ICC mean of the area. ICC: intra-class correlation coefficient; CT: cortical thickness; GMV: gray matter volume; CC: cortical mean curvature; CSA: cortical surface area.

2.2 多动症症状指标预测

       多动症症状指标预测后进行相关分析的结果显示，左脑额上回、左脑颞上回与左脑海马旁回三个区域与多动症症状显著相关（|r|=0.33~0.52，P<0.05），其中左脑海马旁回的CSA指标表现出与两种症状指标预测的显著负相关，如表3所示。

点击查看表格

表3  多动症症状指标预测中的异常脑区
Tab. 3  Abnormal brain regions in the prediction of ADHD symptom indicators

2.3 其他投影方法与本方法的对比

       PCA降维作为一种依赖于特征的线性组合来构造主子空间的常用降维手段^[30]在小型脑区的效果尚可以接受，但是一旦脑区变大、在球面所占比例增加（脑区内顶点数增加），重合点占比也显著增加（r=0.94，P<0.01），如图6所示，使得PCA的映射效果变差。

       由于PCA是一种无中心的降维方法，可能使投影方向发生反向折叠（如图7圈示部分，PCA对已存在映射值的二维区域重复赋值，使本不相邻的顶点混杂在一起），让随后将二维坐标嵌入图像矩阵的工作变得十分复杂。

点击查看大图

图6  主成分分析（PCA）投影与球面局部投影在不同顶点数量下分别产生的重合点占比。
图7  PCA降维投影与球面局部投影效果。7A：左脑第16区PCA投影；7B：左脑第16区球面局部投影（插值前）。
图8  最优插值球面映射对脑区指标均值的影响。
Fig. 6  The proportion of overlapping points generated by principal components analysis (PCA) projection and spherical local projection at different numbers of vertices.
Fig. 7  The effect of PCA dimensionality reduction projection and spherical local projection. 7A: PCA projection of area 16 in the left brain; 7B: Spherical local projection of area 16 in the left brain (before interpolation).
Fig. 8  Effects of optimal interpolation spherical mapping on the mean values of brain region indices.

2.4 映射准确性

       对于所有被试，对比DK脑图谱进行最优插值球面映射前后每个脑区CT均值的变化，发现映射前每个脑区的CT均值与映射后呈显著正相关（P<0.01），且相关系数接近于1（r=0.97），该结果显示该映射方法较好地保持了脑区原本的形态特征（图8），从一个方面证明了映射的有效性与准确性。

3 讨论

       本研究首次提出了一种基于球面局部映射的最优插值方法实现了DK脑图谱从三维球面到二维图像矩阵的完整、可靠映射，并结合基于表面的形态学分析与影像组学方法完成了对脑区形态特性的表征。通过评估不同形态指数、不同类型、不同脑区下影像组学特征的复测信度，发现基于不同模式提取的特征在可复测性上体现出较大差别。此外，结合多动症数据集实现了对多动症症状指标的预测，表明了该映射方法的临床应用价值。本研究旨在为脑影像的研究提供新的方法和客观参考。

3.1 脑影像组学特征的复测信度差异

       本研究通过评估形态指数、组学特征、脑区三个角度的脑影像组学特征复测信度的差异，发现不同形态指数提取的特征在可复测性上体现出较大差别，其中以GMV与CSA为较好，与基于CT和CC的提取结果存在显著差异，这或许是由于GMV与CSA表现的特征较为简单，而涉及到CT、CC等更为复杂的测量时，特征的可重复性则出现下降^[31]。虽然GMV与CSA的优秀率较高，但在ICC的差级和低级中的占比同样较高（相较CT），这可能是由于这两种形态指数的提取在大脑不同区域存在空间异质性，如在颞下皮质、额极等区域甚至表现出更差的可重复性所致^[32]。

       对于同一形态指数所提取的影像组学特征，以原始图像经过指数变换所提取的效果最差，这可能是因为指数曲线在形态指数值较低的区域斜率较大，而四种形态指数中的数据往往梯度很大且拥有大量高梯度值，由此导致了指数变换的值变异，使得可复测性降低^[33]。而在关于CT的同种变换中，一阶特征与GLCM特征则表现出较好的可复测性，这与JO等^[34]的研究有相近的结果。

       此外，脑影像组学特征在内嗅皮层、额叶、颞叶三个区域的复测信度下降与KNUSSMANN等^[32]的研究高度一致。总体而言，虽然不同表面形态指数、不同类型影像组学特征、不同脑区的划分下脑影像组学特征的复测信度存在差异，本研究设计的表面重建方法所提取的脑影像组学特征均具有较高的可重测性。一方面，这可能是由于脑形态学特征本身具有良好的可重复性，使得脑影像组学特征的可重复性普遍较高^[35]；另一方面，被试个体差异导致形态学特征的异质性，使得脑影像组学特征的可重复性存在差异^[32]。

3.2 潜在临床价值与意义

       本研究将设计的脑影像组学特征提取方法与多动症症状指标预测的临床应用相结合，发现海马旁回、额上回与颞上回与多动症症状显著相关。同类研究中，多动症被认为是一种与认知控制缺陷相关的神经疾病^[36]，在研究与该功能相关的区域时，NORMAN等^[37]发现相较于正常人，多动症患者在额叶、颞叶、顶叶存在大脑结构异常。对于这些区域所涉及的形态指数，皮质表面积的细微差异在儿童中普遍存在，有研究发现多动症患儿的脑表面积值较低主要集中在额叶、扣带回和颞叶区域^[38]。较高水平的多动症症状则与尾中额叶、颞叶和枕叶区域较薄的皮质有关^[39]。此外多种研究共同表明，多动症执行功能损伤与认知表现下降这些行为问题的神经性原因来自于额纹状体和颞顶叶相关区域的活动性降低^{[40, 41, 42, 43]}。上述研究均与本文的结果相吻合，也证明了本研究设计的方法在临床疾病研究中的可行性。同时，本研究设计的映射方法也可迁移至其他脑图谱，为脑影像的研究提供新的方法与思路。

3.3 其他投影方法与本方法的对比

       本研究对比PCA方法与球面局部投影基于DK脑图谱的映射效果，发现随着脑区顶点数目的增多，PCA投影的重合点占比显著增大，对球面局部投影的重合点占比则影响不大。而重合点占比的增加会从很大程度上降低投影坐标嵌入图像矩阵的效率，增加该工作的难度。另外，除了PCA投影，基于近邻结构的流形学习方法局部线性嵌入（locally linear embedding, LLE）也常可作为从三维到二维的映射方法，其由于从样本近邻点关系出发的特点^[44]，对完全定义了近邻关系的DK脑图谱是适用的。局部线性嵌入的方法本质上是求解近邻点的局部线性表示和低维空间中重构的优化问题^[45]，它首先整理一个近邻表，我们可以用DK脑图谱中定义的近邻关系与最大近邻数来代替K-邻近，然后求解表示局部线性结构的权重矩阵^[46]。此处当近邻点数量大于输入维数，定义邻点的矩阵会发生秩亏^[47]，因此为解决正则化问题，改进的局部线性嵌入（modified locally linear embedding，MLLE）加入一个和局部权重矩阵协方差阵的迹相关的任意正则化参数，再作特征值分解求解最小特征值的特征向量实现降维^[48]。即便如此，MLLE的性能也仅在复杂空间拓扑结构上表现得更好，对于相较规则的球面则容易在内部产生形变，拉伸投影区域。

3.4 映射准确性

       为提取2D影像组学特征必须将三维脑表面降维至二维。而在顶点从高维空间投影到低维空间的过程中，特征本身发生了变化，并伴随信息损失。本研究对所有被试每个脑区映射前后CT均值的改变计算了相关系数并验证了显著性，发现二者相关程度很高，说明在投影过程中形态特征的丢失较少。不过，这种方式对映射准确性的检验并不深入，本研究所设计的降维与插值方法对影像特征造成的改变尚需进一步的研究与讨论。

3.5 本研究的不足之处

       本研究存在一定的局限性。首先，本研究纳入的样本量较少，未来将纳入更多被试，进一步验证本研究结果。此外，DK脑图谱中一个点的邻点数最大为6，它们之间的关系以三角面片的形式表述，因此在原本的定义中一个点最多表现为6连通。而为了提取影像组学特征，将有效值点嵌入图像矩阵不可避免地扩展了点的近邻关系，使其最多扩展为8连通。这种连通性的改变既会增加特征的维度，也可能会导致无效影像特征的产生。因此在后续的研究中可以在群论的理论基础上，对这种拓扑变化的同伦群展开讨论和相关研究，评估它的具体后果。

4 结论

       本研究利用大脑表面形态指数和DK脑图谱建立了一个从三维球面到二维图像矩阵完整有效的映射，为提取大脑表面2D影像组学特征提供了一种新的方法。借助组内相关系数ICC评价提取特征的复测性作为筛选冗余特征的指标，并发现不同形态指数在数据一致性上表现出一定的差异。本研究还结合多动症数据集进行了多动症症状指标预测，预测结果与已有研究的成果高度吻合，这也证明了本研究设计的方法的临床应用价值与可行性。这些工作将为理解与预测脑疾病、研究脑结构提供方法学上的新建议。