具体实施方式
除非另外限定,否则本文中所用的全部技术与科学术语具有如本发明所属领域的普通技术人员通常理解的相同含义。本文所提及的全部出版物、专利申请、专利和其他参考文献通过引用的方式完整地并入。此外,本文中所述的材料、方法和例子仅是说明性的并且不意在是限制性的。本发明的其他特征、目的和优点将从本说明书及附图并且从后附的权利要求书中显而易见。
定义
为了解释本说明书,将使用以下定义,并且只要适当,以单数形式使用的术语也可以包括复数,并且反之亦然。要理解,本文所用的术语仅是为了描述具体的实施方案,并且不意欲是限制性的。
术语“约”在与数字数值联合使用时意为涵盖具有比指定数字数值小5%的下限和比指定数字数值大5%的上限的范围内的数字数值。
如本文所用,术语“和/或”意指可选项中的任一项或可选项的两项或多项。
在本文中,当使用术语“包含”或“包括”时,除非另有指明,否则也涵盖由所述及的要素、整数或步骤组成的情形。例如,当提及“包含”某个具体序列的引物时,也旨在涵盖由该具体序列组成的引物。
除非另有说明,否则“一个”、“一种”、“该”和“至少一个”可互换使用并且表示一个或多于一个。
“核酸序列”指核苷酸(即,例 如,腺嘌呤(A),胸腺嘧啶(T),胞嘧啶(C),鸟嘌呤(G),和/或尿嘧啶(U))的任何聚合物。将包含脱氧核糖核苷的核酸序列称为脱氧核糖核酸(DNA)。将包含核糖核苷的核酸序 列被称为核糖核酸(RNA)。RNA可进一步表征为以下类型,(1)蛋白质编码RNA,也称为信使RNA (mRNA);(2)非编码RNA(non-codingRNA,ncRNA),ncRNA的核苷酸一般不长,共同特点是转录自基因组,但是不翻译成蛋白质,在RNA水平上行使各自的生物学功能,可以从长度上来划分为小于50nt的ncRNA,例如微量RNA(miRNA),小干扰RNA(siRNA)等;50nt到500nt的ncRNA,例如,转运RNA(tRNA)、核糖体RNA(rRNA)、小核(snRNA)和小核仁RNA(snoRNA)等;大于500nt的ncRNA,包括长的mRNA-like的非编码RNA、长的不带polyA尾巴的非编码RNA等。
如本文所用,术语“测序”通常是指用于确定一种或多种多核苷酸中的核苷酸碱基的序列的方法和技术。多核苷酸可以是核酸分子例如脱氧核糖核酸(DNA)或核糖核酸(RNA),包括它们的变体或衍生物(例如,单链DNA)。可以通过目前可用的各种装置执行测序,例如但不限于通过Illumina®、Pacific Biosciences(PacBio®)、Oxford Nanopore®或Life Technologies(IonTorrent®)的测序装置。可以使用核酸扩增、聚合酶链式反应(PCR)或等温扩增来执行测序。在一些实施方案中,此类装置提供测序读长(本文中也称为“read”)。读长(read)包括一串核酸碱基,其对应于已经测序的核酸分子的序列。
“二代测序”与“下一代测序(Next-generation sequencing,NGS)”可互换地使用,是一种基于PCR和基因芯片发展而来的测序技术。与一代测序通过合成终止测序不同的是:二代测序开创性的引入了可逆终止末端,从而实现边合成边测序(Sequencing bySynthesis)。二代测序在DNA复制过程中通过捕捉新添加的碱基所携带的特殊标记(一般为荧光分子标记)来确定核酸的序列。二代测序的出现将分子生物学研究推向了一个高通量发展的时代,利用NGS产生了大量基因组和转录组数据。但是,NGS测序平台的读长基本不超过500个碱基,是短读长测序,主要适用于检测单个细胞中的单核苷酸变异(SNV)、小的***缺失(<50bp,Indel)、以及拷贝数变异(CNV)等。通过NGS测序获得的序列准确性高,但覆盖度低。
三代测序(third-generation sequencing, TGS)是指单分子实时测序技术,也叫“从头测序技术”。与一代测序技术和二代测序技术相比,三代测序技术最大的特点就是单分子测序,测序过程中无需进行PCR扩增,实现了对每一条DNA分子的单独测序,克服了NGS测序平台读长短的缺点。目前TGS技术主要分为两类,一类是由美国太平洋生物科学公司(PacBio)推出的SMRT测序(Single Molecule Real-Time Sequencing, SMRT),即单分子实时测序,其中,核酸分子在被合成过程中,通过读取荧光信号识别碱基。另一类是由英国Oxford Nanopore公司推出的基于核苷酸通过纳米孔时合成不同碱基所产生的电信号差异以实现单分子测序得方法。三代测序是一种长读长(long-read)测序技术,现有技术中存在利用例如Nanopore MinION测序仪或PacBio平台实施三代测序。
术语“单细胞测序”是指对单个细胞的基因组进行测序。单个细胞是生物体的基本结构和功能单位。在相同的外部刺激或生理条件下,源自同一类型细胞的细胞可能会表现出细胞间差异,即,异质性(heteroplasmy)。构成例如大脑、血液***、免疫***的细胞之间存在异质性。如果能够从组织***中挑选多个不同的单细胞进行研究,则可以提供该组织中多个不同单细胞的有价值的信息,甚至能够重建出所述***。
术语“测序深度(Sequencing Depth)”是指测序得到的总碱基数与待测基因组或转录组或染色体区段大小的比值。假设一个基因组大小为2M,测序深度为10X,那么获得的总数据量为20M。测序深度能够用测序得到的碱基总量(bp)与待测基因组或转录组或染色体区段大小的比值来表示。
术语“覆盖度(coverage)”是指,基因组或转录组或染色体区段上获知序列信息的序列部分占整个组或区段的比例。在一些实施方案中,覆盖度是指,(例如通过序列检测手段,如测序)测到序列信息的碱基数占所检测区域的总碱基数的比例。例如,在测序检测全基因组序列时,由于存在大片段拼接的缺口(gap)、测序读长有限、重复序列等问题,测序后得到的基因组序列通常无法完全覆盖基因组的所有区域,此时,覆盖度就是最终得到的测序碱基数占整个基因组碱基数的比例。例如,对人的基因组测序获得的覆盖度为98.5%,这表明该基因组还有1.5%的区域的序列没有得到。在另一些实施方案中,覆盖度是指,就所检测的区域而言,例如通过测序分析测到序列信息的基因位点(例如SNP位点或基因变异位点)的数目,占该区域中所检测的总基因位点数的比例。所检测区域可以是全基因组、特定染色体、或特定染色体区段、或转录物组、或特定转录区域。
“文库构建”是指在DNA的一端或两端加上可以用于测序以及文库拆分的接头(adaptor)结构。在一个实施方案中,所构建的文库由P5/Truseq Read/TSO/barcode/UMI和/或DNA***片段(DNA Insert)组成,其中,P5与流动池(flow cell)上的接头引物序列互补和相同,进行簇(cluster)生成;Truseq Read与全长Illumina TruSeq Read 1 测序引物互补和相同。
单细胞测序的基础流程是:获得单细胞,单细胞裂解,对单细胞内的核酸进行扩增,以达到测序的最小上样量的要求,然后上机,读取数据并分析。
可以通过流式细胞术(FACS)、激光捕获显微切割技术(Laser capturemicrodissection,LCM)、微流控技术(Microfluidics)获得单细胞,其中微流控技术是一种用于精确控制微量液体的技术,微流控芯片是实施该技术的平台。微流控芯片通过细微的管道对液体实施操控,微流控对液体的操控程度刚好适合于单细胞样品的处理操作。
细胞裂解方法可以是物理方法、化学方法或酶法。存在三种主要的物理细胞裂解形式:机械裂解、热裂解和电裂解。化学裂解应用裂解缓冲液并诱导高效裂解以破坏细胞。酶促细胞裂解是减少DNA断裂的最温和的方法。可以使用蛋白酶,如胃蛋白酶和胰蛋白酶。
可以通过全基因组扩增(Whole Genome Amplification,WGA)和/或全转录组扩增(Whole Transcriptome Amplification,WTA)对细胞内的核酸进行扩增。全基因组扩增方法包括例如简并寡核苷酸引发的PCR(Degenerate Oligonucleotide Primed PCR,DOC-PCR)、多重置换扩增(Multiple Displacement Amplification,MAD)、基于多重退火和成环的扩增循环(MALBAC)。全转录组扩增包括对单细胞中提取的RNA逆转录出cDNA,然后对逆转录生成的cDNA进行扩增。
单细胞基因组测序,即scDNA-seq,用来获得细胞基因组的突变和结构变化。
单细胞转录组测序,即 scRNA-seq,可以以单细胞分辨率来测量转录组中的基因表达,识别细胞簇中的生物学相关差异。单细胞RNA测序方法包括例如Smart-seq2、CEL-seq2、sci-RNA-seq、10x Chromium、Drop-seq、Seq-Well、inDrops,其中Smart-seq2和CEL-seq2是基于微孔板的低通量的方法,其他5种方法是高通量方法。
如本文所用,术语“条形码(barcode)”是指能够提供分析物的信息的标记或标识符。条形码可以是分析物的一部分。条形码也可以独立于分析物。除了分析物的内源特征(例如,分析物的大小或一个或多个末端序列)之外,条形码可以是连接至分析物(例如,核酸分子)的一个标签或者是多个标签的组合。条形码是独特的。条形码可以有多种不同的形式。例如,条形码可包括多核苷酸条形码;随机核酸和/或氨基酸序列;和合成核酸和/或氨基酸序列。条形码可以以可逆或不可逆的方式连接于分析物。在样品测序之前、期间和/或之后,可以将条形码添加到例如脱氧核糖核酸(DNA)或核糖核酸(RNA)样品的片段中。条形码可以允许确定和/或定量各个测序读长。在一些实施方案中,当条形码是一段核苷酸序列时,任意两个条形码之间至少存在两个或者两个以上碱基的不同,因为测序存在对碱基的误读,这样的设计可以避免将两个条形码混淆。对于一段长16bp的核苷酸序列,大约存在350万种条形码(barcode)。在一个实施方案中,barcode用来区别GEMs,是对一个单细胞分配的标记。
术语“固定”是指固定细胞的任何方法或过程。在本文中与“交联”可互换地使用。许多化学品能够用来固定细胞,包括但不限于,甲醛、***、或戊二醛。
术语“富集”指对线粒体核酸的任何分离或目的线粒体核酸区段的比例相对于核酸组合物中的其他核酸区段增加。例如,富集或分离目的线粒体核酸核酸或核酸区段可以包括正向方法,如“分出”目的线粒体核酸核酸或核酸区段,或可以包括负向方法,如排除不作为目的线粒体核酸的核酸或排除不包含线粒体核酸区段的核酸区段。备选地,富集包括选择性或定向扩增目的线粒体核酸或核酸区段。这类选择性或定向扩增目的线粒体核酸将增加线粒体核酸在核酸组合物中的比例(即富集所述线粒体核酸)。
“单核苷酸多态性”或“SNP”是基因组内的单核苷酸变化(即A,C,G或T),其在生物物种的成员之间或成对染色体之间不同。
本文提及的“纯化”,可以指核酸已经历处理(即,例如分级)以除去各种其它组分。当使用术语“基本上纯化的”时,表示目的核酸构成组合物中的主要组分,诸如构成组合物的约50%、约60%、约70%、 约80%、约90%、约95%或更多(即,例如,重量/重量和/或重量/体积)。
多重 PCR (Multiplex polymerase chain reaction,MPCR) 是将不同的引物对和模板分散于相对独立的空间中进行PCR扩增。例如,在同一个反应体系中加入不同的引物对,针对不同的模板或同一模板的不同区域进行特异性的PCR扩增,从而得到多个目的扩增产物,结合一定的检测手段能够实现同时对多个靶标进行检测。MPCR 具有高效率、高通量、低成本的特性。
术语“通用引物结合位点”是指存在于(通常,通过人工添加至)不同的许多目标核酸中的、与通用引物结合的位点。“通用引物”通过结合“通用引物结合位点”指导引物延伸。测序时所使用的通用引物也称为通用测序引物。PCR时所使用的通用引物也称为通用PCR引物。通用PCR引物和通用测序引物可以相同或不同。
术语受试者的“样本”是指例如来自所述受试者个体的血液、唾液、口腔拭子、尿液、指甲、毛囊、皮屑、组织、体液样本。
术语“混合样本(pooling)”是指在进行高通量测序时,对多个样本进行混合,以最大化利用测序仪实施高通量测序。在本发明的方法中可以混合任意数目的样本,不局限于混合几个样本。
术语“样本索引(sample index)”是在三代测序文库构建中使用的样本标记物,用于区分不同样本的三代测序数据。样本索引(sample index)通常是一段碱基组成的DNA序列。在三代测序文库构建过程中,每个样本的文库都要分别带上不同的sample index序列,像样本的“身份证号”一样,能将流动池(flow cell)里所测得的海量序列正确地拆分。选择合适的样本索引集合(sample index sets)以确保在多重测序(multiplexed sequencing)运行中多个样本索引之间没有重叠。
术语“流动池(flow cell)”是指固体表面的腔室,一种或多种流体物质可以流过该固体表面。流动池以及相关的流体***和检测平台的实例例如在US 7,329,492;US 7,211,414中描述。
术语线粒体“核酸文库”是本发明创建的线粒体核酸的集合,可以通过生物合成来制备。
术语“Cell Ranger”是10X genomics公司提供的一套针对单细胞 RNA 测序输出结果进行比对、定量、聚类及基因表达分析的分析流程,它包含有与单细胞基因表达分析相关的四个pipelines,分别是:
cellranger mkfastq 流程:其功能为将 Illumina 测序仪产生的raw base call(BCL) 文件解析成 FASTQ 文件;
cellranger count 流程:其功能为将 cellranger mkfastq 产生的或其他来源的FASTQ 文件进行比对、过滤、barcode 计数以及 UMI 计数,并可以生成 feature-barcode 定量矩阵,随后确定细胞群并进行基因表达分析;
cellranger aggr 流程:其功能为将多个cellranger count 产生的数据进行整合、标准化,并可以对整合后的数据进行分析;
cellranger reanalyze 流程:其功能为使用 cellranger count 或 cellrangeraggr 产生的表达矩阵重新进行降维、聚类等后续分析。
以上四个管线(pipeline)均将转录组常用比对软件STAR封装其中,所述STAR用来判断reads是否比对到了基因组上。
术语“GATK”是 Genome Analysis ToolKit的缩写,是一款从高通量测序数据中分析变异信息的软件(下载链接:https://software.broadinstitute.org/gatk/download/),可以用于寻找SNP和Indel,一个标准的分析流程称为GATK Best Practise,主要包括以下几个步骤:
数据预处理:对从测序仪下机后的数据进行质控,去除低质量的reads,将过滤后的reads比对到参考基因组上,产生BAM格式的比对文件。
寻找变异:进行变异召回(variant calling),寻找例如SNP和Indel,将比对数据存储在VCF格式的文件中。
使用寻找出的变异位点进行后续的分析。
术语“Fastq”是序列数据存储的标准格式之一,每4行为一条read的信息,包含测序read名、序列、正反链标示、序列质量值。
本发明的确定线粒体转录组中的变异和确定多细胞真核生物中的单细胞谱系的方法
线粒体在细胞功能中起着重要的作用,不仅提供细胞所需90%以上的能量,还是细胞凋亡启动和执行的主要细胞器之一。人类线粒体基因组包括16569个碱基对,其中包含2个转录为核糖体RNA的基因、13个转录为线粒体信使RNA并编码蛋白质的基因、以及22个tRNA基因。人类线粒体基因组中的突变会导致多种人类疾病,例如,癌症、心脏病、糖尿病等。
本发明提供了使用线粒体基因组中的体细胞突变来推断细胞谱系的方法,所述线粒体基因组中的体细胞突变也可以用作细胞的遗传条形码(genetic barcode)。在本发明的方法中,线粒体基因组中的体细胞突变是通过对单细胞转录组富集线粒体cDNA并进行二代测序和三代测序确定的,基于细胞中突变的存在对细胞进行聚类,由此推断细胞的谱系。本发明的方法可以确定线粒体DNA中的缺失(deletion)、***(insertion)、重复(duplication)和倒位(inversion)等突变。
本发明提供了确定线粒体转录组中的变异和确定多细胞真核生物中的细胞谱系的方法,所述方法的一个示例性流程图如图1所示,其涉及以下方面。
1) 利用酶法或者物理方法制备多细胞真核生物的单细胞悬液。
2) 捕获单细胞,然后裂解细胞,使细胞mRNA自细胞中游离。
在一个实施方案中,使用10x Genomics公司的仪器捕获细胞;将获得的细胞按照10x Genomics 操作说明书进行上机,可捕获的细胞数为500个细胞-10000个细胞。该方案也同样适用于其他单细胞WTA(Whole transcriptome analysis)捕获平台如Drop-seq、InDrop-seq、BD Rhapsody及smart-seq等。
3) 逆转录细胞mRNA,获得3’端或5’端细胞条形码标记的全长cDNA。
在一些实施方案中,获得如图2所示的3’/5’细胞条形码标记的全长cDNA。
在一些实施方案中,使用了来自10x Genomics公司的产品制备单细胞转录组,利用微流控技术将单细胞与反应试剂包裹形成GEMs(Gel bead in emulsion),每一个GEM里面除了包含一个细胞和试剂外,还具备带有双层标签的Gel beads oligo,第一层标签是10x barcode,用来标记各个细胞,即,每个GEM里面的barcode是唯一的;第二层标签是用来进行表达定量的UMI,用来标记逆转录的cDNA分子。例如,使用10x Chromium Next GEMSingle Cell 5' Kit v2试剂盒(10x Genomics公司,目录号:PN-1000263),在10xChromium Controller上机结束后会形成GEMs,按照相应的说明书操作流程将GEMs进行RT反应程序孵育及cDNA扩增,即可获得带有单细胞barcode的全长cDNA。具体地,在GEM生成后,将凝胶珠(Gel Bead)溶解并使任何共同分配在该凝胶珠中的细胞裂解。释放出包含(i)Illumina R1序列(结合read 1测序引物)、(ii) 16 nt 10x Barcode、(iii) 10 nt UMI和(iv) 13 nt TSO的寡核苷酸,将其与细胞裂解物和含有逆转录 (RT)试剂和poly(dT) RT引物的Master Mix混合。孵育GEM,获得反转录自聚腺苷酸化mRNA的包含10x Barcoded的全长cDNA。
通过所述步骤1)至步骤3)可以获得细胞cDNA,其包含核基因转录组和线粒体基因转录组的cDNA。在一个实施方案中,获得的细胞全长转录组cDNA如图2所示,其包含条形码(barcode)标记。该部分适配于市面上大多数单细胞测序技术,包括10x genomics 3’RNAkit、10x genomics 5’RNA kit、10x genomics multiome kit、smart-seq、smart-seq2、smart-seq3、drop-seq、indrop、microwell-seq等可产生单细胞经标记的全长cDNA的产品或者技术。
4) 设计线粒体特异性引物,通过PCR富集线粒体cDNA。
在一些实施方案中,当步骤3)中的细胞条形码序列标记在cDNA的5’端时,对于每个线粒体基因设计一条或者多条位于cDNA 3’端的特异性引物,利用3’端特异性引物与5’端通用PCR引物进行该线粒体基因的PCR反应。在特异性引物结合位置方面要求尽量靠近基因的3’端,同时所述引物具有较高的特异性,能够区分线粒体基因组和核基因组序列。
优选地,设计并合成针对线粒体DNA(mtDNA)上各线粒体基因的引物。线粒体DNA包含如下线粒体基因:13个编码基因,2个rRNA基因,以及22个tRNA基因。本发明的方法由于仅需要对每个线粒体基因合成1条下游反向引物,由此降低了成本。由于所述引物是特异性靶向线粒体转录组的引物,因此具有富集线粒体基因的作用。
又在一些实施方案中,当步骤3)中的细胞条形码序列标记在cDNA的3’端时,对于每个基因设计一条或者多条位于cDNA 5’端的特异性引物,利用5’端特异性引物与3’端通用PCR引物进行该基因的PCR反应。在特异性引物结合位置方面要求尽量靠近基因的5’端,同时所述引物具有较高的特异性,能够区分线粒体基因组和核基因组序列。
优选地,设计并合成针对线粒体DNA(mtDNA)上各线粒体基因的引物。线粒体DNA包含如下线粒体基因:13个编码基因,2个rRNA基因,以及22个tRNA基因。由于仅需要对每个线粒体基因合成1条上游反向引物,由此降低了成本。由于所述引物是特异性靶向线粒体转录组的引物,因此具有富集线粒体基因的作用。
在一些实施方案中,当步骤3)中的细胞条形码序列标记在cDNA的5’端且获得的全长转录组cDNA如图2所示时,为获得带有单细胞条形码(Single cell barcode)标记的线粒体全长cDNA,线粒体基因的靶向富集用上游引物使用的是通用PCR引物。在一个实施方案中,所述通用PCR引物为Read 1引物序列(5’-CTACACGACGCTCTTCCGATCT-3’(SEQ ID NO:28))。因此仅需要设计每个线粒体基因的下游引物以特异性富集线粒体基因的全长cDNA。在一个实施方案中,设计的下游反向引物的位置尽量靠近mtDNA基因的3’端200bp范围内,例如,尽量靠近mtDNA基因的3’端约200bp范围内,例如,约150bp范围内、约100bp范围内、约50bp范围内、约30bp范围内,避开SNP位点,且能够与核基因组进行区分,以获得mtDNA基因的全长转录组。
在一些实施方案中,本发明设计了如表1所示的具有线粒体基因特异性的下游反向引物。对于每种mtRNA反转录获得的cDNA,仅需设计一条针对其的线粒体特异性反向引物(reverse primer)即可特异性扩增该cDNA。鉴于线粒体基因与细胞核基因具有较高的同源性,因此,需要设计仅对特定的线粒体基因序列具有更强的特异性,从而能够与核基因组进行区分的引物。
现有技术中没有公开对线粒体基因的全长cDNA进行特异性富集,然后用于测序的技术方案。例如,MEASTER技术(参见 Miller, T. E.等人, "Mitochondrial variantenrichment from high-throughput single-cell RNA-seq resolves clonalpopulations", Nat Biotechnol (2022). https://doi.org/10.1038/s41587-022-01210-8)需要设计65对引物(参见 Miller, T. E.等人,第14页图1)以满足短序列测序时覆盖全长线粒体cDNA的要求。与MEASTER技术针对每个线粒体基因设计多个引物相比,本发明的引物能够特异性扩增出全长线粒体cDNA,而不是如MEASTER技术那样仅扩增约250 bp长的线粒体cDNA片段。
进一步地,在富集线粒体cDNA时,本发明人为了解决构建的测序文库中各线粒体转录本信息可能存在覆盖度不均一的问题,对步骤4)中用于PCR扩增各线粒体cDNA的富集用特异性引物的比例进行调整,旨在得到数量均等的线粒体各转录本来源信息,以防止引入PCR偏差(PCR bias),从而有利于后续分析的准确性。
在一些实施方案中,将用于PCR扩增各线粒体cDNA的富集用特异性引物的比例参照HPA数据库https://www.proteinatlas.org/中线粒体各基因的表达水平进行调整。
在一些实施方案中,当通过PCR扩增特异性富集各线粒体cDNA,例如,富集选自线粒体13个编码基因、2个rRNA基因的cDNA时,将各特异性引物依照各基因表达水平来调整所述特异性引物的使用浓度,例如通过HPA数据库https://www.proteinatlas.org/调取线粒体各基因的表达水平,由此设定针对不同基因的特异性引物的添加量,具体如下表A所示。
表A. 线粒体基因表达量与扩增用引物的终浓度范围
待扩增的线粒体cDNA所对应的基因名称 |
线粒体基因表达量(覆盖深度)占比(%) |
扩增用引物的终浓度范围uM |
MT-RNR1 |
0.07~0.55 |
0.07~1.02 |
MT-RNR2 |
0.25~2.62 |
0.05~1.7 |
MT-ND1 |
0.19~2.12 |
0.3~2.68 |
MT-ND2 |
0.26~1.81 |
0.27~0.48 |
MT-CO1 |
0.83~4.86 |
0.05~0.12 |
MT-CO2 |
0.44~3.09 |
0.06~0.47 |
MT-ATP8 |
0.42~4.15 |
0.16~0.5 |
MT-ATP6 |
0.5~3.72 |
0.08~0.34 |
MT-CO3 |
0.57~3.93 |
0.07~0.2 |
MT-ND3 |
0.12~1.84 |
0.09~0.4 |
MT-ND4L |
0.23~2.32 |
0.17~0.57 |
MT-ND4 |
0.39~3.07 |
0.3~0.76 |
MT-ND5 |
0.12~1.23 |
0.3~1.55 |
MT-ND6 |
0.13~1.63 |
0.5~2.82 |
MT-CYB |
0.24~2.11 |
0.2~0.6 |
使用表A所示的针对各线粒体基因的扩增用引物的添加量对样本进行线粒体cDNA富集。经过测序分析不同线粒体基因的覆盖深度,由此对各特异性引物的添加量优化而获得各特异性引物的经调整的添加量。对其他样本同样地获得各特异性引物的经调整的添加量,也能实现稳定且相对均一的线粒体基因富集。
按所获得的针对线粒体cDNA的各特异性引物经调整的添加量将多条特异性引物(例如,针对线粒体13个编码基因、2个rRNA基因的cDNA的15条特异性下游反向引物)混合在1个PCR管中,补加无核酸酶的ddH2O至总体积为100ul,获得下游反向引物Mix。在一个实施方案中,如下进行线粒体转录组扩增PCR1,为1管反应/样本,PCR1扩增程序为:
预变性:95℃ 3min
6个循环(98℃ 30s;65℃ 30s;72℃ 3min)
终延伸:72℃ 5min
保存温度:4℃。
线粒体转录组扩增PCR1反应结束之后,对PCR1反应产物进行纯化,纯化产物作为后续反应的模板。在一些实施方案中,使用1.0XSPRI磁珠(Solid Phase ReversibleImmobilization,Beckman,B23317)或1.0X Ampure XP磁珠(Beckman,A63880)等进行纯化,得到靶向富集的线粒体cDNA。
与本发明可以便捷、省时地得到引物Mix相比,在使用MAESTER技术(参见Miller,T. E.等人,见上文)时,将65条引物溶解后,按固定比例配制成12管不同的下游反向引物Mix(参见Miller, T. E.等人,第14页,图1),通过补加不同体积的缓冲液,保证每条引物在扩增用混合物中的浓度为1μM。随后,每个样本需要配制12管扩增用混合物进行线粒体转录组扩增的PCR1反应。MAESTER技术的PCR1反应为:
预变性:95℃ 3min
6个循环(98℃ 20s;65℃ 15s;72℃ 3min)
终延伸:72℃ 5min
保存温度:4℃。
MAESTER技术在PCR1反应结束后,基于获自测序的覆盖率而凭经验确定需要从12管PCR产物分别取不同体积的产物量混合为1管(参见Miller, T. E.等人,第14页,图1)进行下一步纯化,纯化产物作为PCR2的反应模板,进行PCR2扩增反应。MAESTER技术的PCR2反应为:
预变性:95℃ 3min
8个循环(98℃ 20s;60℃ 30s;72℃ 3min)
终延伸:72℃ 5min
保存温度:4℃。
因此,与采用MAESTER技术相比较,本发明的步骤4)具有耗材用量少,节省成本,配制时间短的优点;以及在样本数量较多时,操作更快,且更节省时间的优点。
5) 对经纯化的步骤4)的线粒体转录组富集产物实施下述i)和ii),进行二代Bulk线粒体转录组片段测序和三代单细胞线粒体全长转录组测序。
i) 取一部分步骤4)的线粒体转录组富集产物进行三代单细胞线粒体全长转录组测序。
在一个实施方案中,如下实施步骤5)的i):取一部分步骤4)的线粒体转录组富集产物进行单一样本长读长三代测序建库并实施三代测序。例如,按照ONT连接测序试剂盒(SQK-LSK110)操作说明书执行末端修复和测序接头连接,实施单一样本三代测序文库构建,然后上机测序。
在一个实施方案中,如下实施步骤5)的i):取一部分步骤4)的线粒体转录组富集产物作为来自目的样本的线粒体转录组富集产物,与作为来自其他非目的样本的线粒体转录组富集产物混合,进行多个样本(即,包括目的样本和非目的样本)的线粒体转录组长读长三代测序建库并实施三代测序。具体地,通过PCR扩增的方式向线粒体转录组富集产物的混合物(即,此轮PCR的模板)添加样本索引(Sample index)序列,由此,获得带样本索引序列的多个样本单细胞线粒体全长转录组扩增产物,然后,进行长读长混合样本三代测序建库并实施混合样本三代测序。在一些实施方案中,混合多个样本,例如,混合1、2、3…、N个样本的线粒体转录组富集产物,其中N是正整数。例如,按照ONT连接测序试剂盒(SQK-LSK110)操作说明书执行末端修复和测序接头连接,实施文库构建,然后上机测序。
ii) 取另一部分步骤4)的线粒体转录组富集产物进行Bulk线粒体转录组片段构建,例如,使用NEBNext® Ultra™ II FS DNA Library Prep Kit for Illumina(NEB,E7805)按照产品说明书分别执行片段化、末端修复和加A尾、接头(adaptor)连接、及PCR扩增,构建二代测序文库,然后使用Illumina测序平台对文库进行测序,以获得Bulk线粒体转录组短读长二代测序数据。
具体而言,将步骤4)获得的线粒体cDNA富集产物进行片段化,筛选满足二代测序读长的片段范围进行片段选择,然后通过接头连接及PCR扩增的方式得到所有的Bulk线粒体转录组二代文库,经二代测序获得的数据中是单个样本的一群细胞的线粒体转录组短读长的数据。在本发明中,将对单个样本中多个细胞的具有细胞条形码(barcode)序列和UMI序列的线粒体cDNA富集产物进行二代建库和测序也简称为Bulk线粒体转录组二代建库和测序,所获得的测序数据也称为“Bulk线粒体转录组二代测序数据”,与之相关的二代测序也称为“短读长Bulk线粒体转录组二代测序”、“二代Bulk线粒体转录组测序” 、“Bulk线粒体转录组短读长测序”或‘“Bulk线粒体转录组短读长测序”。
在一个实施方案中,步骤5)的i)采用的是ONT(Oxford Nanopore Technology公司)三代长读长的测序技术,可以用于分析单细胞全长线粒体cDNA,测序数据准确度也随着测序技术的发展和分析方法的优化得到改善,结合步骤5)的ii)Illumina Novaseq6000平台的二代测序数据共同分析,能够得到更多更准确的数据。
现有技术例如MAESTER采用的是二代测序技术分析线粒体基因组突变,虽然二代测序数据具有较高的测序数据准确度,成本低,但同时也存在读长短,信息丢失的弊端,因此,本发明的方法克服了现有技术中的不足。
6) 任选地,以步骤3)中制备的单细胞全长cDNA转录组为模板,进行二代短读长单细胞转录组测序。
在一个实施方案中,按照10x Genomics CG000331 Rev B中的Gene Expression(GEX) Library Construction制备单细胞转录组文库,其流程如图4所示,将单细胞全长转录组cDNA进行片段化,然后经过末端修复、加A碱基及接头连接,即可在后续PCR过程中富集含有cell barcode及UMI信息的基因表达片段。也可以使用其他能够制备单细胞转录组二代测序文库的试剂。单细胞转录组文库随后进行二代短读长测序,与步骤5)的ii) 短读长Bulk线粒体转录组二代测序数据及步骤5)的i)三代长读长测序数据共同分析。
所述步骤6)是任选实施的步骤,用于获得靠近转录组3’端或5’端序列信息 ,以进一步校正三代测序reads中的barcode和UMI的信息。
另外,所述步骤步骤5)的i)和步骤5)的ii)可以以任何顺序实施;所述步骤6)也可以在步骤4)之前实施,并且所述步骤4)至步骤6)可以以不同的组合二代测序和三代测序的方式进行实施。
在一些实施方案中,如图3和图4所示实施所述步骤4)至步骤6)。在一个具体实施方案中,图3中的Bulk线粒体短读长序列(Short reads sequence)二代测序文库是以线粒体转录组富集产物为模板,使用NEBNext® Ultra™ II FS DNA Library Prep Kit forIllumina(NEB,E7805)试剂盒构建的,因此,所获得的Bulk线粒体转录组二代测序数据中仅含有线粒体基因的数据,并且通过本发明的线粒体转录组富集步骤,可以确保样本的线粒体转录组文库中各片段的二代测序深度一致,更利于生物信息的分析。在一个具体实施方案中,图4中的单细胞转录组短读长序列(Short reads sequence)二代测序文库是以单细胞全长cDNA转录组为模板,使用10x Genomics CG000331 Rev B中的Gene Expression(GEX) Library Construction Kit, PN-1000190构建的。由于单细胞转录组二代测序数据中主要包含的是核基因表达的信息,其中也包含部分线粒体数据,但通常是否包含所有线粒体基因数据是未知的,且样本的线粒体各基因测序深度也难以一致。在常规单细胞转录组分析流程中,在并不需要线粒体数据时,线粒体数据往往作为背景噪音被去除。因此,如果仅采用单细胞转录组二代测序数据进行线粒体信息分析,可能会存在有效的线粒体转录组二代测序数据极少且不全面等技术问题。将本发明的所述步骤4)至步骤6)组合实施克服了该现有技术中存在的技术问题。
7) 对步骤5)和任选地步骤6)生成的测序数据进行生物信息分析。
本发明利用长读长三代测序获得的读长更长、且得到全长转录组信息的优点,对单细胞线粒体全长转录组进行了三代测序(所获得的reads中也包含barcode/UMI信息)。另一方面,由于长读长三代测序,特别是Oxford Nanopore测序的单碱基准确率低,在barcode、UMI及全长转录组都会出现测序错误,因此,在一个实施方案中,使用步骤5)的ii)中得到的Bulk线粒体转录组二代测序数据,对步骤5)的i)的长读长三代测序单细胞线粒体转录组的barcode/UMI和变异信息(SNP/InDel)进行校正,以得到携带有高准确度的barcode/UMI、SNP/InDel的全长单细胞线粒体转录组序列。
又在一个实施方案中,使用步骤5)的ii)中得到的Bulk线粒体转录组二代测序数据和进一步地使用步骤6)中得到的单细胞转录组短读长二代测序数据,对步骤5)的i)的长读长三代测序单细胞线粒体转录组的barcode/UMI和变异信息(SNP/InDel)进行校正,以得到携带有高准确度的barcode/UMI、SNP/InDel的全长单细胞线粒体转录组序列。图5例示了该技术方案中的一种生物信息分析框架图,具体如下所述。
a模块. 长读长测序单细胞线粒体转录组数据分析,以期获得单细胞线粒体全长转录组水平的变异信息(Nanopore 测序因为单碱基准确性低,所以会携带低准确性barcode/UMI序列、低准确性 SNP/InDel位点):
(a) 采用Smith-Waterman局部比对算法(https://github.com/mengyao/Complete-Striped-Smith-Waterman-Library) 将样本index序列与read两端 300bp的序列进行比对,允许4 个碱基错配,将 read分配到具有最佳比对得分的样本;
(b) 单个样本采用 minimap2软件(https://github.com/lh3/minimap2) 将reads比对到参考基因组(例如,GRCh38);
(c) 提取read 5’端未比对到参考基因组的序列(定义为softclip5, 包含例如P5/Truseq Read/TSO/barcode/UMI等通过10X Genomics技术添加的序列);
(d) 在比对完成的bam文件中提取转录本比对到参考基因组的起始位置plrs。
b模块. 短读长测序单细胞转录组数据提供高准确性barcode/UMI序列:
(a) 使用 fastp软件(https://github.com/OpenGene/fastp)选择默认参数对单细胞转录组短读长测序下机的原始数据进行质控,去除低质量 reads;
(b) 质控后的数据使用 Cell Ranger进行处理,所述Cell Ranger 是 10Xgenomics公司提供的一套针对单细胞 RNA 测序输出结果进行比对、定量、聚类及基因表达分析的分析流程,由此获得高质量的barcode/UMI序列;
(c) 质控后的序列使用 STAR软件(https://github.com/alexdobin/STAR)比对到参考基因组(例如,GRCh38),提取转录本比对到参考基因组的起始位置psrs。
c模块. 短读长测序Bulk 线粒体转录组数据提供高准确性的barcode/UMI序列及基因水平的变异:
(a) 使用 fastp软件,选择默认参数对短读长测序Bulk 线粒体转录组的原始数据进行质控,去除低质量 reads;
(b) 采用STAR软件将序列比对到参考基因组(例如,GRCh38),随后采用GATK4Mitochondria Pipeline (https://github.com/gatk-workflows/gatk4-mitochondria-pipeline)进行变异检测,得到高准确性的SNP/InDel;
(c) 质控后的数据使用 Cell Ranger进行处理,由此获得高质量的barcode/UMI序列。
d模块. 长读长测序单细胞线粒体转录组数据与短读长测序单细胞转录组数据整合,用来校正线粒体全长转录组的barcode/UMI序列:
(a) 根据基因名称将单细胞短读长测序的barcode/UMI 及长读长测序softclip5进行分组。同一基因的barcode/UMI与softclip5具有相同的转录本来源,进行比对时,可以减少比对时间,提高比对效率;
(b) 将长读长测序单细胞线粒体转录本数据的softclip5序列分为P5/TruseqRead/TSO等模板序列(有效 reads应携带这些序列,且所有 reads的这些序列相同)与barcode/UMI序列。将模板序列作为锚点序列,锚点序列允许更多的错配(15%的错配率),可以获得更多有效 reads,同时定位到 barcode/UMI序列在softclip5的位置。根据此位置,获得softclip5的barcode/UMI序列,进而与短读长测序单细胞转录本的相同基因来源的barcode/UMI序列进行两两比对;
(c) 整合长读长、短读长barcode/ UMI 比对的结果(mismatches_gaps)与转录本在参考基因组的起始位置(pos_diff =|psrs-plrs|), 比较不同mismatches_gaps与pos_diff下的长读长测序reads数目,最终选择(pos_diff<30 & mismatches_gaps<=3) |(pos_diff<20 & mismatches_gaps<=4) |(pos_diff<10 & mismatches_gaps<=5)作为长读长测序与短读长测序 barcode/UMI 映射的标准。
e模块. 对d模块得到的barcode/UMI校正后的单细胞全长转录本序列与c模块得到的高准确性的SNP/InDel进行整合,预期得到高准确性的SNP/InDel的全长转录本:
(a) 以短读长测序 Bulk 线粒体转录本数据得到的SNP/InDel为依据,确定 SNP/InDel 在长读长测序线粒体转录本单条reads 中是否存在;
(b) 携带相同 barcode/UMI的不同 reads 来自于相同细胞的相同转录本,对这样的reads进行SNP/InDel的相互校正,得到高准确性的变异结果。
f模块. 根据上述各模块得到的携带高准确度 barcode/UMI、高准确度 SNP/InDel的全长转录本计算单个细胞(barcode)水平 SNP/InDel的变异频率,构建细胞-变异(SNP/InDel) 矩阵,并根据细胞间的欧式距离 (Euclidean distance),进行层次聚类分析,据此推断细胞的进化谱系。
本发明中采用的线粒体测序技术是二代测序与三代测序的结合。在线粒体转录组三代测序数据中需要单细胞条形码信息,以便进行三代单细胞线粒体全长转录组的序列分析。在线粒体三代测序数据中,单细胞条形码(Single cell barcode)能够将三代测序数据归类到单细胞的水平。另外,在Bulk线粒体二代测序数据中能够提供高准确性单细胞变异位点(例如,SNP/InDel位点)和单细胞条形码(Single cell barcode)信息,因此,本发明的单细胞线粒体转录组的三代测序数据可以被Bulk线粒体转录组二代测序数据校正。任选地,对单细胞转录组靠近3’端或5’端的序列进行二代测序,也可以进一步校正三代测序reads中的barcode和UMI的信息。
本发明中将步骤5)的ii)的Bulk线粒体转录组二代短读长测序与将步骤5)的i)单细胞线粒体全长转录组的长读长三代测序(例如,ONT纳米孔测序技术、Pacbio单分子测序技术或者其他长读长测序技术)组合,任选地,还与步骤6)的单细胞转录组靠近3’端或5’端的序列的短读长二代测序组合,能够将短读长的二代测序技术与长读长的三代测序技术的各自优势结合起来。
对于同一样本同时执行Bulk线粒体转录组二代测序和单细胞线粒体全长转录组三代测序,以及任选地,还执行单细胞转录组靠近3’端或5’端的序列的短读长二代测序,具有如下优势。其一,能够准确地获得线粒体全长转录组信息;其二,利用同一样本中多个单细胞的Bulk线粒体转录组进行的二代测序结果作为线粒体变异的参照,可显著地降低三代测序变异位点的假阳性,能够降低三代测序所需测序数据量和测序成本。
现有的MAESTER技术仅利用了二代测序平台。二代测序准确度高,但读长短,往往会丢失全长信息,MAESTER技术可能会造成全长线粒体转录组覆盖不完全、测序深度不均一的问题,进而漏检部分真实变异。
本发明的方法将二代测序平台和三代测序平台的测序结果结合起来进行分析,二代测序数据可以用于纠正三代测序数据,降低线粒体变异位点的假阳性;二代测序数据结合三代测序数据的全长特点,可以获得更多线粒体转录组的信息。
因此,本发明的二代测序平台和三代测序平台的组合能获得更加准确的线粒体全长转录组序列信息。
本发明的试剂盒
本发明提供了试剂盒。试剂盒可以包括一种或多种实施本发明方法的物品和/或试剂。
在一个实施方案中,试剂盒用于制备线粒体转录组测序文库,其包含含有通用PCR引物结合位点、用于确定转录组来源细胞的细胞条形码(barcode)序列和用于确定所述来源细胞中的各转录本的唯一分子标识符(UMI)序列的核苷酸序列,用于添加至线粒体RNA反转录获得的cDNA上,以及与线粒体RNA反转录获得的cDNA退火的一个或多个引物组,且每个引物组用于一条线粒体核酸链的靶向扩增,由此扩增一条或多条线粒体核酸链。在一个实施方案中,所述每个引物组中的一条引物是针对线粒体转录本的3’端特异性引物或5’端特异性引物,例如,表1所示的针对线粒体转录本的一条或多条3’端特异性引物;所述每个引物组中的另一条引物是通用PCR引物,用于结合线粒体RNA反转录获得的cDNA上的通用PCR引物结合位点。
在进一步的实施方案中,试剂盒还包含向扩增的一条或多条线粒体核酸添加至少一个样本索引的核苷酸序列(例如,所述添加至少一个样本索引的核苷酸序列选自至少一对SEQ ID NO: 16和17、SEQ ID NO: 18和19、SEQ ID NO: 20和21、SEQ ID NO: 22和23、SEQID NO: 24和25、SEQ ID NO: 26和27样本索引引物),用于多样本三代测序文库的构建。
试剂盒还可以包含用于产生线粒体转录组测序文库所需的一种或多种其他组分。例如,引物延伸或扩增的酶、缓冲剂等。
试剂盒的组分通常存在于合适的包装材料中。
术语“包装材料”是指用于容纳试剂盒内容物的一种或多种物理结构。包装材料通过常规方法构造,一般提供无菌的、无污染的环境。例如,保护内容物免于外部环境的合适容器如小瓶。
包装材料可具有标签,其表明可以使用产生线粒体转录组测序文库的组分。
另外,包装材料包含说明如何使用试剂盒中的材料的说明。“使用说明”通常包括描述试剂浓度或至少一种测定方法参数的有形表达,如待混合的试剂和样品的相对量、试剂/样品混合物的维持时间、温度、缓冲条件等等。
识别线粒体测序数据中线粒体变异和确定多细胞真核生物中细胞谱系的装置
本发明提供了基于二代测序和三代测序来识别单细胞线粒体测序数据中线粒体变异和确定多细胞真核生物中细胞谱系的装置,所述装置被配置成能够实现本发明的确定单细胞线粒体转录组中的变异的方法和/或确定多细胞真核生物中细胞谱系的方法,例如,所述装置被配置成:
-接收输入,所述输入包括二代测序数据和三代测序数据,例如,单细胞线粒体全长三代测序数据和Bulk线粒体转录组二代测序数据,以及任选地,单细胞转录组二代测序数据,
-对输入的各测序数据进行质控和过滤,并进行生物信息学分析。
在一些实施方案中,所述装置包含一个或多个以下模块:
模块1. 接收单细胞线粒体全长三代测序数据,并映射到靶线粒体DNA样本;
模块2. 接收Bulk线粒体转录组二代测序数据,获得高准确性的基因水平的变异和barcode/UMI序列,例如,获得高准确性单细胞变异位点(例如,SNP/InDel位点)和barcode/UMI序列;
模块3. 用模块2中获得的高准确性barcode/UMI序列校正模块1的barcode/UMI序列, 再与模块2得到的高准确性的变异位点(例如,SNP/InDel位点)进行整合,获得高准确性的线粒体变异位点(例如,SNP/InDel位点);以及
模块4. 基于模块3得到的高准确性的线粒体变异位点(例如,SNP/InDel位点)的存在对细胞进行聚类,并推断细胞的谱系;优选地,通过计算单细胞变异位点(例如,SNP/InDel位点)的变异频率,构建细胞-变异(例如,SNP/InDel变异) 矩阵,由此获得细胞的进化谱系。
任一所述模块1-4可以通过计算机执行。因此,本发明还提供经编程可以执行任一所述模块1-4的计算机。计算机典型地包括:与计算机通讯界面连接的CPU,***记忆(RAM),非暂时性存储器(ROM),和一个或多个其他存储装置如硬板、软盘、CD ROM drive。计算机还可以包括展示装置,例如打印机、CRT监测器或LCD展示器,以及输入装置,例如键盘、鼠标、笔、触摸屏或语音激活***。输入装置可以接收数据,例如通过界面直接从测序仪接收。
在一些实施方案中,所述装置包含一个或多个以下模块:
模块1. 接收单细胞线粒体全长三代测序数据,并映射到靶线粒体DNA样本;
模块2. 接收单细胞转录组二代测序数据并分析,获得高准确性barcode/UMI序列;
模块3. 接收Bulk线粒体转录组二代测序数据并分析,获得高准确性的基因水平的变异,从而获得高准确性单细胞变异位点(例如,SNP/InDel位点);
模块4. 将所述模块1和模块2中分别获得的长、短读长单细胞数据barcode/UMI映射,获得高准确性barcode/UMI序列及较高准确性单细胞变异位点(例如,SNP/InDel位点);
模块5. 对模块4得到的barcode/UMI校正后的单细胞全长转录组序列与模块3得到的高准确性的变异位点(例如,SNP/InDel位点)进行整合,获得高准确性的线粒体变异位点(例如,SNP/InDel位点);以及
模块6. 基于模块5得到的高准确性的线粒体变异位点(例如,SNP/InDel位点)的存在对细胞进行聚类,并推断细胞的谱系;优选地,通过计算单细胞变异位点(例如,SNP/InDel位点)的变异频率,构建细胞-变异(例如,SNP/InDel变异) 矩阵,由此获得细胞的进化谱系。
任一所述模块1-6可以通过计算机执行。因此,本发明还提供经编程可以执行任一所述模块1-6的计算机。计算机典型地包括:与计算机通讯界面连接的CPU,***记忆(RAM),非暂时性存储器(ROM),和一个或多个其他存储装置如硬板、软盘、CD ROM drive。计算机还可以包括展示装置,例如打印机、CRT监测器或LCD展示器,以及输入装置,例如键盘、鼠标、笔、触摸屏或语音激活***。输入装置可以接收数据,例如通过界面直接从测序仪接收。
本发明的方法、试剂盒和装置的应用
本发明的方法、试剂盒和装置能够用于单细胞线粒体相关突变的高通量检测。
常规的线粒体相关突变的检测包括:PCR-RFLP和DNA一代测序、芯片方法等,但这些方法通常集中在对少数常见的线粒体基因位点进行突变和缺失筛查,如对线粒体相关性糖尿病常见突变A3243G的筛查,但阳性检出率低,这是因为线粒体疾病具有量效现象,即小量的线粒体DNA突变往往不出现临床症状,但随着突变线粒体比例增高,才会出现临床表现,且临床严重程度可能和突变比例成正相关,而现有的常规检测方法对于单细胞的线粒体突变都不能有效地检测。本发明的方法、试剂盒和装置能够克服现有技术中存在的对线粒体相关突变的检测结果假阳性高、无法高通量检测的技术问题,故本发明对单细胞中线粒体转录组的准确检测具有非常重要的临床指导意义。
进一步地,本发明的方法、试剂盒和装置能够用于追踪细胞谱系。细胞谱系追踪在肿瘤和免疫***研究中尤其具有重要意义。通过肿瘤细胞谱系追踪可以对肿瘤克隆的性质精确定义并为寻找靶点提供信息。
实施例
描述以下实施例以辅助对本发明的理解。不意在且不应当以任何方式将实施例解释成限制本发明的保护范围。
实施例1. 细胞条形码标记的全长cDNA的制备
将来自健康个体的5mL全血收集在含肝素的采血管中。通过聚蔗糖密度梯度离心分离外周血单个核细胞(PBMC)(供应商:迈顺生物,货号:P121111505C)。计数细胞,获得约10000个单细胞的细胞稀释液。
使用10x Chromium Next GEM Single Cell 5' Kit v2试剂盒(10x Genomics公司,目录号:PN-1000263),芯片为Chromium Next GEM Chip K(10x Genomics公司,目录号:PN-1000287),按照生产商的操作说明书所述,将获得的细胞稀释液加载到Chromium NextGEM Chip K芯片进样孔内。
将加载有细胞稀释液的芯片装载到10× Chromium Controller仪器(购自10×Genomics公司)上,运行程序,形成油包水乳液微滴(GEMs,Gel Beads-in-emulsion)。在所产生的GEMs中大部分(约90–99%)不含有细胞,剩余的GEMs中大部分含有仅一个细胞,GEMs中的凝胶微珠上的寡核苷酸包含(i) 一个Illumina R1序列(结合read 1测序引物),(ii)10x条形码(10x Barcode),具有16 nt长度,一个微珠对应于一种条形码;(iii) 一个10 nt唯一分子标识符(UMI),所述UMI是一段随机序列,10个碱基长的UMI有100多万种序列的变化(410 = 1,048,576);和(iv) 13 nt模板转换寡核苷酸(TSO)。
在GEMs产生后,通过10x Chromium Next GEM Single Cell 5' Kit v2试剂盒中的含有逆转录(RT)试剂和多聚(dT)逆转录引物的Master Mix逆转录细胞mRNA,并通过模板转换和转录物延伸获得如图2所示的10x 细胞条形码标记的全长cDNA。
对GEMs进行破油和纯化,然后对经纯化的10x 细胞条形码标记的全长cDNA使用10x Chromium Next GEM Single Cell 5' Kit v2试剂盒中的引物进行全转录组扩增PCR(whole transcriptome amplification (WTA) PCR),其中所用引物结合单细胞转录组逆转录期间在各全长cDNA两端添加的共有5'和3'末端引物结合位点,正向引物为:5’-CTACACGACGCTCTTCCGATCT-3’(SEQ ID NO:28),Poly-dT反向引物为:5’- Poly-dTAAGCAGTGGTATCAACGCAGAG-3’(SEQ ID NO:29),由此获得各单细胞的全转录组扩增产物。
实施例2. 线粒体cDNA的富集
对一部分实施例1中制备的细胞条形码标记的全转录组全长cDNA进行了线粒体cDNA的富集。
对线粒体如下基因的cDNA进行了富集:MT-RNR1、MT-RNR2、MT-ND1、MT-ND2、MT-CO1、MT-CO2、MT-ATP8、MT-ATP6、MT-CO3、MT-ND3、MT-ND4L、MT-ND4、MT-ND5、MT-ND6、MT-CYTB。具体而言,设计并合成了下表1所示的15条下游反向引物。将15条下游反向引物各溶解至100μM,然后按表1所示引物终浓度条件混合在1管PCR管中,补加无核酸酶的ddH2O至总体积为100ul,即得到下游反向引物Mix;向该PCR管中添加上游引物,为10x Chromium NextGEM Single Cell 5' Kit v2试剂盒中的Read1引物序列(5’-CTACACGACGCTCTTCCGATCT-3’(SEQ ID NO:28)),并添加实施例1获得的WTA产物,实施线粒体转录组的PCR扩增。
表1. 用于特异性扩增线粒体cDNA的下游反向引物序列信息:
SEQ ID NO: |
核苷酸序列(5’- 3’) |
扩增的转录本的来源基因 |
引物终浓度uM |
1 |
CCATGTTACGACTTGTCTCCTCTATAT |
MT-RNR1 |
0.1 |
2 |
GTATAATACTAAGTTGAGATGATATCATTTACGG |
MT-RNR2 |
0.05 |
3 |
GATTGTAATGGGTATGGAGACATATCATATAAG |
MT-ND1 |
0.3 |
4 |
CGTGGTAAGGGCGATGAGT |
MT-ND2 |
0.3 |
5 |
GTTCTTCGAATGTGTGGTAGGGTG |
MT-CO1 |
0.1 |
6 |
AGATTTTTAGGGGAATTAATTCTAGGACGA |
MT-CO2 |
0.2 |
7 |
TGAAGCGAACAGATTTTCGTTCA |
MT-ATP8 |
0.5 |
8 |
CTAGAAGTGTGAAAACGTAGGCTTG |
MT-ATP6 |
0.2 |
9 |
AGACATACAGAAATAGTCAAACCACATCTA |
MT-CO3 |
0.2 |
10 |
CTCATAGGCCAGACTTAGGGCTAG |
MT-ND3 |
0.3 |
11 |
GTCTAGGCCATATGTGTTGGAGAT |
MT-ND4L |
0.4 |
12 |
GGATAGGAGGAGAATGGGGGATAG |
MT-ND4 |
0.3 |
13 |
GAGTAGGGTTAGGATGAGTGGG |
MT-ND5 |
0.3 |
14 |
TGTGGGTTTAGTAATGGGGTTTGT |
MT-ND6 |
0.5 |
15 |
TAGGGAGATAGTTGGTATTAGGATTAGGAT |
MT-CYTB |
0.3 |
扩增线粒体转录组的PCR反应以1个PCR管/1个全血个体样本的线粒体cDNA实施,具体PCR扩增程序为:
预变性:95℃ 3min
6个循环(98℃ 30s;65℃ 30s;72℃ 3min)
终延伸:72℃ 5min
保存温度:4℃。
使用1.0X SPRI磁珠(Solid Phase Reversible Immobilization,Beckman,B23317)纯化PCR扩增产物,获得富集的线粒体转录组cDNA产物。
实施例3. 单细胞线粒体转录组的三代测序
将一部分实施例2获得的线粒体转录组富集产物进行长读长的高通量三代测序。
具体地,采用ONT(Oxford Nanopore Technology公司)的长读长三代测序技术,对实施例2获得的单细胞线粒体转录组构建Nanopore文库,使用Nanopore MinION测序仪(Oxford Nanopore Technology公司)对构建的Nanopore文库进行测序,获得如下结果。
长读长三代测序线粒体全长转录组的原始下机数据reads数目为95,123,378,采用 minimap2软件(https://github.com/lh3/minimap2) 将reads比对到参考基因组GRCh38,结果表明比对到参考基因组的reads数目为94,930,502,比对率为99.08%,目标区域reads数目为50,892,494,占比53.61%;在bam文件中提取softclip5,并提取plrs。
实施例4. 线粒体转录组的Bulk二代测序
将一部分实施例2获得的线粒体转录组富集产物进行短读长的高通量Bulk二代测序。
具体地,按照NEBNext® Ultra™ II FS DNA Library Prep Kit for Illumina(NEB,目录号:E7805)产品说明书对实施例2的单细胞线粒体转录组富集产物执行快速可靠的线粒体cDNA片段化、末端修复、接头连接及Bulk二代测序文库构建。
Bulk二代测序文库构建后,使用Illumina Novaseq6000平台上机测序,获得Bulk线粒体转录组短读长的二代测序数据。
Bulk线粒体转录组二代测序的原始下机reads数目为83,164,940,使用 fastp软件,选择默认参数对该下机reads进行质控, 经质控过滤后的剩余reads数目为78,310,624, 占原始下机reads数目的94.16%。
采用STAR软件将质控后的reads序列比对到参考基因组GRCh38,随后使用GATK4Mitochondria Pipeline进行线粒体变异检测,得到高准确性的SNP/InDel,结果如表2所示。
表2. Bulk线粒体转录组短读长二代测序后的变异信息统计
项目 |
对项目的描述 |
数值 |
raw_reads |
原始reads数目 |
83,164,940 |
clean_reads |
质控后reads数目 |
78,310,624 |
primary mapped |
比对到参考基因组的reads数目 |
65,894,825 |
mean depth of variant |
目标基因变异的平均覆盖深度 |
3847 |
total_variant_count |
目标基因的变异数目 |
335 |
snp_count |
目标基因SNP数目 |
196 |
indel_count |
目标基因InDel数目 |
96 |
mean depth of variant[PASS] |
目标基因变异的平均覆盖深度[高质量] |
4110.24 |
total_variant_count[PASS] |
目标基因的变异数目[高质量] |
125 |
snp_count[PASS] |
目标基因SNP数目[高质量] |
78 |
indel_count[PASS] |
目标基因InDel数目[高质量] |
28 |
实施例5. 单细胞全转录组的短读长二代测序
对一部分实施例1中制备的细胞条形码标记的单细胞全转录组全长cDNA进行了短读长高通量二代测序。
使用10x Chromium Next GEM Single Cell 5' Kit v2试剂盒(10x Genomics公司,CG000331 Rev B)中的Gene Expression (GEX) Library Construction如图4所示制备单细胞转录组的短读长测序文库。
具体而言,将实施例1中制备的单细胞全长转录组cDNA进行片段化,然后经过末端修复、加A碱基及接头连接,即可在PCR过程中富集含有细胞条形码和UMI信息的基因表达片段。将制备的单细胞转录组二代测序文库随后使用Illumina Novaseq6000平台进行短读长二代测序,获得单细胞转录组的短读长二代测序数据。
单细胞转录组的短读长二代测序原始下机reads数目为672,118,966,使用 fastp软件,选择默认参数对该下机reads进行质控, 经质控过滤后剩余reads数目为672,118,966, 占原始下机reads数目的100%。
将该质控后的数据使用 Cell Ranger进行处理,统计结果如下表3所示。
表3. Cell Ranger 统计结果
项目 |
对项目的描述 |
数值 |
Estimated Number of Cells |
捕获的细胞数目 |
10,929 |
Mean Reads per Cell |
每个细胞的平均reads数目 |
30,749 |
Median Genes per Cell |
每个细胞中位基因数目 |
1,614 |
Number of Reads |
reads数目 |
336,059,483 |
Valid Barcodes |
有效 barcode 占比 |
89.70% |
Reads Mapped to Genome |
比对到参考基因组的比例 |
91.20% |
Fraction Reads in Cells |
细胞内 reads的比例 |
90.30% |
Total Genes Detected |
检测到的基因总数 |
22,447 |
Median UMI Counts per Cell |
每个细胞的中位 UMI数目 |
4,107 |
将该质控后的reads序列使用 STAR软件比对到参考基因组GRCh38,提取转录组reads序列比对到参考基因组的起始位置psrs,由此,提供了高准确性的单细胞条形码/UMI序列和其后接的基因表达片段序列。
实施例6. 对测序数据的生物信息学分析
本实施例对通过二代测序和三代测序获得的测序数据进行了生物信息学分析。
实施例6.1. 利用线粒体转录组Bulk二代测序数据和单细胞线粒体全长转录组三代测序数据获得高准确度 barcode/UMI、高准确度 SNP/InDel的线粒体全长转录组序列
将实施例3获得的长读长三代测序单细胞线粒体转录组数据与实施例4获得的Bulk线粒体转录组短读长二代测序数据进行整合,旨在进一步校正所述长读长三代测序单细胞线粒体转录组的序列。
具体而言,根据线粒体基因名称将实施例4中的Bulk线粒体转录组短读长二代测序数据通过Cellranger软件获得的barcode/UMI和实施例3中的长读长三代测序单细胞线粒体转录组数据的softclip5进行分组,获得下表4。
表4. Bulk线粒体转录组短读长二代测序和单细胞线粒体转录组长读长三代测序的reads数目统计
目标线粒体基因 |
Bulk线粒体转录组短读长二代测序barcode/UMI数目 |
单细胞线粒体转录组长读长三代测序softclip5数目 |
MT-ATP6 |
20,434 |
7,722,652 |
MT-ATP8 |
95,735 |
3,079,504 |
MT-CO1 |
46,962 |
5,280,040 |
MT-CO2 |
149,669 |
9,610,192 |
MT-CO3 |
96,974 |
5,546,072 |
MT-CYB |
82,807 |
4,978,790 |
MT-ND1 |
12,291 |
1,467,608 |
MT-ND2 |
31,765 |
2,756,783 |
MT-ND3 |
31,103 |
2,179,502 |
MT-ND4 |
26,954 |
2,075,683 |
MT-ND4L |
58,752 |
2,113,493 |
MT-ND5 |
14,624 |
1,130,849 |
MT-ND6 |
7,983 |
790,002 |
MT-RNR1 |
191,186 |
4,702,531 |
MT-RNR2 |
408,461 |
3,811,227 |
采用blast软件,将长读长三代测序的单细胞线粒体转录本softclip5比对到短读长二代测序的Bulk线粒体转录本barcode/UMI,选择最佳比对的barcode/UMI。在进行该最佳比对时,不考虑错配(mismatches)和空位(gaps)。
整合短读长二代测序的Bulk线粒体转录本、长读长三代测序的单细胞线粒体转录本barcode/ UMI比对的结果(mismatches_gaps)与转录本在参考基因组的起始位置(pos_diff =|psrs-plrs|),比较不同mismatches_gaps与pos_diff下的长读长测序reads数目(见表5所示),最终选择(pos_diff<30 & mismatches_gaps<=3) |(pos_diff<20 &mismatches_gaps<=4) |(pos_diff<10 & mismatches_gaps<=5)作为长读长三代测序单细胞线粒体转录组与短读长二代测序Bulk线粒体转录组barcode/UMI映射的标准。
表5. 不同比对结果与位置差异条件下的短读长Bulk线粒体转录本、长读长测序单细胞线粒体转录本barcode/UMI匹配情况
过滤条件 |
长读长测序单细胞线粒体转录本softclip5数目 |
匹配到短读长Bulk线粒体转录本barcode/UMI的长读长测序单细胞线粒体转录本reads数目 |
占比 |
mismatches_gaps<=3 |
47,751,920 |
19,563,044 |
0.4097 |
pos_diff<10 & mismatches_gaps<=3 |
47,751,920 |
18,351,810 |
0.3843 |
pos_diff<20 & mismatches_gaps<=3 |
47,751,920 |
18,967,567 |
0.3972 |
pos_diff<30 & mismatches_gaps<=3 |
47,751,920 |
19,119,653 |
0.4004 |
mismatches_gaps<=4 |
47,751,920 |
21,538,332 |
0.4510 |
pos_diff<10 & mismatches_gaps<=4 |
47,751,920 |
19,528,885 |
0.4090 |
pos_diff<20 & mismatches_gaps<=4 |
47,751,920 |
20,226,372 |
0.4236 |
pos_diff<30 & mismatches_gaps<=4 |
47,751,920 |
20,403,109 |
0.4273 |
mismatches_gaps<=5 |
47,751,920 |
22,671,477 |
0.4748 |
pos_diff<10 & mismatches_gaps<=5 |
47,751,920 |
20,228,268 |
0.4236 |
pos_diff<20 & mismatches_gaps<=5 |
47,751,920 |
20,977,926 |
0.4393 |
pos_diff<30 & mismatches_gaps<=5 |
47,751,920 |
21,175,084 |
0.4434 |
(pos_diff<30 & mismatches_gaps<=3)| (pos_diff<20 &mismatches_gaps<=4) | (pos_diff<10 & mismatches_gaps<=5) |
47,751,920 |
30,016,535 |
0.4414 |
将匹配到短读长测序Bulk线粒体转录本barcode/UMI的长读长测序单细胞线粒体转录本reads序列(即,经短读长测序Bulk线粒体barcode/UMI校正后的单细胞线粒体全长转录本序列)与实施例4获得的短读长Bulk测序线粒体转录本reads序列进行整合,并以实施例4获得的短读长Bulk测序线粒体转录本reads序列中的SNP/InDel为依据,获得包含高准确性SNP/InDel的长读长测序单细胞线粒体转录本序列。
由于携带相同barcode/UMI的不同reads来自于相同细胞的相同转录本,通过对这样的reads进行SNP/InDel的相互校正,最终获得28个高质量线粒体变异。
根据以上生物信息学分析获得的携带高准确度barcode/UMI、高准确度SNP/InDel的线粒体全长转录本序列计算单个细胞(即,barcode标记的单细胞)水平的变异频率,构建了10795个细胞-28个变异矩阵,并根据细胞间欧式距离进行层次聚类分析(图6),据此推断出细胞的进化谱系(图7,不同颜色表示不同的细胞群分类)。
实施例6.2. 利用线粒体转录组Bulk二代测序数据、单细胞转录组二代测序数据和线粒体转录组三代测序数据获得高准确度barcode/UMI、高准确度SNP/InDel的线粒体全长转录组序列
将实施例3获得的长读长测序单细胞线粒体转录组数据与实施例5获得的短读长测序单细胞转录组数据进行整合,旨在进一步校正长读长测序单细胞线粒体转录组的barcode/UMI序列。
根据线粒体基因名称将实施例5中的短读长测序单细胞转录组数据中的barcode/UMI和实施例3中的长读长测序单细胞线粒体转录组数据的softclip5进行分组,获得下表6。
表6. 短读长测序单细胞转录本与长读长测序单细胞线粒体转录本reads数目统计
目标线粒体基因 |
短读长测序单细胞转录本barcode/UMI数目 |
长读长测序单细胞线粒体转录本 softclip5数目 |
MT-ATP6 |
325,189 |
7,722,652 |
MT-ATP8 |
164,661 |
3,079,504 |
MT-CO1 |
558,951 |
5,280,040 |
MT-CO2 |
246,860 |
9,610,192 |
MT-CO3 |
188,977 |
5,546,072 |
MT-CYB |
158,675 |
4,978,790 |
MT-ND1 |
86,467 |
1,467,608 |
MT-ND2 |
72,650 |
2,756,783 |
MT-ND3 |
62,762 |
2,179,502 |
MT-ND4 |
85,837 |
2,075,683 |
MT-ND4L |
78,658 |
2,113,493 |
MT-ND5 |
88,469 |
1,130,849 |
MT-ND6 |
37,381 |
790,002 |
MT-RNR1 |
574,902 |
4,702,531 |
MT-RNR2 |
2,069,266 |
3,811,227 |
采用 blast软件,将长读长测序的单细胞线粒体转录本softclip5比对到短读长测序的单细胞转录本barcode/UMI,选择最佳比对的barcode/UMI。在进行该最佳比对时,不考虑错配(mismatches)和空位(gaps)。
整合短读长测序单细胞转录本、长读长测序单细胞线粒体转录本barcode/ UMI比对的结果(mismatches_gaps)与转录本在参考基因组的起始位置(pos_diff =|psrs-plrs|),比较不同mismatches_gaps与pos_diff下的长读长测序reads数目(见表7所示),最终选择(pos_diff<30 & mismatches_gaps<=3) |(pos_diff<20 & mismatches_gaps<=4) |(pos_diff<10 & mismatches_gaps<=5)作为长读长测序单细胞线粒体转录本与短读长测序单细胞转录本barcode/UMI映射的标准。
表7. 不同比对结果与位置差异条件下的短读长测序单细胞转录本、长读长测序单细胞线粒体转录本barcode/UMI匹配情况
过滤条件 |
长读长测序单细胞线粒体转录本softclip5数目 |
匹配到短读长测序单细胞转录本barcode/UMI的长读长测序单细胞线粒体转录本reads数目 |
占比 |
mismatches_gaps<=3 |
47,751,920 |
28,008,938 |
0.5866 |
pos_diff<10 & mismatches_gaps<=3 |
47,751,920 |
26,047,223 |
0.5455 |
pos_diff<20 & mismatches_gaps<=3 |
47,751,920 |
26,991,453 |
0.5652 |
pos_diff<30 & mismatches_gaps<=3 |
47,751,920 |
27,225,384 |
0.5701 |
mismatches_gaps<=4 |
47,751,920 |
30,725,553 |
0.6434 |
pos_diff<10 & mismatches_gaps<=4 |
47,751,920 |
27,695,975 |
0.5800 |
pos_diff<20 & mismatches_gaps<=4 |
47,751,920 |
28,755,991 |
0.6022 |
pos_diff<30 & mismatches_gaps<=4 |
47,751,920 |
29,026,108 |
0.6079 |
mismatches_gaps<=5 |
47,751,920 |
32,475,194 |
0.6801 |
pos_diff<10 & mismatches_gaps<=5 |
47,751,920 |
28,722,588 |
0.6015 |
pos_diff<20 & mismatches_gaps<=5 |
47,751,920 |
29,860,299 |
0.6253 |
pos_diff<30 & mismatches_gaps<=5 |
47,751,920 |
30,160,527 |
0.6316 |
(pos_diff<30 & mismatches_gaps<=3) | (pos_diff<20 &mismatches_gaps<=4) | (pos_diff<10 & mismatches_gaps<=5) |
47,751,920 |
30,016,535 |
0.6286 |
将匹配到短读长测序单细胞转录本barcode/UMI的长读长测序单细胞线粒体转录本reads序列(即,经短读长测序单细胞转录本barcode/UMI校正后的单细胞线粒体全长转录本序列)与实施例4获得的短读长Bulk测序线粒体转录本reads序列进行整合,并以实施例4获得的短读长Bulk测序线粒体转录本reads序列中的SNP/InDel为依据,获得包含高准确性 SNP/InDel的长读长测序单细胞线粒体转录本序列。
由于携带相同barcode/UMI的不同reads来自于相同细胞的相同转录本,通过对这样的reads进行SNP/InDel的相互校正,最终获得63个高质量线粒体变异。
根据以上生物信息学分析获得的携带高准确度 barcode/UMI、高准确度 SNP/InDel的线粒体全长转录本序列计算单个细胞(即,barcode标记的单细胞)水平的变异频率,构建了10871个细胞-63个变异矩阵,并根据细胞间欧式距离,进行层次聚类分析(图8),据此推断出细胞的进化谱系(图9,不同颜色表示不同的细胞群分类)。
实施例7. 人源肝细胞肝癌的细胞谱系的确定
将来自人源肝细胞肝癌的细胞(中美冠科 Liver-Cancer,批号LI0050)制备单细胞悬液,计数细胞。除了所使用的细胞样本不同之外,与实施例1同样地制备细胞条形码标记的全长cDNA;与实施例2同样地进行线粒体cDNA的富集。
实施例7.1. 混合细胞样本长读长测序线粒体全长转录组
采用ONT(Oxford Nanopore Technology公司)三代长读长的测序技术,使用Nanopore MinION测序仪(Oxford Nanopore Technology公司)对本实施例的肝细胞肝癌的经富集的单细胞线粒体全长转录组进行测序。
将本实施例上述富集的单细胞线粒体全长转录组的肝癌样本(即,目标肝细胞肝癌样本)与自两例人源PBMC样本(迈顺生物,货号:P121111505C(注:该公司此货号仅表示提供的是人源PBMC样本,所述人源PBMC样本可以为自不同的人类个体(donor ID)获得的))和三例人源肝细胞肝癌样本(中美冠科,批号:LI0050)(注:该公司此货号仅表示提供的是人源肝细胞肝癌样本,所述人源肝细胞肝癌样本可以为自不同的人类个体(donor ID)获得的)分别富集的单细胞线粒体全长转录组血液和肝癌样本(共6个样本)混合在1个PCR管中,添加下述样本索引(sample index)引物,进行进一步的PCR扩增并建库,由此实施混合样本三代测序。样本索引的设计原理:首先,要满足引物设计的基本要求,碱基平衡,无发卡结构,无串联重复序列等,其次,在序列上需要满足与线粒体基因无互补序列,避免造成非特异性扩增。使用了下述样本索引(sample index)正向(F)和反向(R)引物:
SIF-1:ATAGTGATACTGACCACCGAGATCTACACATATACGCAGTCGACAACTTTCTTGCGCCATGGGACTACACGACGCTC (SEQ ID NO: 16);
SIF-2:ATAGTGATACTGACCACCGAGATCTACACCGCACGACGTACAAACGGAATCGAGCGCCATGGGACTACACGACGCTC (SEQ ID NO: 17);
SIF-3:ATAGTGATACTGACCACCGAGATCTACACGGGACAGAAGGGGACACAAGACTCGCGCCATGGGACTACACGACGCTC (SEQ ID NO: 18);
SIF-4:ATAGTGATACTGACCACCGAGATCTACACTTAACAAGTATGATAGAATCCGAAGCGCCATGGGACTACACGACGCTC (SEQ ID NO: 19);
SIF-5:ATAGTGATACTGACCACCGAGATCTACACACTCCTGGGTAAACCCTGGACAAGGCGCCATGGGACTACACGACGCTC (SEQ ID NO: 20);
SIF-6:ATAGTGATACTGACCACCGAGATCTACACCAGAAAACCCCGCCTTTGCGAGAAGCGCCATGGGACTACACGACGCTC (SEQ ID NO:21);
SIR-1:CATCAGCGAACCGGCATACGAGATGCACGAGCAGTCCCACGGTAACACTGTTGAACTCCAGTTCAGACGTGT (SEQ ID NO: 22);
SIR-2:CATCAGCGAACCGGCATACGAGATTATTATGCTTCTTTGTTCCCTGAATGTTGAACTCCAGTTCAGACGTGT (SEQ ID NO: 23);
SIR-3:CATCAGCGAACCGGCATACGAGATAGCGACTGTATGCTGTGCCTAGTTTGTTGAACTCCAGTTCAGACGTGT (SEQ ID NO: 24);
SIR-4:CATCAGCGAACCGGCATACGAGATCCTTGCGACAGGGTTTCAACGCTTTGTTGAACTCCAGTTCAGACGTGT (SEQ ID NO: 25);
SIR-5:CATCAGCGAACCGGCATACGAGATGTTAGATTGCACGATAGATGAAACTGTTGAACTCCAGTTCAGACGTGT (SEQ ID NO: 26);
SIR-6:CATCAGCGAACCGGCATACGAGATTAAATACCATCTTCTTTCTACCTGTGTTGAACTCCAGTTCAGACGTGT (SEQ ID NO: 27);
PCR反应体系
试剂 |
体积 |
终浓度 |
2X KAPA HiFi HotStart ReadyMix |
25 ul |
1X |
SIF各引物(分别为10uM) |
5 ul |
1uM |
SIR各引物(分别为10uM) |
5 ul |
1uM |
DNA模板 (为线粒体富集产物) |
15 µL |
|
总体积 |
50 µL |
|
PCR反应程序:
预变性:98℃ 3min
8个循环(98℃ 20s;54℃ 30s;72℃ 3min)
终延伸:72℃ 5min
保存温度:4℃。
将PCR产物使用0.8X SPRI磁珠进行纯化。使用Nanopore MinION测序仪(OxfordNanopore Technology公司)进行测序。获得如下结果:长读长测序混合样本线粒体转录组的原始下机数据reads数目为100,532,000。采用Smith-Waterman局部比对算法将样本index序列与read两端 300bp的序列进行比对,允许 4 个碱基错配,将read分配到具有最佳比对得分的样本,从而将本实施例中的目标肝细胞肝癌样本从混合的6个样本中拆分得到reads数目为16,375,154(占混合样本reads数目比为16.29%)。
采用minimap2软件(https://github.com/lh3/minimap2) 将拆分完成后的reads比对到参考基因组GRCh38,结果表明:比对到参考基因组的reads数目为16,317,110,比对率为99.64%,目标区域reads数目为16,212,816 ,占比99.36%。
提取bam文件的softclip5序列,并提取plrs。
实施例7.2. 短读长Bulk线粒体转录组二代测序
按照NEBNext® Ultra™ II FS DNA Library Prep Kit for Illumina(NEB,目录号:E7805)产品说明书对本实施例的肝细胞肝癌的经富集的单细胞线粒体全长转录组执行快速可靠的线粒体cDNA片段化、末端修复、接头连接及Bulk二代测序文库构建。
Bulk二代测序文库构建后,使用Illumina Novaseq6000平台上机测序,获得Bulk线粒体转录组短读长的二代测序数据。
Bulk线粒体转录组二代测序原始下机reads数目为66,946,990, 使用 fastp软件,选择默认参数对该下机reads进行质控, 经质控过滤后的剩余reads数目为59,857,212, 占原始下机reads数目的89.41%。
采用STAR软件将质控后的reads序列比对到参考基因组GRCh38,随后使用GATK4Mitochondria Pipeline进行线粒体变异检测,得到高准确性的SNP/InDel,结果如表8所示。
表8. Bulk线粒体转录组短读长二代测序后的变异信息统计
项目 |
对项目的描述 |
数值 |
raw_reads |
原始reads数目 |
66,946,990 |
clean_reads |
质控后reads数目 |
59,857,212 |
primary mapped |
比对到参考基因组的reads数目 |
53,204,720 |
mean depth of variant |
目标基因变异的平均覆盖深度 |
4,958 |
total_variant_count |
目标基因的变异数目 |
321 |
snp_count |
目标基因SNP数目 |
175 |
indel_count |
目标基因InDel数目 |
115 |
mean depth of variant[PASS] |
目标基因变异的平均覆盖深度[高质量] |
4832.144 |
total_variant_count[PASS] |
目标基因的变异数目[高质量] |
104 |
snp_count[PASS] |
目标基因SNP数目[高质量] |
43 |
indel_count[PASS] |
目标基因InDel数目[高质量] |
40 |
实施例7.3. 单细胞全转录组的短读长二代测序
对本实施例制备的细胞条形码标记的单细胞全转录组全长cDNA进行了短读长高通量二代测序。
使用10x Chromium Next GEM Single Cell 5' Kit v2试剂盒(10x Genomics公司,CG000331 Rev B)中的Gene Expression (GEX) Library Construction如图4所示制备单细胞转录组文库。
具体而言,将本实施例制备的单细胞全长转录组cDNA进行片段化,然后经过末端修复、加A碱基及接头连接,即可在PCR过程中富集含有细胞条形码和UMI信息的基因表达片段。将制备的单细胞转录组二代测序文库随后使用Illumina Novaseq6000平台进行短读长二代测序,获得单细胞转录组的短读长二代测序数据。
单细胞转录组短读长二代测序的原始下机reads数目为423,788,290,使用 fastp软件,选择默认参数对该下机reads进行质控, 经质控过滤后剩余reads数目为423,788,290, 占原始下机reads数目的100%。
将该质控后的数据使用 Cell Ranger进行处理,统计结果如下表9所示。
表9. Cell Ranger 统计结果
项目 |
对项目的描述 |
数值 |
Estimated Number of Cells |
捕获的细胞数目 |
18,187 |
Mean Reads per Cell |
每个细胞的平均reads数目 |
11,651 |
Median Genes per Cell |
每个细胞中位基因数目 |
929 |
Number of Reads |
reads数目 |
211,894,145 |
Valid Barcodes |
有效 barcode 占比 |
91.70% |
Reads Mapped to Genome |
比对到参考基因组的比例 |
78.6% |
Fraction Reads in Cells |
细胞内 reads的比例 |
73.80% |
Total Genes Detected |
检测到的基因总数 |
24,037 |
Median UMI Counts per Cell |
每个细胞的中位 UMI数目 |
2,001 |
将该质控后的序列使用 STAR软件(https://github.com/alexdobin/STAR)比对到参考基因组GRCh38,提取转录本比对到参考基因组的起始位置psrs,由此,提供了高准确性的单细胞条形码/UMI序列和其后接的基因表达片段序列。
实施例7.4. 对测序数据的生物信息学分析
将实施例7.1获得的长读长三代测序单细胞线粒体转录组数据与实施例7.3获得的短读长二代测序单细胞转录组数据进行整合,旨在进一步校正长读长测序单细胞线粒体转录组的barcode/UMI序列。
根据线粒体基因名称将实施例7.3中的短读长二代测序单细胞转录组数据中的barcode/UMI和实施例7.1中的长读长三代测序单细胞线粒体转录组数据的softclip5进行分组,获得下表10。
表10. 短读长测序单细胞转录本与长读长测序单细胞线粒体转录本reads数目统计
目标线粒体基因 |
短读长测序单细胞转录本barcode/UMI数目 |
长读长测序单细胞线粒体转录本 softclip5数目 |
MT-ATP6 |
325,189 |
1,862,460 |
MT-ATP8 |
164,661 |
1,623,516 |
MT-CO1 |
558,951 |
1,628,099 |
MT-CO2 |
246,860 |
824,904 |
MT-CO3 |
188,977 |
1,540,157 |
MT-CYB |
158,675 |
742,738 |
MT-ND1 |
86,467 |
195,213 |
MT-ND2 |
72,650 |
1,135,344 |
MT-ND3 |
62,762 |
3,106,842 |
MT-ND4 |
85,837 |
1,329,685 |
MT-ND4L |
78,658 |
1,687,007 |
MT-ND5 |
88,469 |
121,260 |
MT-ND6 |
37,381 |
168,240 |
MT-RNR1 |
574,902 |
22,065 |
MT-RNR2 |
2,069,266 |
44,434 |
采用 blast软件,将长读长测序单细胞线粒体转录本softclip5比对到短读长测序单细胞转录本barcode/UMI,选择最佳比对的barcode/UMI。在进行该最佳比对时,不考虑错配(mismatches)和空位(gaps)。
整合短读长测序单细胞转录本、长读长测序单细胞线粒体转录本barcode/ UMI比对的结果(mismatches_gaps)与转录本在参考基因组的起始位置(pos_diff =|psrs-plrs|), 比较不同mismatches_gaps与pos_diff下的长读长测序单细胞线粒体转录本reads数目(表11所示),最终选择(pos_diff<30 & mismatches_gaps<=3) |(pos_diff<20 &mismatches_gaps<=4) |(pos_diff<10 & mismatches_gaps<=5)作为长读长测序单细胞线粒体转录本与短读长测序单细胞转录本barcode/UMI 映射的标准。
表11. 不同比对结果与位置差异条件下的短读长测序单细胞转录本、长读长测序单细胞线粒体转录本barcode/UMI 匹配情况
过滤条件 |
长读长测序单细胞线粒体转录本 softclip5数目 |
匹配到短读长测序单细胞转录本barcode/UMI的长读长测序单细胞线粒体转录本reads数目 |
占比 |
mismatches_gaps<=3 |
15,184,700 |
2,819,913 |
0.1857 |
pos_diff<10 & mismatches_gaps<=3 |
15,184,700 |
2,240,071 |
0.1475 |
pos_diff<20 & mismatches_gaps<=3 |
15,184,700 |
2,320,856 |
0.1528 |
pos_diff<30 & mismatches_gaps<=3 |
15,184,700 |
2,347,948 |
0.1546 |
mismatches_gaps<=4 |
15,184,700 |
3,620,135 |
0.2384 |
pos_diff<10 & mismatches_gaps<=4 |
15,184,700 |
2,641,501 |
0.1740 |
pos_diff<20 & mismatches_gaps<=4 |
15,184,700 |
2,741,979 |
0.1806 |
pos_diff<30 & mismatches_gaps<=4 |
15,184,700 |
2,782,715 |
0.1833 |
mismatches_gaps<=5 |
15,184,700 |
4,713,299 |
0.3104 |
pos_diff<10 & mismatches_gaps<=5 |
15,184,700 |
3,149,672 |
0.2074 |
pos_diff<20 & mismatches_gaps<=5 |
15,184,700 |
3,277,754 |
0.2159 |
pos_diff<30 & mismatches_gaps<=5 |
15,184,700 |
3,339,542 |
0.2199 |
(pos_diff<30 & mismatches_gaps<=3) | (pos_diff<20 & mismatches_gaps<=4) | (pos_diff<10 &mismatches_gaps<=5) |
15,184,700 |
3,277,242 |
0.2158 |
将匹配到短读长测序单细胞转录本barcode/UMI的长读长测序单细胞线粒体转录本reads序列(即,经短读长测序单细胞转录本barcode/UMI校正后的单细胞线粒体全长转录本序列)与实施例7.2获得的短读长Bulk测序线粒体转录本reads序列进行整合,并以实施例7.2获得的短读长Bulk测序线粒体转录本reads序列中的SNP/InDel为依据,获得包含高准确性 SNP/InDel的长读长测序单细胞线粒体转录本序列。
由于携带相同barcode/UMI的不同reads来自于相同细胞的相同转录本,通过对这样的reads进行SNP/InDel的相互校正,最终获得37个高质量线粒体变异。
根据以上生物信息学分析获得的携带高准确度 barcode/UMI、高准确度 SNP/InDel的线粒体全长转录本序列计算单个细胞(即,barcode标记的单细胞)水平的变异频率,构建了17700个细胞-37个变异的矩阵,并根据细胞间欧式距离,进行层次聚类分析(图10),据此推断出细胞的进化谱系(图11,不同颜色表示不同的细胞群分类)。
实施例8. 使用商业化细胞系检验本发明的方法
在本实施例中,使用3种具有已知线粒体基因突变的细胞系验证了本发明的方法可以有效识别线粒体基因突变。与实施例1-6类似地,对3种细胞系(HepG2(供应商:ATCC,货号:HB-8065)、THP-1(供应商:ATCC,货号:TIB-202)、AGS(供应商:ATCC,货号:CRL-1739))实施了单细胞线粒体全长转录组三代测序、单细胞全转录组二代测序和Bulk线粒体转录组二代测序,并对获得的测序数据进行了分析。
实施例8.1. 单细胞线粒体全长转录组三代测序
与实施例3类似地,对3种细胞系(HepG2、THP-1、AGS)实施了单细胞线粒体全长转录组三代测序。该长读长三代测序线粒体转录组的原始下机数据reads数目为63,334,135。采用 minimap2软件(https://github.com/lh3/minimap2) 将reads比对到参考基因组GRCh38,结果表明:比对到参考基因组的reads数目为56,414,636,比对率为89.07%,目标区域reads数目为38,774,532,占比68.73%。
提取bam 文件的softclip5序列,并提取plrs。
实施例8.2. 短读长Bulk线粒体转录组二代测序
与实施例4类似地,按照NEBNext® Ultra™ II FS DNA Library Prep Kit forIllumina(NEB,目录号:E7805)产品说明书对3种细胞系(HepG2、THP-1、AGS)的经富集的单细胞线粒体全长转录组执行快速可靠的线粒体cDNA片段化、末端修复、接头连接及Bulk二代测序文库构建,并测序。
Bulk线粒体转录组二代测序原始下机reads数目为184,685,770, 使用 fastp软件,选择默认参数对该下机reads进行质控, 过滤后剩余reads数目为184,589,764, 占原始下机reads数目的99.94%。
采用STAR软件将质控后的reads序列比对到参考基因组GRCh38,随后使用GATK4Mitochondria Pipeline进行线粒体变异检测,得到11个高准确性的SNP(表12),用于后续分析。
表12. Bulk线粒体转录组短读长二代测序后的变异信息统计
细胞系名称 |
疾病名称 |
基因 |
基因组改变 |
覆盖深度 |
变异频率 |
AGS |
胃癌 |
MT-ATP6 |
g.chrM:9055G>A |
5125 |
0.233 |
AGS |
胃癌 |
MT-ND4L |
g.chrM:10550A>G |
2261 |
0.238 |
AGS |
胃癌 |
MT-CO3 |
g.chrM:9698T>C |
3019 |
0.298 |
AGS |
胃癌 |
MT-ND2 |
g.chrM:4561T>C |
2803 |
0.232 |
HepG2 |
肝癌 |
MT-CO3 |
g.chrM:9950T>C |
110 |
0.607 |
HepG2 |
肝癌 |
MT-ND3 |
g.chrM:10373G>A |
1116 |
0.114 |
HepG2 |
肝癌 |
MT-CYB |
g.chrM:14757T>C |
1015 |
0.602 |
THP-1 |
白血病 |
MT-ND1 |
g.chrM:3970C>T |
5455 |
0.199 |
THP-1 |
白血病 |
MT-ND2 |
g.chrM:4732A>G |
4062 |
0.087 |
THP-1 |
白血病 |
MT-ND5 |
g.chrM:13928G>C |
1160 |
0.176 |
THP-1 |
白血病 |
MT-CO1 |
g.chrM:6962G>A |
2565 |
0.208 |
实施例8.3. 单细胞全转录组的短读长二代测序
与实施例5类似地,对3种细胞系(HepG2、THP-1、AGS)分别制备了细胞条形码标记的单细胞全转录组全长cDNA,并进行了短读长高通量二代测序。
单细胞转录组短读长二代测序的原始下机reads数目为483,079,508,使用 fastp软件,选择默认参数对该下机reads进行质控, 经质控过滤后剩余reads数目为483,079,508, 占原始下机reads数目的100%。
将该质控后的数据使用 Cell Ranger进行处理,统计结果如下表13所示。
表13. Cell Ranger统计结果
项目 |
对项目的描述 |
数值 |
Estimated Number of Cells |
捕获的细胞数目 |
16,670 |
Mean Reads per Cell |
每个细胞的平均reads数目 |
14,489 |
Median Genes per Cell |
每个细胞中位基因数目 |
1,164 |
Number of Reads |
reads数目 |
241,539,754 |
Valid Barcodes |
有效 barcode 占比 |
92.30% |
Reads Mapped to Genome |
比对到参考基因组的比例 |
93.1% |
Fraction Reads in Cells |
细胞内 reads的比例 |
97.20% |
Total Genes Detected |
检测到的基因总数 |
25,401 |
Median UMI Counts per Cell |
每个细胞的中位 UMI数目 |
1,833 |
将该质控后的序列使用 STAR软件(https://github.com/alexdobin/STAR)比对到参考基因组GRCh38,提取转录本比对到参考基因组的起始位置psrs,由此,提供了高准确性的单细胞条形码/UMI序列和其后接的基因表达片段序列。
实施例8.4. 对测序数据的生物信息学分析
将实施例8.1获得的长读长三代测序单细胞线粒体转录组数据与实施例8.3获得的短读长二代测序单细胞转录组数据进行整合,旨在进一步校正长读长测序单细胞线粒体转录组的barcode/UMI序列。
根据线粒体基因名称将实施例8.3中的短读长二代测序单细胞转录组数据中的barcode/UMI和实施例8.1中的长读长三代测序单细胞线粒体转录组数据的softclip5进行分组,获得下表14。
表14. 短读长测序单细胞转录本与长读长测序单细胞线粒体转录本reads数目统计
目标线粒体基因 |
短读长测序单细胞转录本barcode/UMI数目 |
长读长测序单细胞线粒体转录本 softclip5数目 |
MT-ATP6 |
1,151,299 |
1,427,969 |
MT-ATP8 |
505,344 |
2,133,519 |
MT-CO1 |
1,404,862 |
1,689,369 |
MT-CO2 |
1,127,308 |
4,319,140 |
MT-CO3 |
723,077 |
3,714,381 |
MT-CYB |
388,633 |
2,068,604 |
MT-ND1 |
333,244 |
2,238,233 |
MT-ND2 |
196,785 |
1,367,499 |
MT-ND3 |
129,169 |
1,210,243 |
MT-ND4 |
310,451 |
1,682,299 |
MT-ND4L |
263,600 |
1,466,838 |
MT-ND5 |
230,940 |
683,674 |
MT-ND6 |
124,823 |
511,105 |
MT-RNR1 |
2,089,213 |
5,099,172 |
MT-RNR2 |
5,968,084 |
5,504,621 |
采用 blast软件,将长读长测序单细胞线粒体转录本softclip5比对到短读长测序单细胞转录本barcode/UMI,选择最佳比对的barcode/UMI。在进行该最佳比对时,不考虑错配(mismatches)和空位(gaps)。
整合短读长测序单细胞转录本、长读长测序单细胞线粒体转录本barcode/ UMI比对的结果(mismatches_gaps)与转录本在参考基因组的起始位置(pos_diff =|psrs-plrs|), 比较不同mismatches_gaps与pos_diff下的长读长测序单细胞线粒体转录本reads数目,最终选择(pos_diff<30 & mismatches_gaps<=3) |(pos_diff<20 & mismatches_gaps<=4) |(pos_diff<10 & mismatches_gaps<=5)作为长读长测序单细胞线粒体转录本与短读长测序单细胞转录本barcode/UMI 映射的标准,获得了8,219,127 reads, 占比21.19%。
由于携带相同barcode/UMI的不同reads来自于相同细胞的相同转录本,通过对这样的reads进行SNP/InDel的相互校正,最终获得7个高质量线粒体变异。
根据以上生物信息学分析获得的携带高准确度 barcode/UMI、高准确度 SNP/InDel的线粒体全长转录本序列计算单个细胞(即,barcode标记的单细胞)水平的变异频率,构建了16,602个细胞-7个变异的矩阵,并根据细胞间欧式距离,进行层次聚类分析,获得如图12所示的细胞-变异频率层次聚类图。
将实施例 8.3短读长二代测序单细胞转录组表达数据进行聚类分析(分析方法参考https://scanpy-tutorials.readthedocs.io/en/latest/pbmc3k.html),并将HepG2细胞系特异性表达的TF基因、THP-1细胞系特异性表达的AZU1基因、AGS细胞系特异性表达的PRSS23基因作为标志基因,对这三种细胞系进行细胞群标记,结果如图13所示,图中所示“Mix”亚群表示因为细胞周期等原因造成的一部分细胞无法通过基因表达水平进行准确分类。将谱系结果映射到短读长二代测序单细胞转录本表达数据,得到基于谱系聚类的细胞亚群分类结果(图14)。
将基于短读长二代测序单细胞转录组表达的细胞亚群结果作为标准,基于谱系聚类的HepG2、THP-1、AGS细胞系的召回率(recall)分别为99.48%、99.37%、98.62%,由此,本发明的方法可以有效识别各细胞中的线粒体基因突变。
以上描述了本发明的示例性实施方案,本领域技术人员应当理解的是,这些公开内容仅是示例性的,在本发明的范围内可以进行各种其它替换、适应和修改。因此,本发明不限于文中列举的具体实施方案。
序列表
<110> 北京百图智检科技服务有限公司
<120> 确定细胞谱系的方法、试剂盒和装置
<130>
<160> 28
<170> PatentIn版本3.3
<210> 1
<211> 27
<212> DNA
<213> 人工序列
<220>
<223> 引物
<400> 1
ccatgttacg acttgtctcc tctatat 27
<210> 2
<211> 34
<212> DNA
<213> 人工序列
<220>
<223> 引物
<400> 2
gtataatact aagttgagat gatatcattt acgg 34
<210> 3
<211> 33
<212> DNA
<213> 人工序列
<220>
<223> 引物
<400> 3
gattgtaatg ggtatggaga catatcatat aag 33
<210> 4
<211> 19
<212> DNA
<213> 人工序列
<220>
<223> 引物
<400> 4
cgtggtaagg gcgatgagt 19
<210> 5
<211> 24
<212> DNA
<213> 人工序列
<220>
<223> 引物
<400> 5
gttcttcgaa tgtgtggtag ggtg 24
<210> 6
<211> 30
<212> DNA
<213> 人工序列
<220>
<223> 引物
<400> 6
agatttttag gggaattaat tctaggacga 30
<210> 7
<211> 23
<212> DNA
<213> 人工序列
<220>
<223> 引物
<400> 7
tgaagcgaac agattttcgt tca 23
<210> 8
<211> 25
<212> DNA
<213> 人工序列
<220>
<223> 引物
<400> 8
ctagaagtgt gaaaacgtag gcttg 25
<210> 9
<211> 30
<212> DNA
<213> 人工序列
<220>
<223> 引物
<400> 9
agacatacag aaatagtcaa accacatcta 30
<210> 10
<211> 24
<212> DNA
<213> 人工序列
<220>
<223> 引物
<400> 10
ctcataggcc agacttaggg ctag 24
<210> 11
<211> 24
<212> DNA
<213> 人工序列
<220>
<223> 引物
<400> 11
gtctaggcca tatgtgttgg agat 24
<210> 12
<211> 24
<212> DNA
<213> 人工序列
<220>
<223> 引物
<400> 12
ggataggagg agaatggggg atag 24
<210> 13
<211> 22
<212> DNA
<213> 人工序列
<220>
<223> 引物
<400> 13
gagtagggtt aggatgagtg gg 22
<210> 14
<211> 24
<212> DNA
<213> 人工序列
<220>
<223> 引物
<400> 14
tgtgggttta gtaatggggt ttgt 24
<210> 15
<211> 30
<212> DNA
<213> 人工序列
<220>
<223> 引物
<400> 15
tagggagata gttggtatta ggattaggat 30
<210> 16
<211> 77
<212> DNA
<213> 人工序列
<220>
<223> 引物
<400> 16
atagtgatac tgaccaccga gatctacaca tatacgcagt cgacaacttt cttgcgccat 60
gggactacac gacgctc 77
<210> 17
<211> 77
<212> DNA
<213> 人工序列
<220>
<223> 引物
<400> 17
atagtgatac tgaccaccga gatctacacc gcacgacgta caaacggaat cgagcgccat 60
gggactacac gacgctc 77
<210> 18
<211> 77
<212> DNA
<213> 人工序列
<220>
<223> 引物
<400> 18
atagtgatac tgaccaccga gatctacacg ggacagaagg ggacacaaga ctcgcgccat 60
gggactacac gacgctc 77
<210> 19
<211> 77
<212> DNA
<213> 人工序列
<220>
<223> 引物
<400> 19
atagtgatac tgaccaccga gatctacact taacaagtat gatagaatcc gaagcgccat 60
gggactacac gacgctc 77
<210> 20
<211> 77
<212> DNA
<213> 人工序列
<220>
<223> 引物
<400> 20
atagtgatac tgaccaccga gatctacaca ctcctgggta aaccctggac aaggcgccat 60
gggactacac gacgctc 77
<210> 21
<211> 77
<212> DNA
<213> 人工序列
<220>
<223> 引物
<400> 21
atagtgatac tgaccaccga gatctacacc agaaaacccc gcctttgcga gaagcgccat 60
gggactacac gacgctc 77
<210> 22
<211> 72
<212> DNA
<213> 人工序列
<220>
<223> 引物
<400> 22
catcagcgaa ccggcatacg agatgcacga gcagtcccac ggtaacactg ttgaactcca 60
gttcagacgt gt 72
<210> 23
<211> 72
<212> DNA
<213> 人工序列
<220>
<223> 引物
<400> 23
catcagcgaa ccggcatacg agattattat gcttctttgt tccctgaatg ttgaactcca 60
gttcagacgt gt 72
<210> 24
<211> 72
<212> DNA
<213> 人工序列
<220>
<223> 引物
<400> 24
catcagcgaa ccggcatacg agatagcgac tgtatgctgt gcctagtttg ttgaactcca 60
gttcagacgt gt 72
<210> 25
<211> 72
<212> DNA
<213> 人工序列
<220>
<223> 引物
<400> 25
catcagcgaa ccggcatacg agatccttgc gacagggttt caacgctttg ttgaactcca 60
gttcagacgt gt 72
<210> 26
<211> 72
<212> DNA
<213> 人工序列
<220>
<223> 引物
<400> 26
catcagcgaa ccggcatacg agatgttaga ttgcacgata gatgaaactg ttgaactcca 60
gttcagacgt gt 72
<210> 27
<211> 72
<212> DNA
<213> 人工序列
<220>
<223> 引物
<400> 27
catcagcgaa ccggcatacg agatgttaga ttgcacgata gatgaaactg ttgaactcca 60
gttcagacgt gt 72
<210> 28
<211> 22
<212> DNA
<213> 人工序列
<220>
<223> 引物
<400> 28
ctacacgacg ctcttccgat ct 22