一种基于三代PacBio和Hi-C技术组装和注释霍巴藏绵羊基因组的方法

技术领域

本发明属于生物信息技术领域，具体涉及一种基于三代PacBio和Hi-C技术组装和注释霍巴藏绵羊基因组的方法。

背景技术

霍巴藏绵羊(Huoba Tibetan sheep)又名霍巴藏羊，是在高寒草原自然生态条件下形成的一个地方原始宝贵遗传资源。其中心产区位于西藏自治区日喀则市仲巴县霍巴乡玉烈村和布穷村，主产区位于雅鲁藏布江源头杰玛央宗冰川区域，海波5000m以上。高寒和低氧是高原地区是主要的生态限制因子，高原土著动物在长期的适应进化过程中形成了独特的低氧适应性策略。藏绵羊是青藏高原的土著反刍家畜和草地生态系统的重要组成部分，是藏族人民重要的生活和生产资料。经过长期的自然选择和进化，藏绵羊对高寒、低氧、强紫外线和冷季营养胁迫等恶劣的自然环境具有极强的适应性能力。虽然功能基因组学研究揭示了青藏高原特有物种适应高海拔缺氧的遗传基础，也鉴定到了一些常见的选择靶标，例如EPAS1、EGLN1和HIF-1α等。但是由于缺乏这些特有物种的参考基因组，其高寒低氧适应的分子机制仍然难以捉摸。

基因组序列组装能够为霍巴藏绵羊的生物学研究提供参考基因序列。染色体水平的基因组在应用层面上具备了更高的实用价值，比较基因组学、转录组学、三维基因组研究工作，都需要高质量的染色体水平参考基因组。目前尚无霍巴藏绵羊的参考基因组，这制约了对霍巴藏绵羊高寒低氧适应的分子机制的研究。PacBio测序是基于光信号的三代测序技术，以单分子测序为特征，可以在在目标DNA分子复制过程中捕获序列信息，因其测序读长较长，对于基因组中复杂区域尤其是复杂结构变异的研究具较大的优势。Hi-C技术源于染色体构象捕获技术，利用高通量测序技术，结合生物信息分析方法，研究全基因组范围内整个染色质DNA在空间位置上的关系，获得高分辨率的染色质三维结构信息。Hi-C技术用于基因组组装，能够将杂乱的基因序列组装到染色体水平。本研究结合三代PacBio和Hi-C技术，并采用二代、三代转录组文库校正拼装序列，首次构建了霍巴藏绵羊染色体水平的高质量参考基因组，推动了霍巴藏绵羊基因组学的研究进展，为霍巴藏绵羊高寒低氧适应的分子机制的研究奠定了基础。

发明内容

针对上述技术问题，本发明提供了一种基于三代PacBio和Hi-C技术组装和注释霍巴藏绵羊基因组的方法，所述的方法包括如下步骤：

(1)采集霍巴藏绵羊血液和组织样本：分别提取霍巴藏绵羊血液和组织样本的DNA和RNA；

(2)构建基因组文库和转录组文库：针对步骤(1)提取得到的DNA分别测序构建二代DNA文库、三代DNA文库和Hi-C文库，获得霍巴藏绵羊基因组文库；针对步骤(1)提取得到的RNA测序构建二代转录组文库和三代全长转录组文库，获得霍巴藏绵羊转录组文库；并对构建所得的基因组文库和转录组文库进行过滤；

(3)采用步骤(2)获得的二代DNA文库评估获得霍巴藏绵羊的基因组大小和杂合率；

(4)基因组组装、纠错和评估：利用步骤(3)获得的霍巴藏绵羊的基因组大小和杂合率结果，选择mecat2软件对步骤(2)所述的过滤后的三代DNA文库进行组装，得到原始组装结果；使用smrtlink 7.0的纠错软件arrow基于步骤(2)得到的三代全长转录组文库对原始组装结果进行纠错，使用pilon软件基于步骤(2)得到的二代转录组文库进行再纠错；

(5)Hi-C辅助组装、纠错和评估；使用ALLHi C软件将步骤(4)获得的纠错后的组装结果进行Hi-C辅助组装，使用Juicer软件构建互作图谱，使用JucieBox软件对其进行可视化纠错，得到霍巴藏绵羊基因组；

(6)基因组注释和评估：重复序列的识别；非编码RNA的预测；编码基因结构预测和功能注释：使用De novo从头预测、homolog同源预测和转录组证据支持3种方法进行编码基因的结构预测，参数设置为默认参数；使用MAKER软件，将上述3种方法预测得到的基因集整合成一个基因集；最后借助于蛋白数据库InterPro、GO、KEGG_ALL、KEGG_KO、SwissProt、TrEMBL和NR对基因集中的蛋白进行功能注释。

优选地，步骤(4)所述的原始组装的mecat2软件的参数设置为：MIN_READ_LENGTH＝10000,CNS_OPTIONS＝"-r 0.6 -a 1000 -c 4 -l 2000"，ASM_OVLP_OPTIONS＝"-n 100-z 10-b 2000-e 0.5-j 1-u 0-a 400"。

优选地，步骤(4)所述的纠错软件arrow的版本为v2.2.2，参数设置为默认参数；pilon软件的版本为v1.22，参数设置为默认参数。

优选地，步骤(1)所述的DNA提取自霍巴藏绵羊血液和肝脏组织；RNA提取自霍巴藏绵羊组织，所述的组织是指心脏、肝脏、肺脏、脾脏、瘤胃和肌肉。

优选地，步骤(6)所述的重复序列的识别是结合基于RepBase库的同源预测方法、基于自身序列比对及重复序列特征的De novo从头预测方法检测重复序列；还利用了TRF软件寻找基因组中串联重复序列；非编码RNA的预测：非编码RNA的注释过程中，根据tRNA的结构特征，利用tRNAscan-SE软件来寻找基因组中的tRNA序列。

优选地，步骤(2)所述的二代DNA文库构建方法为：通过超声波破碎仪将提取获得的霍巴藏绵羊DNA随机打断成长度为300-350bp的片段；DNA片段进行末端修复、加A尾、加测序接头、纯化、PCR扩增；

三代DNA文库构建方法为：使用g-Tubes剪切提取获得的霍巴藏绵羊DNA；DNA片段纯化和浓缩；DNA片段进行末端修复、加测序接头；筛选大小为20kb的目的片段筛选；杂交测序引物和DNA聚合酶绑定；

Hi-C文库构建：使用多聚甲醛固定提取获得的霍巴藏绵羊DNA的构象；限制性内切酶处理交联的DNA，产生粘性末端；DNA末端补平修复，并同时引入生物素，标记寡核苷酸末端；使用DNA连接酶连接DNA片段；蛋白酶消化解除与DNA的交联状态，纯化DNA并随机打断为300～500bp片段；使用亲和素磁珠捕获标记的DNA，对DNA片段进行末端修复、加A尾、加测序接头、纯化、PCR扩增；

二代转录组文库构建：使用带有Oligo的磁珠从提取获得的霍巴藏绵羊RNA中富集mRNA；通过超声波破碎仪将富集到的mRNA随机打断成200bp的片段；以片段化的mRNA为模板，利用随机引物反转录合成一链cDNA，合成第二链cDNA时dNTPs中的dTTP用dUTP代替；cDNA片段进行纯化、末端补平、加A尾、加测序接头；使用USER酶消化二链cDNA，使其文库中只含有一链cDNA；cDNA进行PCR富集，获得霍巴藏绵羊的二代转录组文库；

三代全长转录组文库构建：使用Clonetech SMARTerTM PCR cDNA Synthesis Kit合成mRNA的全长cDNA；使用PB磁珠纯化扩增的全长cDNA，去除1kb以下的小片段cDNA；对全长cDNA进行末端修复，连接SMRT哑铃型接头；进行核酸外切酶消化未连接接头的片段，再次使用PB磁珠进行纯化，获得霍巴藏绵羊的三代转录组文库。

优选地，步骤(2)所述的二代DNA文库测序数据过滤条件为：①去除含有接头序列的reads；②去除重复reads；③当单端测序read中的一端含有的N的含量超过该条read长度比例的10％时，需要去除此对paired reads；④当单端测序read中的一端含有的Q≤5的碱基数超过该条read长度比例的50％时，去除此对reads；

三代DNA文库测序数据过滤条件为：①去除含有接头序列的reads；②去除长度短于1000bp的reads；③去除低质量reads；

Hi-C文库测序数据过滤条件为：①去除含有接头序列的reads；②去除测序read两端连续质量小于20的碱基；③当测序read最终长度小于50bp时，去除此条reads；④仅保留成对reads；

二代转录组测序数据过滤条件为：①去除含有接头序列的reads；②去除3’端；③去除低质量reads；

三代全长转录组文库测序数据过滤条件为：①使用SMRTlink软件，参数设置为最短Subreads长度＝50，最大Subreads长度＝15,000，最小测序循环数＝3，最低预测准确性＝0.99，对原始测序数据进行预处理，对单分子测序的高质量reads进行拆分得到subreads，同一高质量reads得到的subreads经过自我纠错形成环化一致序列；②通过检测嵌合体序列、5’和3’端测序引物，对环化一致序列进行分类，找出全长非嵌合序列用于后续分析。

本发明的第二目的是提供通过所述的方法组装和注释得到的基因组序列。

本发明的有益效果是：1.本发明首次组装了染色体级别的高质量霍巴藏绵羊基因组，并对霍巴藏绵羊的基因结构和基因功能进行了全面注释，这不仅为藏绵羊种群遗传资源保护与利用研究提供了宝贵的基因组资源，而且为进一步研究青藏高原特有畜禽品种的环境适应机制奠定了坚实的基础，也为人类缺氧相关疾病研究提供参考价值。2.本发明的组装方法获得了连续性更好的霍巴藏绵羊参考基因组，为后续进行大规模基因组进化和功能研究提供保障。3.本发明构建组装获得的霍巴藏绵羊基因组质量是现有文献中公开的最高水平，最终确定霍巴藏绵羊基因组顺序及方向确定的染色体长度2.73Gb，contigs N50＝16.66Mb，scaffold N50＝104.79Mb，contig长度锚定率为99.17％，contig数量锚定率为71.83％。PacBio组装和Hi-C辅助组装，并纠错后能完整比对BUSCO的基因占93.50％；4.采用本发明方法鉴定到的重复序列更多，基因注释的完整性更高，注释到的基因数目更加接近于绵羊的平均基因数量。

附图说明

图1霍巴藏绵羊基因组组装路线示意图

图2K-mer深度和K-mer种类数频率分布图

具体实施方式

为使本发明实施例的目的、技术方案和优点更加清楚，下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、完整地说明，应当理解地是，以下实施例仅是本发明的一部分实施例而不是全部的实施例。在下述实施例的基础上，本领域技术人员在没有作出创造性劳动的前提下所获得的所有其他实施例，均属于本发明的保护范围。

名词解释：

以下实施例中所述的“reads”是高通量测序平台产生的序列标签。

以下实施例所述的“contig”是拼接软件基于reads之间的overlap区，拼接获得的序列。

以下实施例所述的“Contig N50”是reads拼接后会获得一些不同长度的Contigs。将所有的Contig长度相加，能获得一个Contig总长度。然后将所有的Contigs按照从长到短进行排序，如获得Contig 1，Contig 2，Contig3，……，Contig 25。

将Contig按照这个顺序依次相加，当相加的长度达到Contig总长度的一半时，最后一个加上的Contig长度即为ContigN50。举例：Contig l Contig 2Contig 3Contig 4＝Contig总长度*1/2时，Contig 4的长度即为Contig N50。ContigN50可以作为基因组拼接的结果好坏的一个判断标准。

以下实施例所述的“Scaffold”是基因组de novo测序，通过reads拼接获得Contigs后，往往还需要构建454Paired-end库或llluminaMate-pair库，以获得一定大小片段(如3Kb、6Kb、10Kb、20Kb)两端的序列。基于这些序列，可以确定—些Contig之间的顺序关系，这些先后顺序已知的Contigs组成Scaffold。

以下实施例所述的“Scaffold N50”，Scaffold N50与Contig N50的定义类似。Contigs拼接组装获得一些不同长度的Scaffolds。将所有的Scaffold长度相加，能获得一个Scaffold总长度。然后将所有的Scaffolds按照从长到短进行排序，如获得Scaffold 1，Scaffold 2.Scaffold 3，……，Scaffold 25。将Scaffold按照这个顺序依次相加，当相加的长度达到Scaffold总长度的一半时，最后一个加上的Scaffold长度即为Scaffold N50。举例：Scaffold l Scaffold 2Scaffold 3Scaffold 4Scaffold 5＝Scaffold总长度*1/2时，Scaffold 5的长度即为Scaffold N50。ScaffoldN50可以作为基因组拼接的结果好坏的一个判断标准。

实施例一、一种基于三代PacBio和Hi-C技术组装和注释霍巴藏绵羊基因组的方法1材料与方法

1.1样品采集

在西藏自治区日喀则市仲巴县霍巴乡玉烈村霍巴藏绵羊中心产区选择健康的成年公羊1只(海拔>5000m)，早晨空腹通过颈静脉采血5ml加入EDTA-K2抗凝剂的采血管中，保存于-20℃冰箱。之后屠宰并立即分割心脏、肝肺、肺脏、脾脏、瘤胃和肌肉组织，迅速用预冷的生理盐水将组织表面的血液冲洗干净，切割成0.5cm

1.2文库构建及测序

采用血液基因组提取试剂盒(天根生化科技北京有限公司，北京)和CTAB(十六烷基三甲基溴化胺)法提取霍巴藏绵羊血液和肝脏组织中的DNA并进行质检，分成三份用于二代、三代和Hi-C文库构建。使用TRlzol Reagent(Invitrogen，美国)提取霍巴藏绵羊心脏、肝脏、肺脏、脾脏、瘤胃和肌肉组织中的总RNA并进行质检，用于二代转录组和三代全长转录组文库构建。

二代DNA文库构建：通过超声波破碎仪(Covaris，美国)随机打断成长度为300-350bp的片段；DNA片段进行末端修复、加A尾、加测序接头、纯化、PCR扩增；使用Qubit 2.0(Invitrogen，美国)和Agilent 2100(Agilent，美国)对构建文库的浓度和片段大小进行质量检测；构建好的文库通过IlluminaHiseq PE150测序。本实施例中的二代DNA文库用于预测基因组大小和杂合率，Illumina Hiseq PE150平台测序结果显示霍巴藏绵羊基因组大小约为2719.60Mb，杂合率为0.27％，基因组杂合率较低，适宜于进行三代DNA文库的构建，并且该结果为后续选择组装策略提供了依据。

三代DNA文库构建：使用g-Tubes(Covaris，美国)将基因组DNA剪切至～20kb大小；DNA片段纯化和浓缩；DNA片段进行末端修复、加测序接头；筛选20kb左右的目的片段；杂交测序引物和DNA聚合酶绑定；使用Qubit 2.0(Invitrogen，美国)和Agilent 2100(Agilent，美国)对构建文库的浓度和片段大小进行质量检测；构建好的文库通过PacBio Sequel II平台进行测序。

Hi-C文库构建：使用多聚甲醛固定DNA构象；限制性内切酶处理交联的DNA，产生粘性末端；DNA末端补平修复，并同时引入生物素，标记寡核苷酸末端；使用DNA连接酶连接DNA片段；蛋白酶消化解除与DNA的交联状态，纯化DNA并随机打断为300～500bp片段；使用亲和素磁珠捕获标记的DNA，对DNA片段进行末端修复、加A尾、加测序接头、纯化、PCR扩增；使用Qubit 2.0(Invitrogen，美国)和Agilent 2100(Agilent，美国)对构建文库的浓度和片段大小进行质量检测；构建好的文库通过Illumina HiSeq PE150测序。

二代转录组文库构建：使用带有Oligo(dT)的磁珠进行mRNA富集；通过超声波破碎仪(Covaris，美国)将富集到的mRNA随机打断成约200bp的片段；以片段化的mRNA为模板，利用随机引物反转录合成一链cDNA，合成第二链cDNA时dNTPs中的dTTP用dUTP代替；cDNA片段进行纯化、末端补平、加A尾、加测序接头；使用USER酶消化二链cDNA，使其文库中只含有一链cDNA；cDNA进行PCR富集；使用Qubit 2.0(Invitrogen，美国)和Agilent 2100(Agilent，美国)对构建文库的浓度和片段大小进行质量检测；构建好的文库通过Illumina Hiseq2500平台进行双端测序。

三代全长转录组文库构建：使用Clonetech SMARTerTM PCR cDNA Synthesis Kit合成mRNA的全长cDNA；使用PB磁珠纯化扩增的全长cDNA，去除部分1kb以下的小片段cDNA；对全长cDNA进行末端修复，连接SMRT哑铃型接头；进行核酸外切酶消化未连接接头的片段，再次使用PB磁珠进行纯化，获得测序文库；使用Qubit 2.0(Invitrogen，美国)和Agilent2100(Agilent，美国)对构建文库的浓度和片段大小进行质量检测；构建好的文库通过PacBio测序仪进行全长转录组测序。

1.3测序数据质控

DNA和RNA测序结束后得到原始数据(Raw reads)，使用Fast QC软件进行过滤获得高质量序列(Clean reads)。

二代DNA测序数据过滤条件为：1)去除含有接头序列的reads；2)去除重复reads；3)当单端测序read中的一端含有的N的含量超过该条read长度比例的10％时，需要去除此对paired reads；4)当单端测序read中的一端含有的低质量(<＝5)碱基数超过该条read长度比例的50％时，需要去除此对paired reads。质控后共获得252761708814bp的cleandata。

三代DNA测序数据过滤条件为：1)去除含有接头序列的reads；2)去除长度短于1000bp的reads；3)去除低质量reads。质控后共获得2752460014bp的clean data。

Hi-C测序数据过滤条件为：1)去除含有接头序列的reads；2)去除测序read两端连续质量小于20的碱基；3)当测序read最终长度小于50bp时，去除此条reads。4)仅保留成对reads。质控后共获得2729684311bp的clean data。

二代转录组测序数据过滤条件为：1)去除含有接头序列的reads；2)去除3’端；3)去除低质量reads。质控后3个肝脏样品平均获得6614583500bp的clean data。

三代全长转录组测序数据过滤条件为：1)使用SMRTlink软件(参数设置为最短Subreads长度＝50，最大Subreads长度＝15,000，最小测序循环数＝3，最低预测准确性＝0.99)对原始测序数据进行预处理，对单分子测序的高质量reads进行拆分得到subreads，同一高质量reads得到的subreads经过自我纠错形成环化一致序列；2)通过检测嵌合体序列、5’和3’端测序引物，对环化一致序列进行分类，找出全长非嵌合序列用于后续分析。质控后具有polyA的全长非嵌合序列为349553903bp。

1.4基因组大小和杂合率评估

获得的高质量序列(Clean reads)，采用基于K-mer的分析方法来估计基因组大小、杂合率和重复序列信息，通过Illumina Hiseq PE150平台测序后获得霍巴藏绵羊二代基因文库，共包括252761708814bp的clean data，测序质量正常，测序错误率正常。选择K＝17，可以产生的K-mer种类数为4

1.5基因组组装、纠错和评估

利用1.4获得的霍巴藏绵羊的基因组大小和杂合率结果，选择mecat2软件得到原始组装结果，参数设置为：MIN_READ_LENGTH＝10000,CNS_OPTIONS＝"-r 0.6-a 1000-c 4-l 2000"，ASM_OVLP_OPTIONS＝"-n 100-z 10-b 2000-e 0.5-j 1-u 0-a 400"，。使用smrtlink 7.0的纠错软件arrow(v2.2.2，默认参数)基于三代转录组文库对原始组装结果进行纠错，使用pilon软件(v1.22，默认参数)基于二代转录组文库进行再纠错。

原始组装并纠错后获得2752460014bp的基因组序列，Contig number为2059，Contigs N50为16585647bp，Scaffoldnumber为606，ScaffoldN50为104790809bp。

组装和纠错后对基因组进行评估，统计A、G、C、T和N在基因组中的占比和GC含量；选择绵羊基因组中CLR(Continuous Long Reads)subreads，使用minimap2软件(默认参数)比对到组装好的基因组，统计reads的比对率、覆盖基因组程度和深度分布情况，由此评估组装的完整性和测序覆盖的均匀性；基于OrthoDB中的单拷贝同源基因集，使用BUSCO预测这些基因并统计其完整度、碎片化程度和可能的丢失率，由此评估整个组装结果中基因区的完整性；用BWA将reads比对到参考基因组，用GATK进行SNP calling并过滤，统计纯合和杂合SNP个数，根据比对结果统计insersize。

1.6Hi-C辅助组装、纠错和评估

利用纠错后的组装结果进行Hi-C辅助组装，得到最终的基因组组装结果。质控得到的clean data使用BWA软件进行比对，使用Lachesis软件将离酶切位点500bp以外的序列去除，得到的数据进行辅助组装。基于顺式互作(同一染色体内的互作)远大于反式互作(不同染色体间的互作)，且顺式互作中线性距离越近则互作越强的原理，将contigs或者scaffolds进行聚类、排序、定向，得到染色体水平基因组。将辅助组装后的基因组，利用Juicer软件构建互作图谱，使用JucieBox软件对其进行可视化纠错。基于OrthoDB中的单拷贝同源基因集，使用BUSCO软件预测这些基因并统计其完整度，碎片化程度及可能的丢失率。Hi-C辅助组装后获得2729684311bp的基因组序列，Contig number为1479，Contigs N50为16656116bp，Scaffoldnumber为27，ScaffoldN50为104790809bp，能完整比对BUSCO的基因占93.50％。

构建流程如图1所示。

1.7基因组注释和评估

重复序列的识别：结合基于RepBase库(http://www.girinst.org/repbase)的同源预测方法(RepeatMasker和RepeatProteinMask)、基于自身序列比对(RepeatModeler)及重复序列特征(LTR-FINDER)的De novo从头预测方法检测重复序列。此外，De novo从头预测方法还利用了TRF软件寻找基因组中串联重复序列。

非编码RNA的预测：非编码RNA的注释过程中，根据tRNA的结构特征，利用tRNAscan-SE软件来寻找基因组中的tRNA序列。由于rRNA具有高度的保守性，因此可以选择近缘物种的rRNA序列作为参考序列，通过BLASTN比对来寻找基因组中的rRNA。利用Rfam家族的协方差模型，采用Rfam自带的INFERNAL软件预测基因组上的miRNA和snRNA序列信息。

基因结构预测和功能注释：使用De novo从头预测(Augustus和Genscan软件)、homolog同源预测(选择C.hircus、H.sapiens、O.aries_rambouillet_v1.0和O.aries_Oar_v4.0作为同源物种)和转录组证据支持(trans.orf/ISOseq)3种方法进行编码基因的结构预测，参数设置为默认参数。使用MAKER软件(默认参数)，将上述3种方法预测得到的基因集整合成一个非冗余的、更加完整的基因集。最后借助于外源蛋白数据库(InterPro、GO、KEGG_ALL、KEGG_KO、SwissProt、TrEMBL和NR)对基因集中的蛋白进行功能注释。

使用BUSCO软件，基于直系同源数据库OrthoDB对基因组组装的完整性进行定量评估。BUSCO抽样了数百个基因组，从中选择单拷贝直系同源＞90％的基因构建了六种主要的系统进化分枝的基因集。

2结果

2.1霍巴藏绵羊基因组大小和杂合率评估

通过IlluminaHiseq PE150平台测序后，共获得252761708814bp的clean data，测序质量正常，测序错误率正常。随机抽取10000对reads数据，通过Blast软件比对NCBI核苷酸数据库(NT库)，比对结果显示文库数据中不含有明显的外源污染，建库测序成功。K-mer分析结果显示，K＝17时可以产生的K-mer总数为225506288297，K-mer深度为81(图2)。因此预估霍巴藏绵羊基因组大小约为2719.60Mb，杂合率为0.27％，重复序列比例为61.57％，基因组GC含量约为42％。

2.2霍巴藏绵羊基因组组装

通过PacBio组装霍巴藏绵羊基因组序列长度为2.75Gb，contigN50＝16.59Mbp，共2059条contigs。通过Hi-C辅助组装，最终确定霍巴藏绵羊基因组顺序及方向确定的染色体长度2.73Gb，contigs N50＝16.66Mb，scaffoldN50＝104.79Mb，contig长度锚定率为99.17％，contig数量锚定率为71.83％(表1)。霍巴藏绵羊基因组GC含量为42.44％(1168317864bp)。PacBio组装和Hi-C辅助组装后能完整比对BUSCO的基因占93.50％。

表1霍巴藏绵羊基因组PacBio组装和Hi-C辅助组装情况

2.3霍巴藏绵羊基因组注释

通过TRF、Repeatmasker、Proteinmask和De novo方法进行霍巴藏绵羊基因组重复序列注释。去掉4种方法的重叠部分结果，重复序列大小为1454718416bp，占霍巴藏绵羊基因组的52.84％。非编码RNA中注释到tRNA、rRNA、miRNA和snRNA个数分别为264021(0.7025％)、185(0.0055％)、528(0.0016％)和2049(0.0084％)。

通过MAKER软件整合De novo从头预测、homolog同源预测和转录组证据支持的结果后，在霍巴藏绵羊基因组中共预测到21072个编码蛋白质的基因，平均基因长度为44224.71bp，平均CDS序列长度为1443.76bp，平均每个基因外显子数目为13.40个，平均外显子长度为214.87bp，平均内含长度为4649.68bp(表2)。

表2霍巴藏绵羊基因组预测结果统计

借助外源蛋白数据库InterPro、GO、KEGG_ALL、KEGG_KO、SwissProt、TrEMBL和NR对霍巴藏绵羊基因集中的蛋白进行功能注释，共注释到20748个编码蛋白质的基因，占上述7种蛋白数据库的98.46％(表3)。使用BUSCO软件进行霍巴藏绵羊基因组注释评估，能完整比对BUSCO的基因有3745个，占比为91.2％。

表3霍巴藏绵羊基因组注释结果统计

2.4与已报道的绵羊基因组比较分析

通过Illumina、PacBio和Hi-C技术组装出了高质量染色体水平的霍巴藏绵羊基因组。本发明所述的霍巴藏绵羊基因组组装中Contigs N50为16.66Mb，ScaffoldN50为104.79Mb，相较于湖羊、朗布依埃羊(Rambouillet sheep)、马可波罗羊(Marco Polosheep)、特克赛尔羊(Texel)、雪羊(Snow sheep)以及摩弗仑羊(mouflon)，序列长度更长，质量明显提高，覆盖率更高。另外，鉴定到的重复序列更多，基因注释的完整性更高，注释到的基因数目更加接近于绵羊的平均基因数量(表4)。

表4各品种绵羊基因组组装和注释效果比较分析