【变异检测】CNV分析和应用
前言
拷贝数变异（copy number variation ,CNV）是指基因组上某些大片段的拷贝数增加或减少，可分为缺失（deletion）和重复（duplication）两种类型。CNV是一种基因组结构变异，可通过改变基因剂量和转录结构等来调节有机体的可塑性，是个体表型多样性和群体适应性进化的主要遗传基础之一。在基因组中，CNV的变异形式主要包括单个片段的倍增、缺失和多次重复，其中，最常见的形式是单个片段重复（segmental duplication, SD）,其在不同拷贝之间的序列同源性大于90%。通过种内比较基因组学分析，可挖掘对表型有显著效应的CNV。目前定义的CNV长度为50 bp~5 Mb的基因组序列的插入或缺失变异。
图1 CNV类型示意图[1]
CNV常见的检测方法主要分为两类：全基因组范围内检测未知CNV和定点检测已知CNV。基因组未知CNV常用的检测方法有芯片法和测序法。芯片法主要包括比较基因组杂交芯片（aCGH）和SNP芯片（single nucleotide polymorphism arrays）。DNA测序法主要包括全基因组测序（WGS）和单分子长读长测序（long reads sequencing）。
研究表明，不少人类复杂疾病，动植物的重要经济性状都和CNV有密切关系。相比单碱基SNP变异，CNV的长度更长，覆盖更大范围的基因序列，可通过改变基因的剂量效应等影响与基因相关的表型。随着全基因组测序准确性的提高、信息量的扩大、成本的下降等，高通量大规模CNV检测得以迅速发展，并且随着GWAS方法的不断完善，为基于CNV-GWAS发掘更多更可靠的重要形状相关的功能CNV提供了条件，可为后续功能CNV的分子机理以及在育种中的应用奠定基础。
下面就对这些分析点进行详细介绍。
常见分析内容
1.CNV检测结果及统计
用柱状图显示不同类型的CNV在不同片段长度的分布情况。如图2所示，滩羊CNV分布大约48.5%的CNV片段在10KB-50KB之间，CNV片段在50KB-100KB之间的占比是21.7%。CNVRs片段在10KB-500KB之间的占比是86.2%。大于1MB的CNVR很少，占比0.8%。
图2 滩羊CNV变异分布统计
2.CNVR全基因组图谱
CNV区域（CNV Region, CNVR）是指不同个体间检测得到的CNV具有一部分重叠区域，把重叠的CNV进行整合，合并成为一个CNVR。如图3所示，将47头牛的CNV进行合并，共确定1043个CNVRs，共覆盖44.63 Mb，约占普通牛基因组序列的2.06%。已知染色体上CNVRs的分布如图3所示，共有702个是缺失类型（Loss），270个是插入类型（Gain），71个是属于复杂类型（Both， CNVR中同时含Gain和Loss）。Loss型CNV的数量约是Gain的2.6倍。就长度而言，最长的CNVR长度为2,111,937 bp，最短的CNVR长度为3,600 bp。Loss型CNV的数量约是Gain的2.6倍。就长度而言，最长的CNVR长度为2,111,937 bp，最短的CNVR长度为3,600 bp。该研究为荷斯坦奶牛绘制了准确的全基因组CNVRs图谱。
图3 CNVRs的全基因组图谱[2]
注：蓝色-插入，红色-缺失，绿色-复杂
图4 CNVRs密度圈图
注：粉色-缺失，紫色-插入，蓝色-复杂，深灰色-SNP，浅灰色-大片段插入
3.群体遗传结构和遗传多样性
为了研究山羊属的系统发生关系，基于全基因组常染色体上的SNP计算了两两个体之间的遗传距离。利用NJ法构建了系统发育树（图5 B），结果显示伊朗野山羊（bezoar）与家羊分成两大枝，并且家羊群体按地理来源聚类。利用ADMIXTURE软件对bezoar和家羊群体进行遗传结构分析（图5 C），当k＝3时，bezoar与家羊群体分为3大谱系；当k＝6时，家羊群体分为4大谱系。为了细化家羊群体内部的遗传结构，利用家羊群体SNP进行PCA分析（图5 D），结果显示：PC1和PC2将家养群体按地理来源分成4个亚群，分别为欧洲家山羊群体（EUR）、非洲家山羊群体（AFR）、西南亚家山羊群体（SWA-SAS）和东亚家羊群体（EAS）。
图5 bezoar和世界家羊样本地理来源分布图和群体遗传结构分析
4.基于CNVs的GWAS分析
基于性状和CNV基因型，应用多种模型进行GWAS分析，完成不同模型的GWAS分析后，通过QQ plot比较不同模型下实际Pvalue与理论Pvalue的分布，判断最优分析结果。获得最优模型后，多重检验校正确定Pvalue的显著性阈值，筛选显著区域，并筛选CNV关联的位点基因。获得候选基因后，进行功能注释和富集分析。
基于仔猪存活的遗传率和CNV进行GWAS分析，如图6，共检测到16个区域，分别在2、3、4、11、12、13、14、15、16和17号染色体上，可解释22.54%的遗传表型。
图6 存活的仔猪数量性状GWAS分析曼哈顿图[3]
在SSC2，SSC3，SSC12和SSC17上确定了GWAS和CNVR分析的显著区域，总共包含56种编码基因。这些基因通过基因剂量的变化影响仔猪存活的数量（表1）。
表1 显著区域的基因注释表
5.选择清除与环境适应性位点挖掘--Vst分析
Vst分析是类似于Fst的一个指标，用来衡量群体间每个CNVR差异大小的统计量，计算方法为Vst＝（Vt－Vs）/Vt，其中Vt表示所有样本该区域拷贝数大小的标准差，Vs表示两个群体各自的标准差根据各自群体大小加权之后的值。Vst的值介于0-1之间，值越大表示群体间该区域拷贝数变异差异越大，反之则越小。如图7所示，CN vs WL的Vst平均值是0.11，CN vs RIR的Vst平均值是0.11，WL vs RIR的Vst平均值是0.15，说明WL vs RIR显示出最高的分化程度，并且筛选Vst＞0.79的CNVRs中的基因，并进行功能注释与富集分析，GO富集分析显示主要涉及器官发育，代谢和免疫调节。该研究有助于了解家禽在CNV水平上的遗传特征，这可能为开发鸡的繁殖提供有用的信息。
图7 不同分群鸡的基因组Vst值曼哈顿图[4]
6.选择清除与环境适应性位点挖掘--Bayescan分析
Bayescan用bayes的方法和Fst的经验分布检测outlier标记的方法，一般应用于地理等分群信息明确的比较组差异位点的选择分析中。本研究中，结合材料的分群信息，使用bayescan进行差异位点的检测，检测完成后使用R脚本进行迭代收敛情况评价和outliers的筛选，并将outliers作为候选CNVRs。Bayescan outliers检测结果如图8所示：
(a)将所有基因座的Fst值与BayeScan的log10（q value）作图。垂直虚线表示log10（q value）的q＝0.05，右边的点表示受选择的显著的CNVRs；
(b)按分群进行BayeScan的差异位点检测和outliers的筛选。
图8 Bayescan outlier检测结果[5]
CNV介绍就分享到这里，后续再分享CNV的综合案例应用。
参考文献
[1]. Clinical implications of copy number variations in autoimmune disorders FAU - Yim, Seon-Hee FAU - Jung, Seung-Hyun FAU - Chung, Boram FAU - Chung, Yeun-Jun. Korean J Intern Med, 2015. 30(3): p. 294-304.
[2] Liu, M., Fang, L., Liu, S. et al. Array CGH-based detection of CNV regions and their potential association with reproduction and other economic traits in Holsteins. BMC Genomics 20, 181 (2019).
[3] Stafuzza, N.B., Silva, R.M.d., Fragomeni, B.d. et al. A genome-wide single nucleotide polymorphism and copy number variation analysis for number of piglets born alive. BMC Genomics 20, 321 (2019).
[4] Seol D, Ko BJ, Kim B, Chai H-H, Lim D, Kim H. Identification of Copy Number Variation in Domestic Chicken Using Whole-Genome Sequencing Reveals Evidence of Selection in the Genome. Animals. 2019; 9(10):809. [5] Kvist L, Honka J, Niskanen M, et al. Selection in the Finnhorse, a native all‐around horse breed[J]. Journal of Animal Breeding and Genetics, 2020.
