CN112349350B - 基于一种杜氏藻核心基因组序列进行品系鉴定的方法 - Google Patents

基于一种杜氏藻核心基因组序列进行品系鉴定的方法 Download PDF

Info

Publication number
CN112349350B
CN112349350B CN202011238521.2A CN202011238521A CN112349350B CN 112349350 B CN112349350 B CN 112349350B CN 202011238521 A CN202011238521 A CN 202011238521A CN 112349350 B CN112349350 B CN 112349350B
Authority
CN
China
Prior art keywords
dunaliella
genome
strain
sequencing
sequence
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Active
Application number
CN202011238521.2A
Other languages
English (en)
Other versions
CN112349350A (zh
Inventor
高帆
宋韡
南芳茹
冯佳
谢树莲
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Qingdao Aixin Biotechnology Co.,Ltd.
Original Assignee
Shanxi University
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Shanxi University filed Critical Shanxi University
Priority to CN202011238521.2A priority Critical patent/CN112349350B/zh
Publication of CN112349350A publication Critical patent/CN112349350A/zh
Application granted granted Critical
Publication of CN112349350B publication Critical patent/CN112349350B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • G—PHYSICS
    • G16—INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR SPECIFIC APPLICATION FIELDS
    • G16B—BIOINFORMATICS, i.e. INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR GENETIC OR PROTEIN-RELATED DATA PROCESSING IN COMPUTATIONAL MOLECULAR BIOLOGY
    • G16B30/00—ICT specially adapted for sequence analysis involving nucleotides or amino acids
    • G16B30/10—Sequence alignment; Homology search
    • C—CHEMISTRY; METALLURGY
    • C12—BIOCHEMISTRY; BEER; SPIRITS; WINE; VINEGAR; MICROBIOLOGY; ENZYMOLOGY; MUTATION OR GENETIC ENGINEERING
    • C12Q—MEASURING OR TESTING PROCESSES INVOLVING ENZYMES, NUCLEIC ACIDS OR MICROORGANISMS; COMPOSITIONS OR TEST PAPERS THEREFOR; PROCESSES OF PREPARING SUCH COMPOSITIONS; CONDITION-RESPONSIVE CONTROL IN MICROBIOLOGICAL OR ENZYMOLOGICAL PROCESSES
    • C12Q1/00—Measuring or testing processes involving enzymes, nucleic acids or microorganisms; Compositions therefor; Processes of preparing such compositions
    • C12Q1/68—Measuring or testing processes involving enzymes, nucleic acids or microorganisms; Compositions therefor; Processes of preparing such compositions involving nucleic acids
    • C12Q1/6869—Methods for sequencing
    • C—CHEMISTRY; METALLURGY
    • C12—BIOCHEMISTRY; BEER; SPIRITS; WINE; VINEGAR; MICROBIOLOGY; ENZYMOLOGY; MUTATION OR GENETIC ENGINEERING
    • C12Q—MEASURING OR TESTING PROCESSES INVOLVING ENZYMES, NUCLEIC ACIDS OR MICROORGANISMS; COMPOSITIONS OR TEST PAPERS THEREFOR; PROCESSES OF PREPARING SUCH COMPOSITIONS; CONDITION-RESPONSIVE CONTROL IN MICROBIOLOGICAL OR ENZYMOLOGICAL PROCESSES
    • C12Q1/00—Measuring or testing processes involving enzymes, nucleic acids or microorganisms; Compositions therefor; Processes of preparing such compositions
    • C12Q1/68—Measuring or testing processes involving enzymes, nucleic acids or microorganisms; Compositions therefor; Processes of preparing such compositions involving nucleic acids
    • C12Q1/6876—Nucleic acid products used in the analysis of nucleic acids, e.g. primers or probes
    • C12Q1/6888—Nucleic acid products used in the analysis of nucleic acids, e.g. primers or probes for detection or identification of organisms
    • C12Q1/6895—Nucleic acid products used in the analysis of nucleic acids, e.g. primers or probes for detection or identification of organisms for plants, fungi or algae
    • G—PHYSICS
    • G16—INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR SPECIFIC APPLICATION FIELDS
    • G16B—BIOINFORMATICS, i.e. INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR GENETIC OR PROTEIN-RELATED DATA PROCESSING IN COMPUTATIONAL MOLECULAR BIOLOGY
    • G16B20/00—ICT specially adapted for functional genomics or proteomics, e.g. genotype-phenotype associations
    • G16B20/30—Detection of binding sites or motifs
    • G—PHYSICS
    • G16—INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR SPECIFIC APPLICATION FIELDS
    • G16B—BIOINFORMATICS, i.e. INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR GENETIC OR PROTEIN-RELATED DATA PROCESSING IN COMPUTATIONAL MOLECULAR BIOLOGY
    • G16B30/00—ICT specially adapted for sequence analysis involving nucleotides or amino acids
    • G16B30/20—Sequence assembly
    • C—CHEMISTRY; METALLURGY
    • C12—BIOCHEMISTRY; BEER; SPIRITS; WINE; VINEGAR; MICROBIOLOGY; ENZYMOLOGY; MUTATION OR GENETIC ENGINEERING
    • C12Q—MEASURING OR TESTING PROCESSES INVOLVING ENZYMES, NUCLEIC ACIDS OR MICROORGANISMS; COMPOSITIONS OR TEST PAPERS THEREFOR; PROCESSES OF PREPARING SUCH COMPOSITIONS; CONDITION-RESPONSIVE CONTROL IN MICROBIOLOGICAL OR ENZYMOLOGICAL PROCESSES
    • C12Q2600/00—Oligonucleotides characterized by their use
    • C12Q2600/156—Polymorphic or mutational markers

Landscapes

  • Life Sciences & Earth Sciences (AREA)
  • Chemical & Material Sciences (AREA)
  • Proteomics, Peptides & Aminoacids (AREA)
  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • Health & Medical Sciences (AREA)
  • Analytical Chemistry (AREA)
  • Organic Chemistry (AREA)
  • Biotechnology (AREA)
  • Bioinformatics & Cheminformatics (AREA)
  • Biophysics (AREA)
  • General Health & Medical Sciences (AREA)
  • Zoology (AREA)
  • Wood Science & Technology (AREA)
  • Medical Informatics (AREA)
  • Theoretical Computer Science (AREA)
  • Spectroscopy & Molecular Physics (AREA)
  • Molecular Biology (AREA)
  • Evolutionary Biology (AREA)
  • Genetics & Genomics (AREA)
  • Bioinformatics & Computational Biology (AREA)
  • Immunology (AREA)
  • General Engineering & Computer Science (AREA)
  • Microbiology (AREA)
  • Biochemistry (AREA)
  • Botany (AREA)
  • Mycology (AREA)
  • Measuring Or Testing Involving Enzymes Or Micro-Organisms (AREA)

Abstract

本发明属于植物分子鉴定技术领域，具体涉及基于一种杜氏藻核心基因序列进行品系鉴定的方法。该方法主要包括：样本搜集、纯化与培养；全基因组DNA提取；构建DNA测序文库；获取待测藻株和杜氏藻Dunaliella quartolecta全基因组测序数据；杜氏藻D.quartolecta核心基因组测序片段筛选与从头组装，对组装的核心基因组序列进行基因组分、蛋白功能注释及基因组重叠群共线性分析；利用单核苷酸多态性构建***进化树，当待测藻株与四叶杜氏藻聚为一簇，且分支的数据支持率在0.99～1.00，遗传相似度百分比大于≥99％，待测藻株即为D.quartolecta。

Description

基于一种杜氏藻核心基因组序列进行品系鉴定的方法
技术领域
本发明属于植物分子鉴定技术领域，具体涉及基于一种杜氏藻核心基因组序列进行品系鉴定的方法。
背景技术
杜氏藻Dunaliella quartolecta是一种生活在海洋、盐湖及其它极端环境的真核单细胞微藻，属绿藻门、绿藻纲、团藻目、盐藻科、杜氏藻属，抗逆性强，无细胞壁，含色素体和蛋白核，细胞顶端具鞭毛。杜氏藻D.quartolecta富含甘油、β-胡萝卜素、藻多糖等生物活性物质，属于特色经济微藻。以杜氏藻D.quartolecta中的特色品系作为生物反应器，对其活性物质进行提取及工业化生产，在食品加工、医疗保健、生物柴油等领域具有重要的应用前景。然而，目前国内外已被鉴定的杜氏藻属共23个种，形态相似且广谱耐盐，单从形态学角度对其中的杜氏藻D.quartolecta进行鉴定难度很大。从DNA标记、基因标记、蛋白标记的角度虽提高了藻株鉴定的效率，然而准确率仍然受到分子标记手段、片段的保守性及扩增或实验程序的非通用性等因素限制，一些近缘藻株的常规分子鉴定常会遇到候选扩增片段少、通用标记特异性差、新型标记及特异性引物开发周期长、PCR扩增程序仍需优化等缺点，获得的鉴定结果也往往存在假阳性。作为杜氏藻属中一类重要的高附加值特色品系，杜氏藻D.quartolecta资源的分子鉴定十分关键。因此，研发一种更精准、快速、通用的杜氏藻D.quartolecta分子鉴定方法很有必要。
由于下一代DNA测序技术的快速发展，基于物种全基因组水平的分子鉴定技术成为可能。较传统的分子鉴定技术，全基因组水平的鉴定遗传信息量更大，检测范围更广，对于近缘物种的鉴定更有效，可获取的遗传变异信息更丰富。目前，很多模式物种的全基因组测序数据已经公布。虽然盐生杜氏藻(D.salina)的参考基因组测序数据已于2017年公布(Dunsal1 v.2)，然而，作为另一种典型的嗜盐杜氏藻D.quartolecta，迄今仍未有该藻株全基因组测序工作的相关报道。利用当前流行的二代联合三代测序技术对物种进行全基因组测序，虽然能获得该物种较为完整的遗传信息，但仍存在以下几方面缺陷：(1)所有测序片段均须进行完全比对，运算耗时较长，数据量产出庞大，会消耗计算机大量的时间和资源，不利于分子鉴定工作的及时开展；(2)基因组组装及生物信息分析不仅高度依赖国内外测序公司的二代和三代高通量测序平台，如Illimina、Nanopore、PacBio等，而且受物种基因组大小及平台计算能力的限制，结果产出周期较长，造价较高，普通实验室往往难以承担；(3)对近缘物种进行分子鉴定，将高度依赖其全基因组重测序质量，而这又与参考物种基因组质量紧密相关，若参考物种基因组测序深度不够、组装质量不高，将影响待测物种基因组的重测序结果，进而导致物种鉴定出现偏差。
因此，如何提供一种准确、高效、经济的从待测藻株中鉴定杜氏藻D.quartolecta的方法是本领域亟待解决的技术问题。
发明内容
针对上述问题本发明提供了基于一种杜氏藻核心基因组序列进行品系鉴定的方法。
为了达到上述目的，本发明采用了下列技术方案：
基于一种杜氏藻核心基因组序列进行品系鉴定的方法，包括以下步骤：
(1)样本搜集、纯化与培养：采集待测藻株及杜氏藻D.quartolecta，将待测藻株经纯化后进行室内扩大培养；
(2)全基因组DNA提取：利用改良的CTAB法分别提取待测藻株及杜氏藻D.quartolecta的全基因组DNA，冷冻保存；
(3)将步骤(2)中的待测藻株及杜氏藻D.quartolecta的全基因组DNA打断、纯化后分别构建DNA测序文库；
(4)采用高通量测序法对步骤(3)中的DNA测序文库分别进行测序，获取待测藻株和杜氏藻D.quartolecta全基因组二代测序数据；
(5)以NCBI已公布的盐生杜氏藻全基因组数据为参考，将步骤(4)中获取的杜氏藻D.quartolecta全基因组测序数据与其进行比对，通过筛选、从头组装、质量评估后获得杜氏藻D.quartolecta核心基因组序列，该核心基因组序列大小为6592916bp，重叠群数量为3000个，最大重叠群长度为1133322bp，重叠群平均长度为2197.64bp，重叠群N50为15270，完整基因占比23.65％，单拷贝基因占比15.18％，多拷贝基因占比13.76％，空位/缺失占比1.89％，不完整片段占比17.45％，构建从头组装的杜氏藻D.quartolecta核心基因组环状图谱，然后对杜氏藻D.quartolecta核心基因组序列进行基因组分、蛋白功能注释及基因组重叠群共线性分析；
(6)以步骤(5)中构建的杜氏藻D.quartolecta核心基因组序列为参考，将步骤(4)获得的待测藻株全基因组测序数据及已公布代表性藻类的基因组测序数据与其进行比对，检测物种间的单核苷酸多态性和***/缺失位点，然后利用单核苷酸多态性构建***进化树，当待测藻株与杜氏藻D.quartolecta聚为一簇，且分支的数据支持率在0.99～1.00，遗传相似度百分比大于≥99％，待测藻株即为杜氏藻D.quartolecta。
进一步，所述步骤(1)中室内扩大培养具体步骤为：对待测藻株藻细胞进行无菌条件下的单克隆挑取，显微镜检合格后在无菌条件下进行室内扩大培养，室内扩大培养条件为：光周期为18h：6h，光照强度为19000lx，温度：23±3℃，保持无菌通风环境，每隔5天对培养皿进行摇动以防藻细胞贴壁，并取0.5～1mL藻液进行镜检，培养周期为28±7天，配制下述培养基溶液对待测藻株进行室内扩大培养，培养基配方如下：
30g/L NaCl，1.5g/L NaNO3，1.4g/L K2HPO4，1.75g/L MgSO4·7H2O，1.36g/LCaCl2·7H2O，1.2g/LNa2CO3，0.006g/L FeC6H5O7，0.005g/LNaH2PO4·2H2O，0.5g/LCo(NO3)2·6H2O，0.8g/LCuSO4·5H2O，2.3g/LZnSO4·7H2O，0.03g/LH3BO3，4.0g/LNa2MoO4·2H2O，0.02g/LMnCl2·4H2O，0.5g/LVB1，0.5g/LVB12，VH 0.5g/L，超纯水定容至1L。
进一步，所述步骤(2)中改良的CTAB法具体步骤为：取600～800mg待测藻株，用超纯水冲洗2～3次，4℃8000r/min离心1.5min，加液氮研磨15sec，加入800μL 20℃预热的2％W/V的CTAB溶液，1μL 1％V/V的β-巯基乙醇，混匀后在60℃水浴1.5h，期间每隔20min摇匀1次，加入800μLTris饱和酚，混匀后4℃12000r/min离心2.5min后取上清，加入体积比为25：24：2的Tris饱和酚、氯仿和异戊醇混合液，涡旋振荡后4℃静置10min，期间混匀2～3次，加入800μL 0.1％V/V DEPC处理的ddH2O，60℃水浴30min，4℃12000r/min离心4min后取上清液，加入150mL 3mol/L的乙酸钠及250mL 4～5℃预冷的无水乙醇，-20℃沉淀50min，4℃10000r/min离心3min后弃上清液，加入1mL 4～5℃预冷的70％V/V乙醇溶液并涡旋振荡20sec，弃上清后在核酸真空干燥***中挥发液体，加入100×TE缓冲液溶解沉淀以保证DNA浓度≥150ng/μL，1％W/V的琼脂糖凝胶电泳联合荧光定量仪检测基因组DNA，确保其电泳条带明亮、无降解，OD260/OD280在1.8～1.9，无污染。
进一步，所述步骤(3)中构建DNA测序文库具体步骤为：全基因组DNA用80～100W的强档超声波段打断6sec，每隔3sec重复1次，共超声5次，设置打断参数为300～400bp；片段经琼脂糖凝胶电泳后，胶回收300～400bp目的片段；用硅基磁珠吸附回收目的片段，再进行荧光定量仪检测吸附回收目的片段的质量；DNA末端修复，3’末端加A；加接头进行连接反应，对连接产物进行纯化、转化、PCR验证；阳性产物经95℃变性20sec后，进行单链DNA环化反应，纯化产物后构建得到可用于上机的全基因组DNA测序文库。
进一步，所述步骤(5)中通过筛选、组装、质量评估后获得杜氏藻D.quartolecta核心基因组序列的具体步骤为：从测序平台筛选获得高质量序列，筛选测序深度在50～80×，平均长度在12～15K，N50长度大于18K的片段作为查询序列，利用SOAPaligner或BWA软件将其回帖到已报道的盐生杜氏藻参考基因组(Dunsal1 v.2)上，进一步筛选序列一致度≥90％，比对结果E值小于1e-10的测序片段作为杜氏藻D.quartolecta核心基因组序列候选数据；所有剩余测序片段与该候选数据集进行比对，获取比对数据间的重叠区；利用Falcon或Pilon软件对比对结果进行纠错及校正运算，利用SOAPde novo 2.04、Mecat、HERA或Canu软件组装重叠群；利用BySS 2.2.3、Velvet 1.2.10或ABySS 2.2.3软件确定各重叠群顺序；利用BAMStats或GATK DepthOfCoverage软件进行全基因组覆盖度测算，筛选参考基因组覆盖度≥50％，重叠群连续排列数≥2000的核心序列；利用BUSCO 2.0或Quast软件对筛选出的重叠群组装质量进行评估，选择完整基因占比≥20％，单拷贝基因占比15％，多拷贝基因占比≥12％,缺失/空位占比≤3％的组装序列作为杜氏藻D.quartolecta核心基因组序列；利用Circos软件构建该物种的核心基因组环状图谱。
进一步，所述步骤(5)中对杜氏藻D.quartolecta核心基因组序列进行基因组分、蛋白功能注释及基因组重叠群共线性分析，具体步骤如下：利用Augustus3.3.3、ESTScan3.0.1、TransDecoder 2.0.1或Prodigal 2.6.1软件对组装数据进行CDS预测，利用Repeatmasker 4.0.9、RepeatProteinMask 3.2.2、LTR-FINDER、Piler 1.0.6或RepeatScout 1.0.5软件对组装数据进行重复序列分析，利用Diamonds 0.9.14或BLASTX软件将CDS编码的蛋白序列比对到NR数据库，并进行功能注释，预测的蛋白序列经BLASTp自身比对后，再利用MCScanX、Last、Mugsy、Spines或Progressivemauve软件进行基因组的共线性分析。
进一步，所述步骤(6)中利用单核苷酸多态性构建***进化树的具体步骤为：利用LASTZ 1.02.00或Mauve 2.3.1软件，将待测藻株及5～6种已在NCBI数据库中报道的代表性藻类基因组数据分别与所述步骤(5)中组装的杜氏藻D.quartolecta核心基因组序列进行比对，根据比对的共线性块结果提取每个物种与杜氏藻D.quartolecta核心基因组的对应基因型，之后以杜氏藻D.quartolecta核心基因组为模板，将所有物种的基因型信息进行合并、提取、过滤，利用BWA0.7.17软件检测其中的单核苷酸多态性数据和***/缺失位点数据；基于单核苷酸多态性数据，利用EasySpeciesTree 1.0、MEGA 5.0、TreeBeST 1.9.2、PHYLIP、Puzzle 5.2或PHYLO-WIN软件中的最大似然算法构建***进化树，进而确定待测藻株与杜氏藻D.quartolecta间的亲缘关系。
进一步，所述过滤的缺失率不高于20％。
本发明方法不完全依赖于已知杜氏藻全基因组测序结果，对未公布基因组测序数据的近缘品系——杜氏藻D.quartolecta的基因组进行测序，利用优化的数据比对方法及序列组装策略，规避并克服了传统基因组测序中数据完全比对运算耗时，高度依赖先进测序***平台，造价昂贵等缺点。操作者从国内测序公司获得二代测序数据后可依据本发明构建的基因组核心序列及程序命令自行进行测序数据处理、组装、信息分析，各步骤可选用软件范围广、实例中的程序设置严谨，易于上机操作，在杜氏藻品系分子鉴定、变异检测及***进化分析等方面具有广阔的应用前景。
在国内外尚未公布杜氏藻D.quartolecta全基因组测序数据的基础上，本发明首次构建了杜氏藻D.quartolecta核心基因组组装序列，该序列包含了目前遗传信息量最丰富、组装质量较高的D.quartolecta核心基因组信息，以其为参考，将为该藻株的遗传定向改良及其工业化应用提供理论与信息支持。
与现有技术相比本发明具有以下优点：
1.本发明利用二代测序联合基因组从头组装技术首次构建了杜氏藻D.quartolecta核心基因组序列，该序列包含目前遗传信息量最丰富、组装质量较高的D.quartolecta核心基因组信息，弥补了该物种基因组信息的空白。
2.本发明构建的杜氏藻D.quartolecta核心基因组序列可应用于该藻类品系的分子鉴定，在大大提高该类藻株精准鉴定效率的同时，还将为国内外杜氏藻属的***发育与进化研究鉴定理论与技术基础。
3.相较于已公布的盐生杜氏藻D.salina全基因组序列，本发明构建的杜氏藻D.quartolecta核心基因组数据量较小，以其为参考序列进行待测藻株基因组测序数据分析，不仅可以大大缩短数据比对时间，提高了待测藻株有效单核苷酸多态性(SNP)数据获得效率，而且对于基因组水平杜氏藻近缘品系的遗传变异分析具有重要的参考价值，为低等藻类，尤其是绿藻的起源与进化的***研究提供了丰富的数据基础。
4.以本发明构建的杜氏藻D.quartolecta核心基因组序列为参考，依据研究者不同的实验目的，设置相应的实验组和对照组，或以该藻株及其近缘品系相比较，从中挖掘差异或特色基因，这为从分子水平改良研究藻株品质，促进其工业化应用奠定了基础。
5.本发明所用到的杜氏藻D.quartolecta及待测藻株的室内扩大培养方法、改良的CTAB法、核心基因组测序数据筛选及测序片段从头组装的方法可广泛应用于藻类，尤其是绿藻植物的人工培养、高质量全基因组DNA提取及基因组测序数据优化处理等方面，较传统方法实验周期更短、效率更高、易操作，是一套可间接复制的技术方法。
附图说明
图1从头组装的杜氏藻D.quartolecta核心基因组环状图谱，图谱最外层为核苷酸序列大小坐标(单位：Mbp)，其内侧为基于序列一致度(相对于参考基因组Dunsal1 v.2)排列的从头组装片段，基因组片段内部的线条代表各类型基因位点，最内侧为对应的重叠群测序丰度图谱，圈图内部为该藻核心基因组基本信息；
图2是一株待鉴定藻株(暂命名为Dunaliella sp.)室内扩大培养30天后的形态学观察结果，上方为宏观条件，下方为显微条件(比例尺：50μm)，从左至右依次为该藻的1～4号样本；
图3是1％琼脂糖凝胶电泳检测待鉴定样本全基因组DNA示意图，M1和M2代表DNAladders；
图4是杜氏藻D.quartolecta核心基因组与待鉴定藻株基因组测序片段间的共线性分析散点图，图中的点代表两物种基因组间的共线性块，图中A和B分别代表杜氏藻D.quartolecta与待鉴定藻株基因组间2个密集分布的共线性区域；
图5是基于单核苷酸多态性(SNP)数据构建的7个不同藻类间的***进化树，进化树构建算法为最大似然法，步长值设为1000，各分支节点间的数据分别代表支持率和遗传相似度百分数；
图6是一株已鉴定杜氏藻株Dq_SX核心基因组内共线性分析圈图，圈图内的各区段间的连接线代表该物种基因组进化过程中可能发生的加倍事件，圈图上的数字代表核心基因组重叠群编号；
图7是已鉴定杜氏藻株Dq_SX的Ka/Ks值频率分布图，柱状图上的数据代表不同区间内的频率值，Ka代表核苷酸非同义替换率，Ks代表核苷酸同义替换率；
图8是一株已鉴定杜氏藻株Dq_SX核心基因组中蛋白COG注释信息统计，COG即直系同源蛋白数据库，该柱状图统计了同源蛋白注释信息排前20位(top20)的功能信息；
图9是已鉴定的杜氏藻株Dq_SX中一种转录调节因子的跨膜结构域预测图，不同的线条分别代表膜所在区域，膜内区域和膜外区域，纵轴代表该区域预测的概率值，横轴代表氨基酸位置；
图10是已鉴定杜氏藻株Dq_SX中一种转录调节因子信号肽结构预测图，C-score、S-score和Y-score分别代表剪切位点评分，信号肽评分及综合评分值；
图11是杜氏藻D.quartolecta与藻株Dq_SX的代谢通路文氏图，交集部分为两藻株间共有的代谢通路，两藻株代谢通路预测均基于KEGG数据库，即日本京都基因和基因组百科全书进行；
图12是杜氏藻D.quartolecta中特有的排名前20(top20)代谢通路富集气泡图，代谢通路信息来自于KEGG，即日本京都基因和基因组百科全书数据库，气泡体积越大代表该通路中参与的基因数越多，气泡颜色越深代表该通路的置信度越高(Q值越低)，富集程度(显著性)用富集比率表示，富集比率＝KEGG通路注释到的基因数/基因总数；
图13是已鉴定藻株Dq_SX中特有的排名前20(top20)代谢通路富集气泡图，代谢通路信息来自于日本京都基因和基因组百科全书数据库(KEGG)，气泡体积越大代表该通路中参与的基因数越多，气泡颜色越深代表该通路的置信度(Q-value越低)越高，富集程度(显著性)用富集比率表示，富集比率＝KEGG通路注释到的基因数/基因总数；
图14是杜氏藻D.quartolecta显著富集代谢通路中排名前20(top20)的GO富集分析柱状图，GO即基因本体联合会建立的数据库，GO条目越多，对应-log10(Q值)越高(置信度越高)，代表基因参与该生物学功能的程度越高；
图15是已鉴定藻株Dq_SX中显著富集代谢通路中排名前20(top20)的GO富集分析柱状图，GO即基因本体联合会建立的数据库，GO条目越多，对应-log10(Q值)越高(置信度越高)，代表基因参与该生物学功能的程度越高；
图16为基于21株杜氏藻ITS基因构建的***进化树，进化树构建算法为最大似然法，步长值设为1000，各分支节点间的数据分别代表支持率和遗传相似度百分数；
图17为基于21株杜氏藻SSR标记构建的***进化树，进化树构建算法为最大似然法，步长值设为1000，各分支节点间的数据分别代表支持率和遗传相似度百分数；
图18为基于21株杜氏藻基因组SNP构建的***进化树，进化树构建算法为最大似然法，步长值设为1000，各分支节点间的数据分别代表支持率和遗传相似度百分数。
具体实施方式
下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例仅仅是本发明一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例，都属于本发明保护的范围。
实施例1
一种杜氏藻D.quartolecta全基因组测序及其核心基因组序列片段从头组装的方法，包括以下步骤：
步骤1，对一株杜氏藻D.quartolecta的藻细胞进行无菌条件下的单克隆挑取，显微镜检合格后在无菌条件下进行室内扩大培养，室内扩大培养条件为：光周期为18h：6h，光照强度为19000lx，温度：23±3℃，保持无菌通风环境，每隔5天对培养皿进行摇动以防藻细胞贴壁，并取0.5～1mL藻液进行镜检，培养周期为28±7天，配制下述培养基溶液对待测藻株进行室内扩大培养，培养基配方如下：
30g/L NaCl，1.5g/L NaNO3，1.4g/L K2HPO4，1.75g/L MgSO4·7H2O，1.36g/LCaCl2·7H2O，1.2g/LNa2CO3，0.006g/L FeC6H5O7，0.005g/LNaH2PO4·2H2O，0.5g/LCo(NO3)2·6H2O，0.8g/LCuSO4·5H2O，2.3g/LZnSO4·7H2O，0.03g/LH3BO3，4.0g/LNa2MoO4·2H2O，0.02g/LMnCl2·4H2O，0.5g/LVB1，0.5g/LVB12，VH 0.5g/L，超纯水定容至1L；
步骤2，利用本发明改良的CTAB法提取该杜氏藻D.quartolecta的全基因组DNA，保证DNA浓度不低于150ng/μL，OD260/OD280在1.8至1.9之间，无蛋白、盐离子和RNA污染；具体程序是：取600～800mg室内扩大培养的藻细胞，4℃8000r/min离心1.5min，加液氮研磨15sec，加入800μL 20℃预热的2％W/V的CTAB溶液，1μL 1％的β-巯基乙醇(V/V)，混匀后在60℃水浴1.5h，期间每隔20min摇匀1次，加入800μLTris饱和酚，混匀后4℃12000r/min离心2.5min后取上清，加入体积比为25：24：2的Tris饱和酚、氯仿和异戊醇混合液，涡旋振荡后4℃静置10min，期间混匀2～3次，加入800μL 0.1％DEPC(V/V)处理的ddH2O，60℃水浴30min，4℃12000r/min离心4min后取上清液，加入150mL3mol/L的乙酸钠及250mL 4-5℃的无水乙醇，-20℃沉淀50min，4℃10000r/min离心3min后弃上清液，加入1mL4～5℃预冷的70％(V/V)乙醇溶液并涡旋振荡20sec，弃上清后在核酸真空干燥***中挥发液体，加入适量的100×TE缓冲液(10mmol/LTris-HCl，1mmol/L EDTA)溶解沉淀；
步骤3，全基因组DNA用非接触式超声破碎仪强档能量(80～100W)打断5次后(6sec/次，On/6s Off，每隔3sec一次)，获得符合长度要求(300～400bp)的短DNA片段；
步骤4，DNA片段经1.5％TBE琼脂糖凝胶回收及磁珠纯化选择(本发明选用AgencourtAMPure XP Beads磁珠)，近一步筛选获得大小在300～400bp的样品，检测样品质量，确保基因组DNA质量符合步骤(1)的质量标准；
步骤5，获得的合格DNA样品在T4 DNA聚合酶和Klenow聚合酶作用下修复末端，制备平末端，并在3’端加A碱基；配制连接反应体系：1μLT4 DNA连接酶，1μLT载体，5μL 1×连接反应缓冲液，5μL接头(10μmol/L)，5μL DNA样品，无菌水定容至20μL；16℃水浴过夜后获得连接反应产物，按照Agencourt AMPure XP试剂盒要求纯化产物；纯化的产物经感受态细胞转化、蓝白筛选后进行菌液PCR验证及测序(此步骤可委托测序公司完成)，选取阳性克隆结果，利用Agilent 2100Bioanalyzer仪检测扩增产物；将阳性扩增产物96℃变性30sec后配制DNA环化扩增体系：2μL DNA样品，4μL 5×Rapid ligation buffer，1μL连接酶，双蒸水定容至20μL；上述扩增体系经25℃水浴15min后，加入线性DNA消化酶消化10min，最后获得DNA测序文库；利用Agilent SureSelectQXT WGS仪检测文库浓度，保证文库的浓度不超过2nmol/L，体积不低于12μL；
步骤6，将步骤5获得的测序文库进行梯度PCR，配制扩增体系：1μL待测文库样本，引物对各1μL(可选用二代测序接头引物试剂盒)，0.5μL DNA聚合酶，2.5μL dNTPs，1.5μLMgCl2，2.5μLbuffer缓冲液，ddH2O定容至25μL；PCR扩增程序为：96℃3min，96℃30sec，循环40次(每0.5sec降低1℃直至56℃，72℃45sec)，72℃8min，4℃保存；扩增后的片段通过联合锚定聚合技术(cPAS)进行高通量测序，该步骤须委托有相关技术资质的测序公司完成；
步骤7，将步骤6获得的杜氏藻D.quartolecta原始测序数据进行过滤，利用ngsQCToolkit 2.3.3过滤掉低质量(长度小于5kb的短序列，平均质量低于8的序列及接头序列)的测序数据，获得的高质量测序数据分别以FASTQ文件格式存储，文件命名为Dq.fq，对杜氏藻D.quartolecta全基因组测序数据(Dq.fq)进行核心片段筛选与组装。
步骤8，核心基因组片段筛选与组装的具体步骤为：从杜氏藻D.quartolecta测序数据(Dq.fq)中筛选测序深度在50～80×，平均长度在12～15K，N50长度大于18K的测序数据集，将其回帖到盐生杜氏藻(D.salinalina)参考基因组(Dunsal1 v.2)上，同时用Picard软件对回帖结果进行质控，设置比对率≥90％，比对参数为1e-10，筛选符合上述条件的序列作为杜氏藻D.quartolecta基因组核心序列候选数据；利用Burrows-WheelerAlignment(BWA)软件将剩余的杜氏藻D.quartolecta基因组测序数据与核心序列候选数据进行BLASTn比对，比对参数设置为1e-8，Falcon软件进行纠错，获取比对数据间的重叠区后，利用SOAPde novo 2.04软件进行重叠群组装，设置的程序命令为：
1)#maximal read length
2)max_rd_len＝100
3)[LIB]
4)#average insert size
5)avg_ins＝300
6)#ifsequence needs to be reversed
7)reverse_seq＝0
8)#in which part(s)the reads are used
9)asm_flags＝3
10)#use only first 100 bps ofeach read
11)rd_len_cutoff＝100
12)#in which order the reads are used while scaffolding
13)rank＝1
14)#cutoffofpair number for a reliable connection(at least 3 forshort insert size)
15)pair_num_cutoff＝3
16)#minimum aligned length to contigs for a reliable read location(atleast 32for short insert size)
17)map_len＝32
18)#a pair offastq file,read 1 file should always be followed by read2 file
19)q1＝/vol3/agis/xiaoyutao_group/wangjingchun/yanzao/Dq_1.fq
20)#SOAPdenovo-63mer all–s config.txt-p 10-K 55-M 3-F-u–o
21)#SOAPdenovo-63mer all-s-config.txt p 40-K 27-D 1-N 500m-o./result/MDCZ_27>MDCZ_27.log
22)SOAPdenovo-63mer all-s/vol3/agis/xiaoyutao_group/wangjingchun/yanzao/soapdenovo/config.txt-p 10-K 55-o
23)/vol3/agis/xiaoyutao_group/wangjingchun/yanzao/soapdenovo/test
24)qsub-l nodes＝1-q queue8./soap.sh
步骤9，用ABySS 2.2.3软件对上述各重叠群进行再次组装，设置的程序命令为：
25)conda install-c conda-forge-c bioconda-c defaults ABySS
26)ABYSS-k 31-o/vol3/agis/xiaoyutao_group/wangjingchun/yanzao/ABySS/31_contigs.fa
27)/vol3/agis/xiaoyutao_group/wangjingchun/yanzao/Dq.fq
28)qsub-l nodes＝1-q queue6./ABySS.sh
步骤10，利用BUSCO 2.0软件对上述杜氏藻D.quartolecta基因组组装序列的质量进行评估，选择完整基因占比≥20％，单拷贝基因占比15％，多拷贝基因占比≥12％，缺失/空位占比≤3％的组装序列作为杜氏藻D.quartolecta核心基因组序列。设置的程序命令为：
29)python/public/home/wangjingchun/miniconda2/bin/run_BUSCO.py-i
30)/vol3/agis/xiaoyutao_group/wangjingchun/yanzao/02busco/Dq_contig.fa-m geno-l
31)/vol3/agis/xiaoyutao_group/wangjingchun/yanzao/02busco/eukaryota_odb10-o results_Dq
步骤11，利用Augustus 3.3.3软件对筛选的核心基因组组装数据进行功能基因CDS预测，设置的程序命令为：
32)augustus--strand＝both--genemodel＝partial--singlestrand＝false--protein＝on--introns＝on--start＝on--stop＝on--cds＝on--codingseq＝on--alternatives-from-evidence＝true--gff3＝on--UTR＝false--outfile＝/vol3/agis/xiaoyutao_group/wangjingchun/yanzao/04gene/Dqaugustus/out.gff--species＝volvox/vol3/agis/xiaoyutao_group/wangjingchun/yanzao/04gene/Dq/Dq_masked.fa
步骤12，利用Circos软件构建该藻的核心基因组环状图谱，设置的程序命令为：
33)#circos.conf
34)karyotype＝data/karyotype/karyotype.Dq.txt
35)<ideogram>
36)<spacing>
37)default＝0.005r
38)</spacing>
39)radius＝0.9r
40)thickness＝20p
41)fill＝yes
42)</ideogram>
43)#The remaining content is standard and required.It is imported
44)#from default files inthe Circos distribution.
45)#These shouldbe present in every Circos configuration file and
46)#overridden as required.To see the content ofthese files,
47)#look in etc/in the Circos distribution.
48)<image>
49)#Included from Circos distribution.
50)<<include etc/image.conf>>
51)</image>
52)#RGB/HSV color definitions,colorlists,location offonts,fillpatterns.
53)#Included from Circos distribution.
54)<<include etc/colors_fonts_patterns.conf>>
55)#Debugging,I/O an dother systemparameters
56)#Included from Circos distribution.
57)<<include etc/housekeeping.conf>>
综合上述基因组组装质量评估结果可知，杜氏藻D.quartolecta中可筛选到一段核心基因组序列，该核心基因组序列大小为6592916bp，重叠群数量为3000个，最大重叠群长度为1133322bp，重叠群平均长度为2197.64bp，重叠群N50为15270，完整基因占比23.65％，单拷贝基因占比15.18％，多拷贝基因占比13.76％，空位/缺失占比1.89％，预测的CDS占比38.03％，核心基因组环状图谱如图1所示。
实施例2
一种利用杜氏藻D.quartolecta核心基因组序列进行品系鉴定的方法，包括以下步骤：
步骤1，样本搜集、纯化与培养：采集待测藻株(暂命名为Dunaliella sp.)，将待测藻株经纯化后进行室内扩大培养，具体步骤为：对待测藻株藻细胞进行无菌条件下的单克隆挑取，显微镜检合格后在无菌条件下进行室内扩大培养，室内扩大培养条件为：光周期为18h：6h，光照强度为19000lx，温度：23±3℃，保持无菌通风环境，每隔5天对培养皿进行摇动以防藻细胞贴壁，并取0.5～1mL藻液进行镜检，培养周期为28±7天，配制下述培养基溶液对待测藻株进行室内扩大培养，培养基配方如下：
30g/L NaCl，1.5g/L NaNO3，1.4g/L K2HPO4，1.75g/L MgSO4·7H2O，1.36g/LCaCl2·7H2O，1.2g/LNa2CO3，0.006g/L FeC6H5O7，0.005g/LNaH2PO4·2H2O，0.5g/LCo(NO3)2·6H2O，0.8g/LCuSO4·5H2O，2.3g/LZnSO4·7H2O，0.03g/LH3BO3，4.0g/LNa2MoO4·2H2O，0.02g/LMnCl2·4H2O，0.5g/LVB1，0.5g/LVB12，VH 0.5g/L，超纯水定容至1L；扩大培养的藻株分为4份样本(编号1～4)。
步骤2，全基因组DNA提取：分别取成熟期(约30天左右)的藻液(图2)，4℃低温离心1.5min(8000r/min)，富集藻细胞后用液氮速冻并快速研磨15sec后，利用改良的CTAB法分别提取全基因组DNA，具体程序是：向研磨粉末中加入800μL 20℃预热的2％(W/V)的CTAB溶液，另加入1μL1％的β-巯基乙醇(V/V)，轻轻混匀后60℃水浴1.5h，加入800μL Tris饱和酚，轻轻混匀，4℃12000r/min离心2.5min后取上清，加入体积比为25：24：2的Tris饱和酚、氯仿和异戊醇混合液，涡旋振荡后4℃静置10min，期间轻轻混匀2～3次，加入800μL 0.1％DEPC(V/V)处理的ddH2O，60℃水浴30min，4℃12000r/min离心4min后取上清液，加入150mL3mol/L的乙酸钠及250mL4～5℃预冷的无水乙醇，-20℃沉淀50min，4℃10000r/min离心3min后弃上清液，加入1mL4～5℃预冷的70％(V/V)的乙醇溶液，涡旋振荡20sec，弃上清后在核酸真空干燥***内挥发液体，加入100μL 100×TE缓冲液(10mmol/L Tris-HCl，1mmol/L EDTA)溶解沉淀，1％(W/V)的琼脂糖凝胶电泳联合荧光定量仪检测基因组DNA质量，保证DNA浓度不低于150ng/μL，OD260/OD280在1.8至1.9之间，无蛋白、盐离子和RNA污染。琼脂糖凝胶电泳检测结果显示(图3)，1号和4号样本的DNA浓度较高，完整性较好；荧光定量检测结果同样显示(表1)，1号和4号样本的DNA浓度较高，污染较少，适合作为下一步文库构建的候选样本。
表1荧光定量检测待鉴定藻类样本全基因组DNA质量
样品编号 稀释倍数(×) 上样量(μL) 检测浓度(ng/μL) OD<sub>260</sub>/OD<sub>280</sub>
1 1 1 204.6 1.85
2 1 1 152.0 1.69
3 1 1 72.2 1.62
4 1 1 384.1 1.89
步骤3，构建DNA测序文库：取约2.0μg全基因组DNA，80～100W非接触式超声破碎仪强档能量打断5次后(6sec/次，On/6s Off，每隔3sec一次)，获得符合长度要求(300～400bp)的短DNA片段，然后进行琼脂糖凝胶电泳(琼脂糖凝胶浓度：1％，电压：150V)，电泳30min后进行EB染色，紫外灯下截取大约300～400bp的片段回收，在溶解的胶回收液中加入10μL吸附范围在300～400bp的硅基磁珠(本发明选用AgencourtAMPure XP Beads磁珠)混匀，置于磁力架进行分离，将分离磁珠用150μL 80％的乙醇清洗2-3次，加入15μL 0.1×TE，混合后在室温静置10min，将离心管置于磁力架上，约8min后收集上清液。荧光定量检测合格后，获得的合格DNA样品在T4 DNA聚合酶和Klenow聚合酶作用下修复末端，制备平末端，3’末端加A；配制连接反应体系：1μL T4 DNA连接酶，1μLT载体，5μL 1×连接反应缓冲液，5μL接头(10μmol/L)，5μL DNA样品，无菌水定容至20μL。16℃水浴过夜后获得连接反应产物，按照Agencourt AMPure XP试剂盒要求纯化产物；纯化的产物经转化、筛选后进行菌液PCR验证及测序(此步骤可委托测序公司完成)，选取阳性克隆结果，利用Agilent2100Bioanalyzer仪检测扩增产物；将扩增产物96℃变性30sec后至于冰上，配制DNA环化扩增体系：2μL DNA样品，4μL 5×Rapid ligationbuffer，1μL连接酶，双蒸水定容至20μL。上述扩增体系经25℃水浴15min后，加入线性DNA消化酶室温消化10min，最后获得DNA测序文库，利用Agilent SureSelectQXT WGS仪检测文库浓度，保证单个文库的浓度不超过2nmol/L，体积不低于12μL。
步骤4，将所述步骤3获得的测序文库进行梯度PCR，配制扩增体系：1μL待测文库样本，引物对各1μL(可选用二代测序接头引物试剂盒)，0.5μL DNA聚合酶，2.5μL dNTPs，1.5μLMgCl2，2.5μLbuffer缓冲液，ddH2O定容至25μL；PCR扩增程序为：96℃3min，96℃30sec，循环40次(每0.5sec降低1℃直至56℃，72℃45sec)，72℃8min，4℃保存；扩增后的片段通过联合锚定聚合技术(cPAS)进行高通量测序，获取待测藻株全基因组测序数据(该步骤可委托有相关技术资质的测序公司完成)。
步骤5，将所述步骤4获得的待测藻株的原始测序数据经质控(Q20>96％，GC含量>45％)后，分别进行数据过滤，利用ngsQCToolkit 2.3.3软件过滤掉低质量(长度小于5kb的短序列，平均质量低于8的序列及接头序列)的测序数据，过滤参数设置为“-l 20-q 0.5-n0.03-A 0.28”，获得的高质量测序数据(表2)以FASTQ文件格式存储，文件命名为Dsp.fq。
表2过滤后的待鉴定藻株测序信息统计表
样本编号 过滤后片段数 过滤后碱基数 读长 Q20(％) GC(％)
1 238,959 23,895,898 100 97.90 49.11
4 155,286 15,528,625 100 95.36 47.47
由表2可知，经质控检测，编号为1的待鉴定藻株样本测序质量更优(Q20和GC含量较高)，可用于下一步的数据比对与分析。
步骤6，将步骤5获得的待鉴定藻株基因组测序数据(Dsp.fq)及NCBI数据库已公布的5种代表藻类即布氏轮藻(Chara braunii)，衣藻(Chlamydomonas eustigma)，铜绿微囊藻(Microcystis aeruginosa)，片状微囊藻(Microcystis panniformis)以及团藻(Volvoxcarteri)的基因组测序数据进行搜集，以实施例1组装构建好的杜氏藻D.quartolecta核心基因组序列为参考，利用LASTZ1.02.00软件将上述几种藻类的基因组数据与杜氏藻D.quartolecta核心基因组数据进行比对，根据比对的共线性块结果(图4中的A和B)，从中提取每个物种与杜氏藻D.quartolecta对应的基因型，将基因型信息合并、提取、过滤(过滤缺失率≤20％)。
步骤7，以杜氏藻D.quartolecta核心基因组序列为参考，利用BWA0.7.17软件检测步骤5所述物种间的单核苷酸多态性(SNP)、***/缺失位点(Indel)，其中对待测藻株数据进行检测的程序命令为：
1)bwa index-abwtsw Dq.fna#建库
2)bwa aln-t 2-f/vol3/agis/xiaoyutao_group/wangjingchun/yanzao/08snp/Dq_results/Dsp_R1.sai/vol3/agis/xiaoyutao_group/wangjingchun/yanzao/08snp/Dq.fna/vol3/agis/xiaoyutao_group/wangjingchun/yanzao/00data/Dsp_1.fq
3)bwa aln-t 2-f/vol3/agis/xiaoyutao_group/wangjingchun/yanzao/08snp/Dsp_results/Dsp_R2.sai/vol3/agis/xiaoyutao_group/wangjingchun/yanzao/08snp/Dq.fna/vol3/agis/xiaoyutao_group/wangjingchun/yanzao/00data/Dsp_2.fq
4)bwa sampe-f/vol3/agis/xiaoyutao_group/wangjingchun/yanzao/08snp/Dsp_results/Dsp.sam/vol3/agis/xiaoyutao_group/wangjingchun/yanzao/08snp/Dq.fna/vol3/agis/xiaoyutao_group/wangjingchun/yanzao/08snp/Dsp_results/Dsp_R1.sai/vol3/agis/xiaoyutao_group/wangjingchun/yanzao/08snp/Dsp_results/Dsp_R2.sai/vol3/agis/xiaoyutao_group/wangjingchun/yanzao/00data/Dsp_1.fq/vol3/agis/xiaoyutao_group/wangjingchun/yanzao/00data/Dsp_2.fq
5)samtools view-@20-b-S/vol3/agis/xiaoyutao_group/wangjingchun/yanzao/08snp/Dsp_results/Dsp.sam-o/vol3/agis/xiaoyutao_group/wangjingchun/yanzao/08snp/Dsp_results/Dsp.bam
6)samtools sort-@20-m 150G/vol3/agis/xiaoyutao_group/wangjingchun/yanzao/08snp/Dsp_results/Dsp.bam-o/vol3/agis/xiaoyutao_group/wangjingchun/yanzao/08snp/Dsp_results/Dsp.sort.bam
7)samtools rmdup-S/vol3/agis/xiaoyutao_group/wangjingchun/yanzao/08snp/Dsp_results/Dsp.sort.bam/vol3/agis/xiaoyutao_group/wangjingchun/yanzao/08snp/Dsp_results/Dsp.rmdup.bam
8)samtools index/vol3/agis/xiaoyutao_group/wangjingchun/yanzao/08snp/Dsp_results/Dsp.rmdup.bam/vol3/agis/xiaoyutao_group/wangjingchun/yanzao/08snp/Dsp_results/Dsp.rmdup.bam.bai
9)samtools mpileup-gf/vol3/agis/xiaoyutao_group/wangjingchun/yanzao/08snp/Dq.fna/vol3/agis/xiaoyutao_group/wangjingchun/yanzao/08snp/Dsp_results/Dsp.rmdup.bam>/vol3/agis/xiaoyutao_group/wangjingchun/yanzao/08snp/Dsp_results/Dsp.bcf
10)bcftools view-A./Dsp.bcf>Dsp.vcf
步骤8，测序片段读长大于100的使用aln mem算法，程序命令如下：
1)samtools view-@20-b-S./result/SRR2602391.sam-o./result/Dsp.fq.bam
2)samtools sort-@20-m150G./result/Dsp.fq.bam-o./result/Dsp.fq.sort.bam
3)samtools rmdup-S./result/Dsp.fq.sort.bam./result/Dsp.fq.rmdup.bam
4)samtools index./result/SRR2602391.rmdup.bam./result/Dsp.fq.rmdup.bam.bai
5)samtools mpileup-gf./database/grape.fa./result/*.rmdup.bam>Vitis_2.bcf
6)bcftools call-Avm Vitis.bcf>Vitis.vcf
步骤9，其它代表性藻类基因组的SNP、InDel检测程序及算法同上述待鉴定藻株步骤。
步骤10，利用BWA0.7.17软件检测其中的有效单核苷酸多态性(SNP)和***/缺失位点(InDel)数据，检测SNP和InDel时，须先标记出重复片段并忽略掉，再对InDel附近区域进行重新比对，最后筛选获得SNP和InDel。由待鉴定藻株Dq_SX的统计结果可知(表3)，该藻株基因组的主要SNP类型以核苷酸转换为主，颠换主要发生在腺嘌呤(A)与胸腺嘧啶(T)之间。
表3待鉴定藻株的SNP和InDel统计
物种 待鉴定藻株Dunaliellasp.
SNP数量 968,450
InDel数量 61,140
SNP类型1 TC转换(数量：167,620)
SNP类型2 AG转换(数量：167,120)
SNP类型3 GA转换(数量：167,060)
SNP类型4 CT转换(数量：266,320)
SNP类型5 AT颠换(数量：200,330)
步骤11，利用EasySpeciesTree 1.0软件，基于上述获得的有效SNP数据进行***进化树构建(图5)，进而确定待测藻株与杜氏藻D.quartolecta间的亲缘关系，采用最大似然算法，步长值为1000，程序命令设置如下：
1)orthofinder-forthsp1-M msa-S diamond-t 16-a 16
2)orthofinder-f/vol3/agis/xiaoyutao_group/wangjingchun/yanzao/06tree-M msa-S diamond-t 10-a 10-o
3)vol3/agis/xiaoyutao_group/wangjingchun/yanzao/06tree/results
4)输入文件-in1 OrthoFinder/Results_Sep25/WorkingDirectory/SpeciesIDs.txt的第二列-in4 cat*.fa>>all.pepe.fa
5)输入文件-in2和-in3来自xxx/orthsp1/OrthoFinder/Results_Sep25/Orthogroups
6)-in2#cp Orthogroups_SingleCopyOrthologues.txt../../../easy/SingleCopyOrthologues.txt
7)-in3#cp Orthogroups.tsv../../../easy/Orthogroups.csv
8)python2.7/vol1/agis/xiaoyutao_group/wangjingchun/software/EasySpeciesTree/EasySpeciesTree.py-in1
9)SpeciesIDs.txt-in2 SingleCopyOrthologues.txt-in3 Orthogroups.csv-in4 all.pep.fa-t 2
步骤12，基于构建好的***进化树各分支间的支持率和遗传相似度百分比值，确定待测藻株是否属于杜氏藻D.quartolecta，即当待测藻株与杜氏藻D.quartolecta间的支持率在0.99-1.00，相似度百分比≥99％时，基因组覆盖度≥55％，可确定该藻为杜氏藻D.quartolecta。由图4可知，待鉴定藻株(Dunaliella sp.)与杜氏藻D.quartolecta间的支持率为1.00，相似度百分数为100％，基因组覆盖度达56.8％，可鉴定该藻株为杜氏藻D.quartolecta。
实施例3
以杜氏藻D.quartolecta核心基因组数据为参考，分析一株已鉴定藻株Dq_SX基因组的遗传变异与进化特征，包括以下步骤：
步骤1，参照本发明构建的杜氏藻D.quartolecta核心基因组测序数据筛选与组装的方法(见实施例1)，利用SOAPde novo 2.04软件对一株已鉴定杜氏藻株(暂命名为Dq_SX)的全基因组测序数据进行核心片段筛选与从头组装(方法可参见实施例1和2)，筛选并组装的该藻株核心基因组序列主要指标见表4。
步骤2，利用LASTZ 1.02.00软件对步骤1中已构建的杜氏藻Dq_SX核心基因组组装数据进行共线性分析，获得该物种基因组内不同区域间发生加倍事件的重复片段(图6)。
步骤3，以本发明构建的杜氏藻D.quartolecta核心基因组序列为参考模板，利用TBtools软件将其与步骤1组装好的Dq_SX核心基因组数据进行比对，利用Orthofinder2.3.11软件从比对结果中筛选两者间的同源基因，筛选条件设置为：p-value<10-50，score>80，程序命名设置为：conda install-c conda-forge-c bioconda-c defaultsorthofinder/#/public/home/wangjingchun/miniconda2/bin/orthofinder/condainstall-c conda-forge-c bioconda-c defaults trimal#/public/home/wangjingchun/miniconda2/bin/trimal/conda install-c conda-forge-c bioconda-c defaultsraxmlHPC/#/public/home/wangjingchun/miniconda2/bin/raxmlHPC-PTHREADS/###/public/home/wangjingchun/miniconda2/bin/mafft/###/vol1/agis/xiaoyutao_group/wangjingchun/software/ASTRAL-master/Astral/astral.5.7.3.jar。
步骤4，以步骤3筛选的同源基因信息作为数据分析集，利用PAML 4.8软件检测其中的同义和非同义突变位点，并计算非同义替换率(Ka)和同义替换率(Ks)值，由Ka/Ks值推测已鉴定藻株Dq_SX的进化选择压力(图7)。
表4一株已鉴定杜氏藻Dq_SX核心基因组组装数据及其质量评估
Figure BDA0002767618840000241
由表4可知，已鉴定杜氏藻株Dq_SX核心基因组组装结果较完整，不完整片段仅占16.12％，仅有1.54％的空位或缺失。由图6可知，已鉴定藻株Dq_SX在进化过程中，其基因组内的不同区域可能发生了大量的加倍事件，涉及到加倍事件的区段对有1007对，这暗示了该物种进化过程的复杂性。由图7可知，相对于本发明已构建的杜氏藻D.quartolecta核心基因组，已鉴定藻株Dq_SX核心基因组中80.52％的基因Ka/Ks比值小于1.0(Ka/Ks平均值为0.47；当Ka/Ks比值在0.35-0.45范围内时，频率最高为0.108)，暗示了该藻株的多数基因在进化过程中受到了纯化选择的压力(图7)。
实施例4
已鉴定杜氏藻株Dq_SX核心基因组的重复片段预测、预测蛋白的功能注释及其结构特征分析，包括以下步骤：
步骤1，利用Repeatmasker 4.0.9软件对实施例3中已鉴定的一株杜氏藻Dq_SX核心基因组组装数据进行重复序列分析，先构建好待测序列数据库(BuildDatabase-nameDq_SX Dq_SX_contig.fa)，设置以下程序命令：
1)RepeatModeler-pa 10-database/vol3/agis/xiaoyutao_group/wangjingchun/yanzao/03repeat/Dq_SX/Dq_SX-engine ncbi-recoverDir/vol3/agis/xiaoyutao_group/wangjingchun/yanzao/03repeat/Dq_SX
2)qsub-l nodes＝1-q queue8./repeatmodeler.sh
步骤2，在家目录下得到consensi.fa consensi.fa.masked families.stk
步骤3，#fasta文件“*-families.fa”为训练得到的共识重复序列，序列id后会标注它属于哪种重复序列家族，若无法归类则用标注为“Unkown”。“*-families.stk”为种子联配(Seed alignments)文件，是Dfam兼容的Stockholm格式，可以使用RepeatModeler安装路径中自带的工具“RepeatModeler/util/dfamConsensusTool.pl”上传到Dfam_consensus数据库中。
步骤4，再搜索杜氏藻Dq_SX核心基因组中的重复序列，设置程序命令为：
1)RepeatMasker-pa 4gff lib/public/home/wangjingchun/RM_Dq_SX/consensi.fa dir/vol3/agis/xiaoyutao_group/wangjingchun/yanzao/03repeat/Dq_SX/Repeatmasker/lib_result/vol3/agis/xiaoyutao_group/wangjingchun/yanzao/03repeat/Dq_SX/Dq_SX_contig.fa
2)qsub-l nodes＝1-q queue8./repeatmasker2.sh
步骤5，利用Diamonds 0.9.14软件将已鉴定藻株核心基因组的CDS编码蛋白序列与非冗余蛋白数据库(NR)进行BLASTp比对，进而获得蛋白的功能注释，设置比对参数为1e-value≤10-5，设置程序命令如下：
1)$diamondmakedb--innr_eukaryon.fasta-d nr_eukaryon_20200805
2)$diamond blastx--db nr_eukaryon_20200805--query reads.fq.gz--outreads.tab
3)$diamond blastp--db nr_eukaryon_20200805--query proteins.fasta--outnr.tab--outfmt 6--sensitive--max-target-seqs 20--evalue 1e-5--id 30--block-size20.0--tmpdir/dev/shm--index-chunks 1
步骤6，利用MCScanX软件进行已鉴定藻株核心基因组的重复片段共线性分析，设置程序命令为：
1)makeblastdb-in Dq_SX.fa-dbtype prot-out Dq_SX
2)Blastp-query/vol3/agis/xiaoyutao_group/wangjingchun/yanzao/07circos/Dq_SX.fa-db/vol3/agis/xiaoyutao_group/wangjingchun/yanzao/07circos/Dq_SXnum_threads 10-evalue 1e outfmt 6out/vol3/agis/xiaoyutao_group/wangjingchun/yanzao/07circos/Dq_SX.blastp
3)MCScanX./Dq_SX
步骤7，由于不同物种间重复序列的保守性相对较低，针对特定物种进行重复序列的预测需要查询特定的重复序列数据库。鉴于此，我们利用RepeatMaskerv4.0.6软件将已鉴定藻株Dq_SX核心基因组测序组装数据与RepBase中的数据进行比对，查询该藻株中可能的散在重复序列。利用RepeatModeler、LTR-FINDER、RepeatScout软件对已鉴定藻株Dq_SX核心基因组数据进行注释，获得串联重复序列(包括：微卫星序列，小卫星序列等)。
步骤8，将上述结果中的重复部分过滤，获得最终非冗余的重复序列注释结果(表5)。
步骤9，将已鉴定藻株Dq_SX核心基因组数据与NR数据库进行比对，并对比对结果进行筛选(e-value<10-5)。
步骤10，利用eggNOG软件对筛选后的同源蛋白序列进行COG功能注释，使用eggNOG中的emapper.py脚本对蛋白序列进行注释，对注释结果中的前20位(top20)蛋白簇进行分类统计(图8)。
步骤11，运行eggNOG软件对其基因编码的同源蛋白进行COG功能注释；设置程序命令如下：
python/public/home/wangjingchun/miniconda2/envs/qiime1/bin/emapper.py-i/vol3/agis/xiaoyutao_group/wangjingchun/yanzao/new/04cog/Dq_SX_protein.fa--output/vol3/agis/xiaoyutao_group/wangjingchun/yanzao/new/04cog/out-mdiamond--data_dir/vol3/agis/xiaoyutao_group/wangjingchun/yanzao/new/04cog/database--cpu20
步骤12，利用TMHMM2.0在线分析工具对上述top20蛋白簇中排名首位的蛋白进行跨膜结构域预测分析(图9)；利用SignalP4.1在线分析工具对上述蛋白进行信号肽预测，设置每条蛋白序列的氨基酸数量阈值不超过6000个(图10)；输出格式选为Extensive,withgraphic，其它参数选择默认。
表5一株已鉴定杜氏藻Dq_SX核心基因组中重复序列分类结果统计
重复序列类型 重复序列大小(bp) 重复序列的基因组占比(％)
LINE 165380 0.26
LTR 118737 0.19
SINE 984126 1.57
其它 1007445 1.60
总数 2275688 3.62
由表5可知，一株已鉴定杜氏藻Dq_SX核心基因组中共搜索到重复序列长度2275688bp，占全基因组比例约3.62％。由图6可知，Dq_SX核心基因组已注释功能蛋白中转录调控因子(88条)和动力蛋白重链(87条)分类数量最多。其中，转录调控因子跨膜结构域预测结果显示(图9)，在该因子60-110位氨基酸的结构很可能处于膜外(概率值约0.8)，130位氨基酸之后的部分大概率(概率值0.82)处于膜内，位于膜上的概率不高于0.4。由该因子的信号肽预测结果可知(图10)，在25-26位氨基酸左右，C值最大，S值陡峭，Y值达最高峰，暗示该处为信号肽剪切位点。
实施例5
基于杜氏藻D.quartolecta与一株已鉴定藻株Dq_SX核心基因组数据的差异代谢通路比较分析及特色基因挖掘，包括以下步骤：
步骤1，将杜氏藻D.quartolecta及实施例3中一株已鉴定藻株Dq_SX核心基因组中预测的蛋白序列(Dq_SX核心基因组测序组装数据获取方法见实施例3，蛋白序列获取方法见实施例4)与KEGG数据库(日本京都基因和基因组百科全书)进行BLASTp比对，获取基因编码产物可能参与的代谢通路，设置的程序命令如下：
1)diamond makedb--in/vol3/agis/xiaoyutao_group/wangjingchun/yanzao/09kegg/ko.pep.fasta-d
2)/vol3/agis/xiaoyutao_group/wangjingchun/yanzao/09kegg/kegg
3)diamond blastp-d/vol3/agis/xiaoyutao_group/wangjingchun/yanzao/09kegg/kegg--query
4)/vol3/agis/xiaoyutao_group/wangjingchun/yanzao/09kegg/Dq_protein.fa-f6-o
5)/vol3/agis/xiaoyutao_group/wangjingchun/yanzao/09kegg/Dq.blastp-p30-e0.00005
6)diamond blastp-d/vol3/agis/xiaoyutao_group/wangjingchun/yanzao/09kegg/kegg--query
7)/vol3/agis/xiaoyutao_group/wangjingchun/yanzao/09kegg/Dq_SX_protein.fa-f6-o
8)/vol3/agis/xiaoyutao_group/wangjingchun/yanzao/09kegg/Dq_SX.blastp-p 30-e 0.00005
步骤2，依据步骤1中各代谢通路分配的KO编号，将杜氏藻D.quartolecta及已鉴定藻株Dq_SX的KEGG通路预测结果进行交集分析，构建文氏图(图11)，从中筛选各自独有的代谢通路。
步骤3，根据步骤2获得的杜氏藻D.quartolecta及已鉴定藻株Dq_SX独有代谢通路top20(步骤2统计中，除交集外的排名前20位KEGG通路)，分别从中筛选富集程度最高，即富集比率最高的特色基因(图12和图13)。
步骤4，分别将上述杜氏藻D.quartolecta与已鉴定藻株Dq_SX中富集程度最高(显著富集)的代谢通路基因在GO(基因本体联合)数据库中进行查询分析，进一步获取其排名前20(top20)的GO功能注释富集结果(图14和图15)，从富集程度较高，即GO条目数较多且对应的-log10(Q值)(置信度)也较高的基因集中筛选研究者感兴趣的特色基因。
由图11可知，基于杜氏藻D.quartolecta和Dq_SX核心基因组测序数据，我们共预测到KEGG通路608条，其中共有代谢通路141条，特色代谢通路467条(其中，杜氏藻D.quartolecta 85条，Dq_SX 382条)。由图12和图13可知，杜氏藻D.quartolecta特有代谢通路中富集程度最高的是剪接体相关代谢，Dq_SX中富集程度最高的是细胞组分合成相关代谢。由图14和图15可知，杜氏藻D.quartolecta剪接体代谢途径中参与基因的功能多与RNA转运、加工及合成紧密相关，Dq_SX膜组分合成代谢中参与基因的功能多与蛋白质结构与加工过程相关。
实施例6
三种不同的杜氏藻D.quartolecta分子鉴定技术比较分析，
搜集20株待测藻株及一株实施例1中已鉴定的杜氏藻D.quartolecta，利用ITS基因、SSR分子标记及基因组测序数据对待测藻株进行分子鉴定，具体包括以下步骤：
步骤1，利用本发明改良的CTAB法(具体可参见实施例1)提取各藻株基因组DNA，设计并合成杜氏藻ITS基因扩增引物如SEQ ID NO.1和SEQ ID NO.2所示：
SEQ ID NO.1：5'-GAAGGAGAAGTCGTAACAAG-3'；
SEQ ID NO.2：5'-CCTCCCTTATTGATATGC-3'；
配制ITS基因PCR扩增体系：2.0μL dNTPs(2mmol/L)，1.0μL Mg2+(25mmol/L)，1.0μLDNA，0.3μLTaq酶(5U/μL)和2.5μL 10×buffer缓冲液，上述引物各1.0μL，ddH2O补齐至25μL；设置PCR反应程序：95℃3min，95℃30sec，52℃40sec，72℃1min，循环35次后72℃延伸10min；1.2％琼脂糖凝胶电泳检测，收集800～1000bp的特异性扩增产物，送测序公司测序。
步骤2，根据测序公司反馈的测序结果，利用MEGA5.0软件，基于最大似然法构建21株藻的ITS基因***进化树，步长值为1000，根据进化树中各分支节点的支持率及遗传相似百分数，从待测藻株中鉴定杜氏藻D.quartolecta(图16)。
步骤3，基于发明人获得的9组杜氏藻转录组测序数据(NCBI数据库编号：SRR8393723、SRR8393722、SRR8393725、SRR8393724、SRR8393727、SRR8393726、SRR8393729、SRR8393728、SRR8393721)，我们从24311条SSR标记中筛选获得特异性标记15条，依据标记信息，设计10对SSR多态性扩增引物，引物信息如下所示：
CL1007：SEQ ID NO.3：5'-CTAAATCCATGCGTTCTTCTTTC-3'；
SEQ ID NO.4：5'-ACAGTACAACCAGAGGCTTTGAA-3'；
CL1008：SEQ ID NO.5：5'-AACAATGTCACCTCTCATTTGCT-3'；
SEQ ID NO.6：5'-TCGTTTTGTTGTTGTTCTTCAAA-3'；
CL102：SEQ ID NO.7：5'-GCCAATTCCAAAAAGTTAAAATCT-3'；
SEQ ID NO.8：5'-ATTGTGGTTTTCTTCCTGGTTTT-3'；
CL1041：SEQ ID NO.9：5'-AGGCAAGCAGTGCATTTGTA-3'；
SEQ ID NO.10：5'-GGCTCTCTATGAGTCGATGTGTC-3'；
CL1047：SEQ ID NO.11：5'-GCAGTGGAAACACACTTCCTTAC-3'；
SEQ ID NO.12：5'-TCTCTCAAATCAAAGGTGCTTTC-3'；
CL1157：SEQ ID NO.13：5'-GAGATCGAACTTGAGGCTTAGAA-3'；
SEQ ID NO.14：5'-AAAATAGAAGCCATCATGAAACG-3'；
CL1160：SEQ ID NO.15：5'-GGATACAGATTTCCACACTGCTC-3'；
SEQ ID NO.16：5'-CTATCTGGCTGAAGGTCATGTTT-3'；
CL1168：SEQ ID NO.17：5'-CGTTTTTGGAACTGATTTCTTTG-3'；
SEQ ID NO.18：5'-TTCTTGTAATACATCGCAGGAAG-3'；
CL1322：SEQ ID NO.19：5'-AACAGAGGAAATTCTGATGATGC-3'；
SEQ ID NO.20：5'-CTTGCAAGAAGGAACAACTCACT-3'；
CL1627：SEQ ID NO.21：5'-GTGGTCACCAGGAAGAGACAG-3'；
SEQ ID NO.22：5'-ACGGTACTGACAGTGGAAACAAT-3'；
上述扩增产物物的大小依次为155bp、131bp、139bp、121bp、158bp、136bp、118bp、149bp、160bp、127bp；
步骤4，送生物公司合成上述SSR引物，配制SSR-PCR扩增体系，即：2.5μL dNTPs(2mmol/L)，1.2μL Mg2+(25mmol/L)，1.0μL DNA(步骤1获得)，0.4μL Taq酶(5U/μL)和2.5μL10×buffer缓冲液，上述引物各0.8μL，ddH2O补齐至25μL；SSR-PCR反应程序为：94℃5min；35个循环(94℃45sec，57℃35sec，72℃1min)；72℃8min；扩增后的SSR产物用4％的变性聚丙烯酰胺电泳分离，银染30min、显色15min、固定20min后，在电泳图谱上进行“1”(有条带)和“0”(无条带)标记；用UPGMA法，利用NTSYSpc 2.2软件进行待测藻株的聚类分析，构建其SSR标记的***进化树(图17)。
步骤5，以本发明构建的杜氏藻D.quartolecta核心基因组组装数据为参考，基于步骤1获得的全基因组DNA，建立测序文库，建库方法可参照实施例1进行；对待测藻株进行基因组测序，测序片段无须进行从头组装，该步骤可委托有资质的测序公司完成。
步骤6，以本发明构建的杜氏藻D.quartolecta核心基因组数据为参考，利用BWA0.7.17软件检测待测藻株间的单核苷酸多态性(SNP)和***缺失(InDel)数据，检测SNP和InDel时，须先标记出重复片段并忽略掉，再对InDel附近区域进行重新比对，最后筛选获得SNP和InDel，程序命令可参照实施例2进行。
步骤7，利用EasySpeciesTree 1.0软件，基于上述获得的SNP数据进行***进化树构建(图18)，采用最大似然算法，设置步长值为1000，程序命令可参照实施例2进行。
步骤8，比较分析上述三种不同的分子鉴定结果，技术优缺点(表6)。
由图16可知，藻株Dsp11与杜氏藻D.quartolecta聚为一簇，支持率0.99，遗传相似度99％，可鉴定为D.quartolecta。由图17可知，藻株Dsp4与杜氏藻D.quartolecta聚为一簇，Dsp4与杜氏藻D.quartolecta的支持率为0.99，遗传相似度99％，可鉴定为D.quartolecta，Dsp11又与Dsp4和D.quartolecta聚为一簇，支持率为1.00，遗传相似度99％，也可鉴定为D.quartolecta；由图18可知，Dsp11、Dsp4可与D.quartolecta共聚为一簇，支持率为1.00，遗传相似度100％，也可鉴定为D.quartolecta。由表6可知，相较于其它两种分子鉴定方法，以本发明构建的杜氏藻D.quartolecta核心基因组数据为参考，对待测藻株进行简化基因组测序并获得SNP数据，可在短期内(7～10天)准确鉴定D.quartolecta，不仅成本较低，还可为后期深入研究提供丰富的生物信息数据。
表6三种分子鉴定杜氏藻D.quartolecta的技术方法比较
Figure BDA0002767618840000331
实施例7
本发明创建的杜氏藻D.quartolecta核心基因组测序组装技术与传统基因组测序技术的比较，包括以下步骤：
步骤1，对实施例3中一株已鉴定杜氏藻D.quartolecta进行基因组DNA提取，方法可参见实施例1。
步骤2，将质控合格的DNA样本(DNA浓度≥150ng/μL，电泳条带明亮、无降解，OD260/OD280在1.8～1.9)送测序公司进行DNA测序文库构建、测序、核心片段筛选，从头组装，测序分析平台分别选择Nanopore、PacBio和HiSeq(该步骤可委托有相关测序平台的公司进行操作)。
步骤3，将本发明自主构建的杜氏藻D.quartolecta核心基因组测序片段组装数据(详见实施例1的操作步骤)与步骤2获得的各测序平台的组装数据进行关键指标比较。
步骤4，以NCBI已公布的盐生杜氏藻参考基因组(Dunsal1 v.2)为参考，将上述各技术平台获得的杜氏藻D.quartolecta核心基因组测序数据与其进行比对，根据比对结果分析各技术间的差异(表7)。
表7核心基因组测序数据组装过程中的比对结果分析
Figure BDA0002767618840000341
步骤5，利用SOAPsnp软件对步骤4获得的唯一比对的测序片段进行单核苷酸多态性(SNP)检测，检测过程中过滤掉重复片段，对***/缺失(InDel)位点附近区域进行重比对，筛选有效的高质量SNP；测序数据中的短序列与参考基因组进行比对及聚类分析，检测InDel，设置gap长度：1～10个碱基。比较分析四种技术获得的有效SNP和InDel平均数(表8)。
表8 SNP和InDel统计结果比较分析
Figure BDA0002767618840000342
步骤6，利用步骤4获得的测序片段结合实施例4中的重复片段预测方法，计算不同技术平台条件下该藻株的重复序列占总测序片段比例(表9)。
表9重复序列占比结果比较分析
技术 重复序列占总测序片段比例(％)
自主技术 1.45％
Nanopore 15.27％
PacBio 12.99％
HiSeq 3.58％
由表7可知，本发明创建的技术条件下，待测藻株测序片段的基因组覆盖率、比对序列和鉴定占比均高于其它三种测序技术。由表8可知，本发明技术条件下检测的有效SNP和InDel数量均高于其它三种技术，且错误率最低。由表9可知，本发明技术条件下检测到的重复序列占比低于其它三种技术。综上可知，本发明创建的杜氏藻D.quartolecta核心基因组测序片段组装技术整体性能优于Nanopore、PacBio和HiSeq。
以上显示和描述了本发明的基本原理和主要特征和本发明的优点,对于本领域技术人员而言，显然本发明不限于上述示范性实施例的细节，而且在不背离本发明的精神或基本特征的情况下，能够以其他的具体形式实现本发明。因此，无论从哪一点来看，均应将实施例看作是示范性的，而且是非限制性的，本发明的范围由所附权利要求而不是上述说明限定，因此旨在将落在权利要求的等同要件的含义和范围内的所有变化囊括在本发明内。不应将权利要求中的任何附图标记视为限制所涉及的权利要求。
此外，应当理解，虽然本说明书按照实施方式加以描述，但并非每个实施方式仅包含一个独立的技术方案，说明书的这种叙述方式仅仅是为清楚起见，本领域技术人员应当将说明书作为一个整体，各实施例中的技术方案也可以经适当组合，形成本领域技术人员可以理解的其他实施方式。
序列表
<110> 山西大学
<120> 基于一种杜氏藻核心基因组序列进行品系鉴定的方法
<160> 22
<170> SIPOSequenceListing 1.0
<210> 16
<211> 20
<212> DNA
<213> ITS基因上游引物(ITS-F)
<400> 16
gaaggagaag tcgtaacaag 20
<210> 17
<211> 18
<212> DNA
<213> ITS基因下游引物(ITS-R)
<400> 17
cctcccttat tgatatgc 18
<210> 18
<211> 23
<212> DNA
<213> CL1007上游引物(CL1007-F)
<400> 18
ctaaatccat gcgttcttct ttc 23
<210> 19
<211> 23
<212> DNA
<213> CL1007下游引物(CL1007-R)
<400> 19
acagtacaac cagaggcttt gaa 23
<210> 20
<211> 23
<212> DNA
<213> CL1008上游引物(CL1008-F)
<400> 20
aacaatgtca cctctcattt gct 23
<210> 21
<211> 23
<212> DNA
<213> CL1008下游引物(CL1008-R)
<400> 21
tcgttttgtt gttgttcttc aaa 23
<210> 22
<211> 24
<212> DNA
<213> CL102上游引物(CL102-F)
<400> 22
gccaattcca aaaagttaaa atct 24
<210> 23
<211> 23
<212> DNA
<213> CL102下游引物(CL102-R)
<400> 23
attgtggttt tcttcctggt ttt 23
<210> 24
<211> 20
<212> DNA
<213> CL1041上游引物(CL1041-F)
<400> 24
aggcaagcag tgcatttgta 20
<210> 25
<211> 23
<212> DNA
<213> CL1041下游引物(CL1041-R)
<400> 25
ggctctctat gagtcgatgt gtc 23
<210> 26
<211> 23
<212> DNA
<213> CL1047上游引物(CL1047-F)
<400> 26
gcagtggaaa cacacttcct tac 23
<210> 27
<211> 23
<212> DNA
<213> CL1047下游引物(CL1047-R)
<400> 27
tctctcaaat caaaggtgct ttc 23
<210> 28
<211> 23
<212> DNA
<213> CL1157上游引物(CL1157-F)
<400> 28
gagatcgaac ttgaggctta gaa 23
<210> 29
<211> 23
<212> DNA
<213> CL1157下游引物(CL1157-R)
<400> 29
aaaatagaag ccatcatgaa acg 23
<210> 30
<211> 23
<212> DNA
<213> CL1160上游引物(CL1160-F)
<400> 30
ggatacagat ttccacactg ctc 23
<210> 31
<211> 23
<212> DNA
<213> CL1160下游引物(CL1160-R)
<400> 31
ctatctggct gaaggtcatg ttt 23
<210> 32
<211> 23
<212> DNA
<213> CL1168上游引物(CL1168-F)
<400> 32
cgtttttgga actgatttct ttg 23
<210> 33
<211> 23
<212> DNA
<213> CL1168下游引物(CL1168-R)
<400> 33
ttcttgtaat acatcgcagg aag 23
<210> 34
<211> 23
<212> DNA
<213> CL1322上游引物(CL1322-F)
<400> 34
aacagaggaa attctgatga tgc 23
<210> 35
<211> 23
<212> DNA
<213> CL1322下游引物(CL1322-R)
<400> 35
cttgcaagaa ggaacaactc act 23
<210> 36
<211> 21
<212> DNA
<213> CL1627上游引物(CL1627-F)
<400> 36
gtggtcacca ggaagagaca g 21
<210> 37
<211> 23
<212> DNA
<213> CL1627下游引物(CL1627-R)
<400> 37
acggtactga cagtggaaac aat 23

Claims (6)

1.基于一种杜氏藻核心基因组序列进行品系鉴定的方法，其特征在于，包括以下步骤：
(1)样本搜集、纯化与培养：采集待测藻株及杜氏藻Dunaliella quartolecta株，将待测藻株经纯化后进行室内扩大培养，具体步骤为：对待测藻株藻细胞进行无菌条件下的单克隆挑取，显微镜检合格后在无菌条件下进行室内扩大培养，室内扩大培养条件为：光周期为18h：6h，光照强度为19000lx，温度：23±3℃，保持无菌通风环境，每隔5天对培养皿进行摇动以防藻细胞贴壁，并取0.5～1mL藻液进行镜检，培养周期为28±7天，配制下述培养基溶液对待测藻株进行室内扩大培养，培养基配方如下：
30g/L NaCl，1.5g/L NaNO3，1.4g/L K2HPO4，1.75g/L MgSO4·7H2O，1.36g/L CaCl2·7H2O，1.2g/L Na2CO3，0.006g/L FeC6H5O7，0.005g/L NaH2PO4·2H2O，0.5g/L Co(NO3)2·6H2O，0.8g/L CuSO4·5H2O，2.3g/L ZnSO4·7H2O，0.03g/L H3BO3，4.0g/L Na2MoO4·2H2O，0.02g/L MnCl2·4H2O，0.5g/LVB1，0.5g/L VB12，VH 0.5g/L，超纯水定容至1L；
(2)全基因组DNA提取：利用改良的CTAB法分别提取待测藻株及杜氏藻株D.quartolecta的全基因组DNA，冷冻保存；所述改良的CTAB法具体步骤为：取600～800mg待测藻株，用超纯水冲洗2～3次，4℃ 8000r/min离心1.5min，加液氮研磨15sec，加入800μL20℃预热的2％W/V的CTAB溶液，1μL 1％V/V的β-巯基乙醇，混匀后在60℃水浴1.5h，期间每隔20min摇匀1次，加入800μLTris饱和酚，混匀后4℃ 12000r/min离心2.5min后取上清，加入体积比为25：24：2的Tris饱和酚、氯仿和异戊醇混合液，涡旋振荡后4℃静置10min，期间混匀2～3次，加入800μL 0.1％V/V DEPC处理的ddH2O，60℃水浴30min，4℃ 12000r/min离心4min后取上清液，加入150mL 3mol/L的乙酸钠及250mL 4～5℃预冷的无水乙醇，-20℃沉淀50min，4℃ 10000r/min离心3min后弃上清液，加入1mL 4～5℃预冷的70％V/V乙醇溶液并涡旋振荡20sec，弃上清后在核酸真空干燥***中挥发液体，加入100×TE缓冲液溶解沉淀以保证DNA浓度≥150ng/μL，1％W/V的琼脂糖凝胶电泳联合荧光定量仪检测基因组DNA，确保其电泳条带明亮、无降解，OD260/OD280在1.8～1.9，无污染；
(3)将步骤(2)中的待测藻株及杜氏藻D.quartolecta的全基因组DNA打断、纯化后分别构建DNA测序文库；
(4)采用高通量测序法对步骤(3)中的DNA测序文库分别进行测序，获取待测藻株和杜氏藻D.quartolecta全基因组二代测序数据；
(5)以NCBI已公布的盐生杜氏藻(D.salina)全基因组数据为参考，将步骤(4)中获取的杜氏藻D.quartolecta全基因组测序数据与其进行比对，通过筛选、从头组装、质量评估后获得杜氏藻D.quartolecta核心基因组序列，该核心基因组序列大小为6592916bp，重叠群数量为3000个，最大重叠群长度为1133322bp，重叠群平均长度为2197.64bp，重叠群N50为15270，完整基因占比23.65％，单拷贝基因占比15.18％，多拷贝基因占比13.76％，空位/缺失占比1.89％，不完整片段占比17.45％，构建从头组装的杜氏藻D.quartolecta核心基因组环状图谱，然后对杜氏藻D.quartolecta核心基因组序列进行基因组分、蛋白功能注释及基因组重叠群共线性分析；
(6)以步骤(5)中构建的杜氏藻D.quartolecta核心基因组序列为参考，将步骤(4)获得的待测藻株全基因组测序数据及已公布代表性藻类的基因组测序数据与其进行比对，检测物种间的单核苷酸多态性和***/缺失位点，然后利用单核苷酸多态性构建***进化树，当待测藻株与杜氏藻D.quartolecta聚为一簇，且分支的数据支持率在0.99～1.00，遗传相似度百分比大于≥99％，待测藻株即为杜氏藻D.quartolecta。
2.根据权利要求1所述的基于一种杜氏藻核心基因组序列进行品系鉴定的方法，其特征在于，所述步骤(3)中构建DNA测序文库具体步骤为：全基因组DNA用80～100W的强档超声波段打断6sec，每隔3sec重复1次，共超声5次，设置打断参数为300～400bp；片段经琼脂糖凝胶电泳后，胶回收300～400bp目的片段；用硅基磁珠吸附回收目的片段，再进行荧光定量仪检测吸附回收目的片段的质量；DNA末端修复，3’末端加A；加接头进行连接反应，对连接产物进行纯化、转化、PCR验证；阳性产物经95℃变性20sec后，进行单链DNA环化反应，纯化产物后构建得到可用于上机的全基因组DNA测序文库。
3.根据权利要求1所述的基于一种杜氏藻核心基因组序列进行品系鉴定的方法，其特征在于，所述步骤(5)中通过筛选、组装、质量评估后获得杜氏藻D.quartolecta核心基因组序列的具体步骤为：从测序平台筛选获得高质量序列，筛选测序深度在50～80×，平均长度在12～15K，N50长度大于18K的片段作为查询序列，利用SOAPaligner或BWA软件将其回帖到已报道的盐生杜氏藻参考基因组上，进一步筛选序列一致度≥90％，比对结果E值小于1e-10的测序片段作为杜氏藻D.quartolecta基因组核心序列候选数据；所有剩余测序片段与该候选数据集进行比对，获取比对数据间的重叠区；利用Falcon或Pilon软件对比对结果进行纠错及校正运算，利用SOAPde novo 2.04、Mecat、HERA或Canu软件组装重叠群；利用BySS2.2.3、Velvet 1.2.10或ABySS 2.2.3软件确定各重叠群顺序；利用BAMStats或GATKDepthOfCoverage软件进行全基因组覆盖度测算，筛选参考基因组覆盖度≥50％，连续排列数≥2000的重叠群；利用BUSCO 2.0或Quast软件对筛选出的重叠群组装质量进行评估，选择完整基因占比≥20％，单拷贝基因占比15％，多拷贝基因占比≥12％，缺失/空位占比≤3％的组装序列作为杜氏藻D.quartolecta核心基因组序列；利用Circos软件构建该物种的核心基因组环状图谱。
4.根据权利要求1所述的基于一种杜氏藻核心基因组序列进行品系鉴定的方法，其特征在于，所述步骤(5)中对杜氏藻D.quartolecta核心基因组序列进行基因组分、蛋白功能注释及基因组重叠群共线性分析，具体步骤如下：利用Augustus 3.3.3、ESTScan3.0.1、TransDecoder 2.0.1或Prodigal 2.6.1软件对组装数据进行CDS预测，利用Repeatmasker4.0.9、RepeatProteinMask 3.2.2、LTR-FINDER、Piler 1.0.6或RepeatScout 1.0.5软件对组装数据进行重复序列分析，利用Diamonds 0.9.14或BLASTX软件将CDS编码的蛋白序列比对到NR数据库，并进行功能注释，预测的蛋白序列经BLASTp自身比对后，再利用MCScanX、Last、Mugsy、Spines或Progressivemauve软件进行基因组的共线性分析。
5.根据权利要求1所述的基于一种杜氏藻核心基因组序列进行品系鉴定的方法，其特征在于，所述步骤(6)中利用单核苷酸多态性构建***进化树的具体步骤为：利用LASTZ1.02.00或Mauve 2.3.1软件，将待测藻株及5～6种已在NCBI数据库中报道的代表性藻类基因组数据分别与所述步骤(5)中组装的杜氏藻D.quartolecta核心基因组序列进行比对，根据比对的共线性块结果提取每个物种与杜氏藻D.quartolecta基因组的对应基因型，之后以杜氏藻D.quartolecta核心基因组为模板，将所有物种的基因型信息进行合并、提取、过滤，利用BWA 0.7.17软件检测其中的单核苷酸多态性数据和***/缺失位点数据；基于单核苷酸多态性数据，利用EasySpeciesTree 1.0、MEGA 5.0、TreeBeST 1.9.2、PHYLIP、Puzzle5.2或PHYLO-WIN软件中的最大似然算法构建***进化树，进而确定待测藻株与杜氏藻D.quartolecta间的亲缘关系。
6.根据权利要求5所述的基于一种杜氏藻核心基因组序列进行品系鉴定的方法，其特征在于，所述过滤的缺失率不高于20％。
CN202011238521.2A 2020-11-09 2020-11-09 基于一种杜氏藻核心基因组序列进行品系鉴定的方法 Active CN112349350B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202011238521.2A CN112349350B (zh) 2020-11-09 2020-11-09 基于一种杜氏藻核心基因组序列进行品系鉴定的方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202011238521.2A CN112349350B (zh) 2020-11-09 2020-11-09 基于一种杜氏藻核心基因组序列进行品系鉴定的方法

Publications (2)

Publication Number Publication Date
CN112349350A CN112349350A (zh) 2021-02-09
CN112349350B true CN112349350B (zh) 2022-07-19

Family

ID=74428639

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202011238521.2A Active CN112349350B (zh) 2020-11-09 2020-11-09 基于一种杜氏藻核心基因组序列进行品系鉴定的方法

Country Status (1)

Country Link
CN (1) CN112349350B (zh)

Families Citing this family (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN113160893B (zh) * 2021-06-09 2022-08-19 中国科学院昆明植物研究所 从二代测序数据挖掘植物ITSs序列并用于鉴别品种家系
CN113549620B (zh) * 2021-07-13 2022-09-23 山西大学 多型杜氏藻盐胁迫响应miRNAs及其应用
CN114664379A (zh) * 2022-04-12 2022-06-24 桂林电子科技大学 一种基于深度学习的第三代测序数据的自校正纠错方法
CN115810393B (zh) * 2022-12-22 2023-08-25 南京普恩瑞生物科技有限公司 一种基于构建人群SNPs库的测序样本同源性检测方法及***
CN116705155A (zh) * 2023-08-03 2023-09-05 海南大学三亚南繁研究院 一种全基因dna数据的定义方法

Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
WO2013177615A1 (en) * 2012-06-01 2013-12-05 Agriculture Victoria Services Pty Ltd Selection of symbiota by screening multiple host-symbiont associations
CN106282330A (zh) * 2015-12-02 2017-01-04 香港中文大学深圳研究院 一种开发沙冬青植物基因组简单重复序列分子标记的方法
WO2018190170A1 (ja) * 2017-04-12 2018-10-18 花王株式会社 微細藻類における硝酸の基質アナログに対する耐性を向上させる方法

Family Cites Families (19)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN101504697B (zh) * 2008-12-12 2010-09-08 深圳华大基因研究院 一种片段连接支架的构建方法和***
WO2011143231A2 (en) * 2010-05-10 2011-11-17 The Broad Institute High throughput paired-end sequencing of large-insert clone libraries
US9506167B2 (en) * 2011-07-29 2016-11-29 Ginkgo Bioworks, Inc. Methods and systems for cell state quantification
US9567597B2 (en) * 2012-05-11 2017-02-14 University Of Hawaii Ultrasound mediated delivery of substances to algae
US10777301B2 (en) * 2012-07-13 2020-09-15 Pacific Biosciences For California, Inc. Hierarchical genome assembly method using single long insert library
WO2016192772A1 (en) * 2015-06-02 2016-12-08 Siemens Healthcare Gmbh Genetic testing for predicting resistance of shigella species against antimicrobial agents
WO2017012659A1 (en) * 2015-07-22 2017-01-26 Curetis Gmbh Genetic testing for predicting resistance of salmonella species against antimicrobial agents
WO2017016600A1 (en) * 2015-07-29 2017-02-02 Curetis Gmbh Genetic testing for predicting resistance of enterobacter species against antimicrobial agents
CA3010724A1 (en) * 2016-01-07 2017-07-13 Commonwealth Scientific And Industrial Research Organisation Plants with modified traits
CN107190003A (zh) * 2017-06-09 2017-09-22 武汉天问生物科技有限公司 一种高效快速分离t‑dna***位点侧翼序列的方法及其用途
WO2019005913A1 (en) * 2017-06-28 2019-01-03 Icahn School Of Medicine At Mount Sinai METHODS OF HIGH RESOLUTION MICROBIOME ANALYSIS
CN110042148B (zh) * 2018-01-16 2023-01-31 深圳华大基因科技有限公司 一种高效获取叶绿体dna测序数据的方法及其应用
CN108034706B (zh) * 2018-01-16 2021-03-26 浙江大学 利用重测序技术快速确定转基因株系***位点的方法
US11913006B2 (en) * 2018-03-16 2024-02-27 Nuseed Global Innovation Ltd. Plants producing modified levels of medium chain fatty acids
CN109295185B (zh) * 2018-09-05 2022-03-22 暨南大学 一种适用于单细胞真核藻类基因组大小的测定方法
CN114807398A (zh) * 2018-10-30 2022-07-29 厦门极元科技有限公司 宏基因组中沙门氏菌的鉴定方法及装置、分型方法及装置
CN111276185B (zh) * 2020-02-18 2023-11-03 上海桑格信息技术有限公司 一种基于二代高通量测序的微生物鉴定分析***及装置
CN111363706A (zh) * 2020-04-13 2020-07-03 天津中医药大学 旱莲草内生细菌、旱莲草组合物及其应用
CN111647680A (zh) * 2020-06-18 2020-09-11 北京市园林科学研究院 基于二代高通量测序的全基因组水平对苔草品种快速鉴定和溯源的方法

Patent Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
WO2013177615A1 (en) * 2012-06-01 2013-12-05 Agriculture Victoria Services Pty Ltd Selection of symbiota by screening multiple host-symbiont associations
CN106282330A (zh) * 2015-12-02 2017-01-04 香港中文大学深圳研究院 一种开发沙冬青植物基因组简单重复序列分子标记的方法
WO2018190170A1 (ja) * 2017-04-12 2018-10-18 花王株式会社 微細藻類における硝酸の基質アナログに対する耐性を向上させる方法

Also Published As

Publication number Publication date
CN112349350A (zh) 2021-02-09

Similar Documents

Publication Publication Date Title
CN112349350B (zh) 基于一种杜氏藻核心基因组序列进行品系鉴定的方法
Peng et al. Tracing the Austronesian footprint in Mainland Southeast Asia: a perspective from mitochondrial DNA
US20180258421A1 (en) Compositions, methods and uses for multiplex protein sequence activity relationship mapping
CN103088120B (zh) 基于SLAFseq技术的大规模样品基因分型方法
CN105112569A (zh) 基于宏基因组学的病毒感染检测及鉴定方法
Mark et al. Barcoding lichen-forming fungi using 454 pyrosequencing is challenged by artifactual and biological sequence variation
CN106868116A (zh) 一种桑树病原菌高通量鉴定及种属分类方法及其应用
CN106947827A (zh) 一种获得鳙性别特异分子标记及其筛选方法和应用
CN108103235A (zh) 一种苹果砧木抗寒性鉴定的snp分子标记、引物及其应用
CN109402241A (zh) 鉴定和分析古dna样本的方法
CN109112217A (zh) 一种与猪体长和***数显著关联的遗传标记及应用
Méndez-García et al. Metagenomic protocols and strategies
Lemoine et al. Assessing the evolutionary rate of positional orthologous genes in prokaryotes using synteny data
Kuster et al. ngsComposer: an automated pipeline for empirically based NGS data quality filtering
Xu et al. Genome reconstruction and haplotype phasing using chromosome conformation capture methodologies
CN111197050A (zh) 桑树拟干枯病病原菌的核糖体rna基因及其应用
CN110438244A (zh) 一种快速提高鸭群体青壳率的分子标记及应用
CN113564266A (zh) Snp分型遗传标记组合、检测试剂盒及用途
CN116083605B (zh) 一种包含67个高效能常染色体微单倍型的遗传标记体系及其检测引物和应用
CN104357563A (zh) 二次dna片段化的基因组单倍型高通量测序方法
US20220243267A1 (en) Compositions and methods related to quantitative reduced representation sequencing
CN107354151A (zh) 基于梅花鹿全基因组开发的str分子标记及其应用
Yang et al. A new perspective on codon usage, selective pressure, and phylogenetic implications of the plastomes in the Telephium clade (Crassulaceae)
CN112359102A (zh) 一种基于基因组学构建烟草核心种质的方法及其应用
Kust et al. Model cyanobacterial consortia reveal a consistent core microbiome independent of inoculation source or cyanobacterial host species

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant
TR01 Transfer of patent right
TR01 Transfer of patent right

Effective date of registration: 20231108

Address after: No. 9 Fulong Road, Shinan District, Qingdao, Shandong Province, 266000, 317

Patentee after: Qingdao Aixin Biotechnology Co.,Ltd.

Address before: 030006 No. 92, Hollywood Road, Taiyuan, Shanxi

Patentee before: SHANXI University