摘要
人类大脑皮层的折叠是一个高度基因调控的过程，它使得我们的颅顶能够容下一个拥有一个更大的表面积的大脑并优化功能组织。脑沟深度是一种可靠但尚未得到充分研究的局部大脑折叠的测量方法，在以往研究中被认为与多种神经发育障碍相关。在这里，我们报告了第一个脑沟深度的全基因组关联研究。通过多变量综合统计（MOSTest）对33,748名英国生物样本库的被试（平均年龄64.3岁；52.0%为女性）进行顶点水平的测量，我们确定了856个于脑沟深度显著相关的全基因组位点(P<5×10−8)。与皮层厚度和表面积的比较表明，脑沟深度具有较高的基因座数量、遗传力和基因有效样本容量。这些指标之间有大量的遗传重叠，基于基因的分析表明它们与神经发育过程有很强的关联。我们的研究结果表明，脑沟深度是一种很有前途的神经影像学表型，可能会增强我们对皮层形态的理解。
一、引言
在早期大脑发育过程中，由于多种内在和外在的机械力的影响，大脑皮层以一种高度可调节的方式折叠成脑回和脑沟。这种皮层折叠不仅允许更大的表面积来适应颅顶，而且还减少了神经元之间的距离，导致更快的信号传输。因此，脑沟形态的测量与认知表现相关，而皮层折叠的缺失（无脑回畸形）伴有严重的智力迟钝。异常的折叠可由神经元增殖、迁移和分化的缺陷引起，并与主要的神经发育和神经退行性疾病有关。
脑沟深度是一种对脑沟形态研究中使用相当少的测量方法，反映了皮层表面任何给顶点的凸性或凹性。因此，这种测量方法非常适合于量化皮层的主要折叠模式，并且它已被证明对小皱纹形式的噪声不敏感，并且在个体中相对稳定。脑沟深度也被认为在整个大脑中具有很高的重测性。使用这一方法的少数研究表明，男性和女性在多个不同区域表现出脑沟深度的差异，并且研究也表明，随着老化，脑沟变得更浅，但这些研究仍有待更完整的描述。
一些研究表明，脑沟深度有很强的遗传基础，这主要是由产前决定的。同卵双胞胎的脑沟比异卵双胞胎更相似，不同被试间的中央沟的平均深度差异有大约56%是由遗传控制。此外，由7号染色体部分缺失引起的威廉姆斯综合征与脑沟深度的广泛减少有关，这介导了其行为症状。然而，目前还没有对这一指标进行大规模的分子遗传学研究。
在这里，我们提供了第一个脑沟深度的全基因组关联研究(GWAS)，将其遗传结构与更常见的大脑皮层厚度和表面积的形态测量进行比较。考虑到基因变异可能在磁共振成像(MRI)表型中具有分布效应，我们针对皮层表面的顶点进行多变量分析，防止了多重比较校正或数据缩减。因此，我们将多元综合统计检验(MOSTest)应用于fsaverage3模板中的1153个顶点的数据，并排除内侧壁的顶点。我们的主要样本包括UK Biobank数据集中的33,748名英国白人被试，平均年龄为64.3岁(SD为7.5岁)，52.0%为女性。我们进一步对另外两个样本进行了分析，其中包括(i)5199名UK Biobank被试（平均年龄62.8岁），由于他们不是英国白人血统，被排除在主要分析数据集之外；(ii)8072名青少年大脑认知发展(ABCD)研究的被试（平均年龄9.9岁）。这些样本在种族和年龄方面与探索样本有很大的不同，这为所报告的关联的广泛性提供了一个强有力的检验。在皮层表面重建之后，我们控制年龄、性别、扫描站点、图像质量和前20个遗传成分对顶点的影响，以控制群体分层效应（群体分层或群体层化指的是在基因研究中，组间或种群间标记位点等位基因频率和易感基因频率差异而导致的假阳性结果）。在应用基于秩的逆正态变换后，对得到的顶点残差测量进行MOSTest，产生与每个单核苷酸多态性(SNP)的多变量关联。此外，我们重复了主要的GWAS分析，同时将所有顶点的平均值作为协变量以消除全局效应，这些分析的结果与主要分析结果高度相似。
二、结果
在探索样本中，MOSTest显示了856个独立的基因位点，全基因组显著性阈值为α=5×10−8（图1A）。相比之下，在皮层表面积和皮层厚度方面，我们分别发现了661个和591个位点。我们发现，已发现的脑沟深度的基因位点集合具有可接受的泛化能力。
我们发现脑沟深度的遗传变异形式遵循脑回和脑沟的空间分布模式，在另外两个样本中遗传效应结果的空间模式高度相似，对15号染色体SNP的效应如图1B所示。
接下来，使用MiXeR工具，我们将高斯混合模型的零模型和非零模型拟合到三个GWAS汇总的统计数据中，估计了多基因性和变异的效应量。结果总结在图1C中，描述了每个指标所发现的SNP所解释的遗传变异的估计比例作为样本容量的函数。不同指标间曲线的水平位移表明，脑沟深度的有效样本容量最高，皮层厚度有效样本容量最低。此外，我们还使用连锁不平衡得分回归(LDSC)应用于每个皮层顶点计算的脑沟深度的平均遗传力，结果表明脑沟深度的遗传力显著高于其他两个指标(见图1D)，即单变量测量也捕获了脑沟深度较高的遗传信号。
图1。基因位点发现。(A)多变量GWAS对沟道深度的曼哈顿图，观察到的每个SNP的−log10(P值)显示在y轴上。x轴表示相关基因位置，按染色体分组，红色虚线表示全基因组显著性阈值（5×10−8）。y轴在−log10(P值)=150处被裁剪。(B)皮层的外侧图，描述了15号染色体在探索样本（顶部）和两个验证样本（中间和底部）中rs4924345的z值。左列对应左半球，右列对应右半球。(C)Power图显示了全基因组显著命中值(y轴)和样本容量(x轴)解释变异之间的关系。图中发现的命中数和解释的遗传变异的百分比显示在图例的括号之间。(D)Bar图每中指标基于SNP的平均遗传力（具有95%的置信区间）。在(C)和(D)中，脑沟深度用绿色表示，皮层表面积用橙色表示，皮层厚度用紫色表示。****P ≤ 0.0001。
此外，我们还通过基因组注释（MAGMA）的多标记分析进行了基于基因的分析。我们发现了2010个在多重比较校正（α=0.05/18,203）后仍然显著的的脑沟深度基因（即所有18,203个基因的11%），而皮层表面积有1486个，皮层厚度有1347个。
接下来，我们在基因位点水平、基因水平和功能通路水平上分析了这三个指标之间的遗传重叠。在基因位点水平上，我们发现脑沟深度和皮层表面积有625个重叠的基因位点（Dice系数为0.82），脑沟深度和皮层厚度有509个重叠的基因位点（Dice系数为0.70），皮层表面积和皮层厚度有450个重叠的基因位点（Dice系数为0.72）。所有3个指标一共有447个重叠的基因位点 (图2A)。我们进一步发现，这些皮层特征的大量基因位点与先前确定的一组大脑障碍和认知表现的基因位点重叠。
在基因水平上，这些指标之间的大量遗传重叠也很明显，如图2B所示。STH基因被认为在tau蛋白的磷酸化中发挥作用，它与这三个指标高度相关。ROBO2、NAV2和SEMA3A是在神经元生长的引导中起关键作用，它们也与这三个指标相关。两个组蛋白基因HIST1H4L和HIST1H2BL都位于6号染色体的大组蛋白基因团簇中，它们与脑沟深度相关；组蛋白活性在调节基因表达模式中起着核心作用，这些基因表达模式决定了大脑皮层形态的增殖和分化过程。
图2C显示了通过MAGMA进行的基因集分析的结果，列出了脑沟深度最重要的10个基因的本体路径及其他两个指标的P值。我们发现这与神经发生和神经元分化途径有很强的关联。神经元切向迁移都与脑沟深度和皮层表面积有关，但与皮层厚度的相关较小，这与神经元的切向迁移在决定皮层折叠中的作用相一致。值得注意的是，与软骨细胞分化和骨骼系统发育相关的功能通路似乎更特定于脑沟深度，这可能暗示了皮层折叠和头盖骨的形成之间的早期生命相互作用。
我们进一步将我们基于基因集分析结果与皮层基因表达模式相结合，这些基因表达模式来自于无显著临床特征的供体的死后大脑。如图2D所示，与这三个指标相关的基因标记探针在整个生命周期中表现出明显的特征，即高的产前表达和低的产后表达。
图2. 基因上的重叠。(A)三个不同指标中发现的基因位点的重叠数量维恩图。(B)基于基因的P值散点图，y轴为脑沟深度P值，x轴为皮层表面积P值，颜色表示皮层厚度P值。请注意，−log10(P值)的最大值为40。(C)脑沟深度的10个最显著的基因功能通路，如y轴所示，x轴上表示三个指标的−log10(P值)。(D)每个指标的基因表达随时间的变化。灰色阴影表示95%的置信区间。
三、讨论
在这里，我们报道了第一个大规模的脑沟深度的分子遗传学研究结果。我们发现了856个基因位点这解释了大约32%的遗传变异，这项研究发现了迄今为止所考虑的所有大脑指标中最多的基因位点。我们进一步表明我们的发现适用于其他群体的数据，并携带有意义的生物学信息。
脑沟深度与皮层表面积和皮层厚度的直接比较表明，脑沟深度更具遗传性。这可能反映了皮层折叠的进化意义，其发展使一个更大的大脑的出现和其功能组织的优化成为可能。研究表明，人类特异性折叠源于机械力和细胞机制之间的相互作用，这是在进化过程中主要通过与细胞周期和神经生长耦合的基因突变产生的。这与已知的神经元增殖和迁移调节中发挥作用的特定基因相一致。因此，我们的研究结果表明脑沟深度指标与这些形成大脑形态的重要特征的遗传过程密切一致。
基因效应对大脑中的一些邻近点有相反的影响方向。这与沿皮层褶皱的神经元和纤维的形态和排列的强烈差别相一致，从沿侧壁的脑回的冠部到脑沟的底部的差异很大。此外，目前的研究结果表明存在广泛的遗传效应说明了其具有遗传结构这一重要特征；其中涉及到许多变异，每一种在整个皮层范围内都有一个复杂的效应模式。正如我们对泛化能力的估计所显示的，这些对脑沟深度的多变量遗传效应似乎在一定程度上独立于人种和年龄，这个结果强调了大脑形态的基本神经生物学相关性。
我们进一步发现，这三种大脑形态指标之间存在大量的遗传重叠，扩展了我们之前的发现，即皮层表面积和皮层厚度共享大部分的遗传决定因素。我们发现脑沟深度的遗传与皮层表面积的重叠大于脑沟深度与皮层厚度的重叠，表明决定皮层折叠程度的神经生物学机制与表面积之间的关系更密切。神经元的切向迁移模式被认为是大脑皮层折叠和其表面积的重要驱动因素，而其厚度受径向迁移的影响更大。功能通路分析的结果似乎证实了这一区别。然而，对重叠的估计以及观察到的表型相关性，也表明脑沟深度和皮层表面积确实在一定程度上捕获了不同的遗传过程。因此，脑沟深度很可能提供关于分子遗传学对大脑形态影响的额外信息，捕捉了决定皮层折叠的机械过程的变化，以补充通过对皮层表面积和皮层厚度的研究所发现的信息。识别的与大脑疾病和认知能力的重叠基因结果还表明，这些大脑指标的遗传学携带临床相关信息。
除了报道的重叠基因位点外，特定的遗传变异、基因和功能通路分析进一步有助于我们了解皮层形态和相关疾病。最重要的功能通路与早期大脑发育特别相关，其中神经生长和神经分化排名最高。这与大量关于驱动皮质折叠的机械力的遗传调控的文献非常吻合。这也与我们的发现相一致，即识别的基因集在胎儿皮层组织中表达最高。此外，皮层折叠已被证明几乎完全发生在产前，出生时的脑沟模式可以预测神经行为结果。我们需要对包括婴儿在内的整个生命周期的神经成像数据进行随访来验证这些发现，并进一步确定脑沟深度遗传效应的时空模式。鉴于脑沟形态的报道与一系列神经发育和神经退行性疾病有关，以及我们的探索发现结果的基因与大脑疾病重叠，探索脑沟深度遗传学与整个生命周期内大脑疾病的发展有何关系将很有意义。
我们注意到一个我们研究中的局限，我们目前缺乏对脑沟深度指标捕获的机制理解。尽管如此，与皮层厚度和皮层表面积相比，皮层折叠的进化和个体生长中的重要性似乎反映在更高的遗传性上。我们对与脑沟深度的多变量遗传关联的探索证实了这一点，为人类皮层的复杂性提供了新的见解，但仍然需要进一步研究以增强我们对大脑和相关疾病的理解。
四、方法
1. 被试
在主要分析中，我们使用了来自UK Biobank数据集的被试数据，这些数据来自注册号为27412的数据存储库。UK Biobank的组成、设置和数据收集协议已经在其他地方被广泛地描述过。UK Biobank已获得美国国家卫生服务机构国家研究伦理服务机构的伦理批准，并获得了其参与者的知情同意与。在本研究中，我们选择了英国白人个体，通过自我报告和遗传主成分分析确定，他们接受了神经成像测量。我们使用了截至2020年3月发布的T1 MRI扫描数据，排除了771个结构扫描质量差的被试，年龄和性别调整的欧拉数比扫描站点的平均值低3个标准差。通过全基因组复杂性状分析(GCTA)，我们进一步排除了每对相关个体中的一个，其阈值为0.0625(n=1138)。最后，分析的样本量为n = 33,748人，平均年龄为64.3岁(SD为7.5岁)，共有52.0%的样本为女性。
对于泛化分析，我们使用了相同的UK Biobank数据和预处理步骤，但将我们的样本限制在那些没有被归类为英国白人的个体。因此，我们共有n=5199个样本，平均年龄为62.8岁(SD为7.7岁)，其中54.1%为女性。
此外，我们还纳入了参与ABCD研究的儿童的数据，完整的遗传数据和基线T1 MRI扫描数据已发布，并通过了ABCD质量控制程序(n=8072)。这些儿童的平均年龄为9.9岁(SD为0.6岁)，其中46.9%为女性。所有程序都得到了加州大学圣地亚哥分校的中央机构审查委员会(IRB)的批准，在某些情况下，还由个别站点的IRBs批准。父母或监护人提供知情同意书，儿童在参与研究前表示同意。
2. 数据预处理
UK Biobank数据集的 T1 MRI扫描是从英国的三个站点收集的，都在相同的配置的西门子Skyra 3T扫描仪上，具有32通道接收线圈。UK Biobank神经成像团队已经发表了大量关于应用扫描协议和程序的信息。我们也使用由西门子Prisma、GE750和飞利浦3T扫描仪从21个站点收集的ABCD数据中获得T1 MRI扫描数据。扫描方案在不同的站点之间进行了协调。
所有扫描都存储在奥斯陆大学的安全计算数据库中。我们应用FreeSurfer v5.3的“recon-all”处理管道，然后在ico3（1284个顶点）和ico4（5124个顶点）分辨率下提取脑沟深度、皮层表面积和皮层厚度的顶点数据，而不使用平滑处理。我们包括了左右半球的测量，并排除了属于内侧壁的皮层顶点。
请注意，我们选择了脑沟深度作为皮层折叠的度量，因为它捕获了顶点方向的局部折叠，提供了与膨胀表面的有方向的距离。
随后，我们从每个顶点测量中回归了年龄、性别、扫描站点、欧拉数和指标的前20个遗传的主成分。在此之后，我们将基于秩的逆正态变换应用于每个指标的残差，从而得到的符合正态分布的指标作为GWAS的输入。
我们在探索性UK Biobank数据中重新运行了如上所述的MOSTest分析，另外回归了三个指标的所有顶点的平均值。
3. MOSTest程序
MOSTest软件可以在https://github.com/precimed/mostest上免费获得。简而言之，假设存在N个变异和M个表型。设zij为第i个变异与第j个表型之间的单变量关联检验的z分数，zi 为M个表型中第i个变异的z分数的向量。设Z={zij}为z分数的矩阵，行为变异，列为表型。对于每个变异，考虑其基因型的随机排列，让={~}为具有排列基因型和表型的变异之间的单变量关联检验的z分数矩阵。对每个变异进行一次基因型的随机排列，并对最终排列后的基因型进行与所有表型的关联测试，从而保持表型之间的相关结构。
设为的相关矩阵，=USVT为它的奇异值分解成分(U和V为正交矩阵，S为对角矩阵，对角线上的奇异值为)。考虑相关矩阵=UVT的正则化形式，其中是通过保持r个最大奇异值用剩余的值替换得到的。然后将第i个变异的MOSTest统计量估计为= zi ~ ziT，其中正则化参数r分别选择皮层表面积和皮层厚度，以最大限度地提高全基因组显著位点的产量。在本研究中，我们观察到r=10对皮层表面积的最大的产量；皮质厚度的最佳选择为r=20，脑沟深度为r=30。
4. 单变量GWAS程序
我们使用了UK Biobank v3的数据。在将BGEN格式转换为PLINK二进制格式后，我们还进行了标准的质量检查程序，包括过滤出缺失超过10%的个体，缺失超过5%的SNPs，以及在P=1×10−9时未通过Hardy-Weinberg平衡检验的SNPs。我们进一步设置了一个次要的等位基因频率阈值为0.005，最终留下了9,061,022个SNPs。
对于ABCD数据库，我们获取了第三次数据发布中的遗传数据，通过精确医学跨组学(TOPMED)估算，并将其从基因组构建hg38映射到hg19。随后，我们对UK Biobank遗传数据应用了相同的计算，留下13,131,314个SNPs。
在排除内侧壁中平均31,153个顶点提供的特征数量后，通过使用MOSTest识别的指标的最大位点数量。我们使用fsaverage3分辨率的1153个顶点的数据输入MOSTest，并计算顶点数据的单变量遗传力。为了可视化结果的区域化，我们以大脑地图的形式展示，另外我们也对fsaverage4的5124个顶点数据进行了单变量GWAS。这样做只是为了提高可视化的分辨率。采用基因型载体和表型载体y之间线性关联的标准加性模型，对每个指标的标准化组进行单变量GWAS。
5. 遗传力
我们通过将LDSC应用于默认设置为的单变量GWAS汇总统计数据中，计算了在fsaverage3分辨率下的每个顶点的基于SNP的遗传力。
6. 基因位点定义
根据精神病学基因组学联盟(PGC)基因位点的定义，从MOSTest汇总统计中确定了显著SNPs和基因位点，这也用于GWAS(FUMA)SNP2GENE的功能定位和注释。首先，我们选择一个通过全基因组显著性阈值5×10−8的SNPs子集，并使用PLINK在连锁不平衡(LD)r2=0.6处执行聚类程序，以确定显著的SNPs列表。其次，我们在LD r2=0.1阈值下的确定显著的SNPs列表，以识别先导SNPs。第三，我们查询了LD的 r2=0.1或更高值的所有候选SNP的参考面板。此外，对于每个先导SNP，其对应的基因组位点被定义为先导SNP染色体的一个相邻区域，包含r2=0.1或更高版本的LD中带有先导SNP的所有候选SNP。最后，如果相邻的基因位点间隔小于250kb，则合并在一起。等位基因LD相关性是从欧洲(EUR)群体1000基因组阶段3数据中计算出来的。此外，我们还根据Enhancing Neuroimaging Genetics联盟使用的定义进行了聚类，以便与以前的成像GWAS研究进行比较。根据这一定义，使用PLINK阈值P=5×10−8和1 Mb的LD截止值为1mb和r2<0.2形成位点。我们利用FUMA在线平台(https://fuma.ctglab.nl/)从MOSTest分析结果中绘制出了显著的SNP。
7. 遗传重叠
不同的MOSTest特征集合之间的遗传重叠被定义为每对汇总统计数据之间在物理上重叠的显著基因位点的数量。如果对基因位点的边界、起始和结束的基因组位置可以通过聚集或重叠来确定，则被认为是物理上重叠的。与全局遗传相关性的指标相反，这种方法可以应用于MOSTest汇总统计，并且对会降低遗传相关性估计的混合效应不敏感。
每对指标间的Dice系数被计算为基因位点重叠数除以两个指标发现的基因位点总数之和。
8. 基因集分析
我们使用默认设置的MAGMA v1.08进行了基于基因的分析，这需要应用基于SNP的平均模型，并使用EUR 1000基因组阶段3的参考。基因集分析是通过将MAGMA应用于基因水平的输出进行的，将正在调查的基因集限制在基因本体生物过程子集(n=7522)的一部分，如分子签名数据库(c5.bp.v7.1)所列出的基因集。
9. 基因表达分析
我们使用了来自56个无临床症状的供体的脑组织的基因表达数据，年龄从受孕后5周到82岁。我们为每个基因选择差异稳定性最高的探针，数量为n=16,660。考虑到大脑皮层样本的表达模式具有相对较高的同质性，我们随后在供体内平均超过13个供体的皮层区域，并对每个探针的基因表达水平值进行跨供体标准化，使得表达值的范围落在0到100之间。
参考文献：The genetic architecture of human cortical folding
