一种基于Nanopore宏基因组RNA-seq的生物信息学检测病原体的方法

技术领域

本发明属于生物信息学领域，具体涉及一种基于Nanopore宏基因组RNA-seq的生物信息学检测病原体的方法。

背景技术

多种病原体(如细菌、真菌、病毒、非典型病原体、寄生虫)可引起人器官感染。临床微生物实验室鉴定病原体的主要方法是培养，但病毒、非典型病原体和寄生虫无法常规鉴定。PCR虽然可以鉴定多种病原体，但仍会忽略一些未知的病原体。宏基因组测序(mNGS)可高通量鉴定样本中的细菌、真菌、病毒、非典型病原体、寄生虫以及新型病原体，弥补了传统微生物检测方法的局限性。

RNA测序比DNA测序更有益。一方面，如果仅提取DNA，则RNA病毒无法检测到。另一方面，进行总RNA提取/测序可捕获DNA和RNA的表达，并且mRNA序列可翻译成蛋白质。氨基酸序列比核苷酸序列更保守，因此可产生更明确的分类信息。

Nanopore单分子实时测序技术也被称为第三代纳米孔测序技术。作为一个新型的测序平台，它具有低成本、高通量、非标记和测序长度长等优势。不同于其他测序技术，Nanopore单分子实时测序技术不需要将基因组提前打断成数百万个几百个碱基长的片段，便能够快速、经济地生产出长达上万个碱基的读长数据。Nanopore单分子实时测序技术读长长、一致性准确度高、单分子实时检测能够突破二代测序技术读长短，存在GC bias等局限。

发明内容

本发明的目的是检测临床样本含有的病原体种类、丰度和/或基因组覆盖度。

本发明首先保护一种高通量检测若干待测样本含有的病原体种类、丰度和/或基因组覆盖度的方法，依次包括如下步骤：

(1)分别获得待测样本的cDNA；

(2)将各个待测样本的cDNA进行高通量测序，得到高通量测序结果；

(3)生物信息学分析，获得若干待测样本含有的病原体种类、丰度和/或基因组覆盖度；

所述步骤(2)中，将各个待测样本的cDNA进行高通量测序的方法包括(2-1)—(2-5)：

(2-1)分别取待测样本的cDNA，进行PCR扩增，得到PCR扩增产物；

(2-2)分别取PCR扩增产物，纯化，得到纯化产物；

(2-3)分别取纯化产物，进行末端修复，得到末段修饰后的cDNA；

(2-4)分别取末段修饰后的cDNA，连接barcode，得到连接barcode的样品；

(2-5)将各个连接barcode的样品混合，高通量测序；

所述步骤(3)中，生物信息学分析，获得若干待测样本含有的病原体种类、丰度和/或基因组覆盖度的方法包括(3-1)—(3-5)：

(3-1)Basecalling；

(3-2)拆分；

(3-3)质控；

(3-4)去宿主；

(3-5)物种比对；

(3-6)与参考基因组比对，确定病原体，获得若干待测样本含有的病原体种类、丰度和/或基因组覆盖度。

所述步骤(1)中，获得待测样本的cDNA的方法可为：将待测样本的RNA进行退火，得到退火产物；之后进行逆转录和二链合成。

所述进行退火时的引物(即引物KN8)的核苷酸序列如SEQ ID NO：1所示。

所述待测样本的RNA可为采用qiagenAllPrepPowerViral DNA/RNA Kit(Catalogno.28000-50)提取待测样本获得。

所述退火反应体系可为13μl，包括1μl浓度为10μM的引物KN8溶液、1μl dNTP(浓度为10mM)和11μl待测样本的RNA。

所述退火程序具体可为：65℃5min；冰上急冷>1min；热盖105℃。

所述逆转录反应体系可为20μl，包括13μl退火产物、4μl super IVase buffer(Thermo Fisher)、1μl 100mM DDT、1μl RNase 007inhibitor和1μl super IVase(ThermoFisher)。

所述逆转录程序可为：42℃50min；70℃10min；4℃∞。

所述二链合成可为向逆转录产物中加入1μl lemon酶，之后混匀，进行二链合成，得到cDNA。

所述二链合成的反应程序可为：37℃3min；75℃15min；4℃保存。

所述(2-1)中，进行PCR扩增的引物(即primer K)的核苷酸序列可如SEQ ID NO：2所示。

所述(2-1)中，进行PCR扩增的反应体系可为50μl，包括5μl cDNA、25μl 2×mix(Q5溶液)、2μl primer K和18μl无核酸水。所述进行PCR扩增的反应程序可为：98℃30s；98℃10s，55℃15s，72℃1min，35cycles；72℃10min。

所述(2-2)中，纯化可采用TakaRaMiniBEST DNA Fragment Purification kit(cat#9761lot#AH70942A)进行。

所述(2-3)中，进行末端修复的反应体系可为15μl，包括Xμl纯化产物(含50ngDNA)、(12.5-X)μl无核酸水、1.75μl ultraⅡEnd-prep reaction buffer和0.75μl ultraⅡEnd-prep enzyme mix。进行末端修复的程序可为：20℃3min；65℃5min。

所述(2-4)中，连接barcode的方法可为：(1)制备反应体系。反应体系为20μl，包括4μl无核酸水、3μl末段修饰后的cDNA、2.5μl barcode、10μl Master mix和0.5μl Enhancer(目的为将barcode加在DNA末端)。(2)取所述反应体系，反应，得到连接barcode的样品。反应程序为：20℃20min；65℃10min。

所述(2-5)中，将各个连接barcode的样品混合后、高通量测序前，还可包括纯化的步骤。其中混合和纯化的步骤具体可为：

(1)将所有连接barcode的样品放入一个1.5ml的LoBind tube中，然后加入磁珠(总cDNA体积：磁珠体积＝1:0.4)，充分混合；

(2)室温孵育10min后，放到磁力架上，直至澄清；

(3)去除上清；

(4)加入500μl 80％乙醇水溶液洗涤，旋转EP管；

(5)去除上清，80％乙醇水溶液再洗一次；去除上清，瞬离，再吸干净，开盖放置30s，晾干；

(6)加入35μl无核酸水，室温放置5min，洗脱DNA；

(7)将LoBind tube管重新放回磁力架上，液体澄清后，将上清转移至新的LoBindtube管中，得到纯化后cDNA样品(已加接头)；

(8)制备反应体系。反应体系为50μl，包括Xμl纯化后cDNA样品(含30-50ngcDNA)(已加接头)、(30-X)μl无核酸水、5μl AMⅡ、10μl 5×NEB Quick Ligation buffer和5μlT

(9)取所述反应体系，21℃20min(开热盖)。

(10)纯化

(10-1)向完成步骤(9)的反应体系汇总加入20μl磁珠，室温孵育10min(10min内反复颠倒混匀)；

(10-2)移到磁力架上，直至澄清，弃上清；

(10-3)加入125μl的SFB到磁珠中，悬浮磁珠，放回磁力架，直至澄清，弃上清；

(10-4)加入125μl的SFB再洗一次，去除上清，瞬离，再吸干净；

(10-5)用15μl的EB悬浮磁珠，放回磁力架，直至澄清，将上清转移至新的LoBindtube管中(不要碰到磁珠)，得到cDNA文库。

所述高通量测序的步骤具体可为：

(1)配制priming mix(测序mix)：30μl FLT直接加到1管FB中。

(2)将1000μl加样枪调到780μl，1μl慢慢的往上调排气泡，枪尖里有液体即可(大约调到800μl)。

(3)用1000μl的加样枪吸取800μl的priming mix，缓慢注入flow cell里，避免气泡，放置5min。

(4)配制文库。文库为75μl，包括37.5μl SQB、25.5μl LB(用前冲打混匀)和12μlcDNA文库。

(5)轻轻的打开SpotON sample port的盖，用1000μl的加样枪吸取200μl的priming mix，从priming port缓慢注入flow cell里，避免气泡。

(6)用200μl的加样枪吸取75μl文库一滴一滴加入到SpotON sample port里。

(7)盖上盖上机。

放10min左右loading后上机。

所述步骤(3)中，生物信息学分析，获得若干待测样本含有的病原体种类、丰度和/或基因组覆盖度的方法可包括(3-1)—(3-5)：

(3-1)Basecalling；

(3-2)拆分；

(3-3)质控；

(3-4)去宿主；

(3-5)物种比对；

(3-6)与参考基因组比对，确定病原体，获得若干待测样本含有的病原体种类、丰度和/或基因组覆盖度。

所述步骤(3)中，Basecalling可为采用如下命令生成fastq文件：guppy_basecaller-i./-s output/--config/raw/chb/lung_tissue_nanopore/ont-guppy-cpu/data/dna_r9.4.1_450bps_hac.cfg-r--num_callers 24--cpu_threads_per_caller 2。拆分可为采用如下命令拆分fastq文件：guppy_barcoder-i./-s barcoder_20200806_1--barcode_kits SQK-RPB004-t 12。质控可为采用nano_pp脚本去除低于500bp的短reads，同时统计reads数、reads长度和分布。去宿主可为采用nano_rm_host脚本去除人宿主序列。物种比对可为采用nano_classify脚本与数据库进行物种比对，得出物种名称、特异性reads数。与参考基因组比对可为采用nano_realign脚本将比对上的reads和参考基因组进行比对。确定病原体可为采用nano_remove_neg脚本将待测样本测得的reads数标化后确定。基因组覆盖度可为采用nano_stats脚本计算每个物种测到的reads覆盖参考基因组的覆盖度，同时以图形显示测得reads覆盖参考基因组情况。采用nano_realign脚本将比对上的reads重新和参考基因组进行直接的比对进一步确定鉴定物种的准确性，得出物种名称、特异性reads数、丰度。采用nano_stats脚本计算每个物种测到的reads覆盖参考基因组的覆盖度，同时以图形显示测得reads覆盖参考基因组情况。采用nano_remove_neg脚本将待测样本测得的reads数标化，同时和阴性对照相比，RPM比率＝RPM(样品)/RPM(阴性对照)大于10的物种为病原体。

所述数据库的构建方法如下：从已有基因组数据中下载细菌基因组库、古菌基因组库、真菌基因组库、病毒基因组库、非脊椎动物基因组库、原生生物基因组库和人全基因组数据，根据物种信息，每个物种选择一个代表序列，构建代表序列数据库；然后使用mummer将同一个物种的基因组进行比对去冗余，得到这个物种的pan-genome序列；去除冗余后的所有数据，使用centrifuge-build构建数据库。

所述已有基因组数据可为NCBI数据库收录的基因组数据。

在本发明的一个实施例中，数据库中共包含细菌5208个、古菌274个、真菌322个、病毒12947个、原生生物91个和非脊椎动物216个。

上述任一所述的方法在检测待测样本含有的病原体种类、丰度和/或基因组覆盖度中的应用也属于本发明的保护范围。

上述任一所述待测样本可为临床样本。所述临床样本可为脑脊液、肺泡灌洗液、穿刺液、血浆、宫颈分泌物或肺组织。

所述临床样本可置于DNA/RNA Shield稳定剂(Zymo，Catalog Code：R1100-50)中保存。临床样本为组织样本时，可先剪成小块再处理。临床样本为痰样本时，可先用痰消化液处理。

本发明将宏基因组测序和Nanopore单分子实时测序有机结合，应用于临床标本的病原体高通量检测。Nanopore单分子实时测序所需起始RNA量大，临床标本提取的RNA的量不能满足其建库要求，本发明采用随机引物扩增将RNA极大的富集，使其能够应用于临床。再次，本发明开发了基于Nanopore宏基因组RNA-seq的生物信息学分析方法，能够迅速、准确的检测临床样本中的病原体，为临床疑难感染性疾病的病原学诊断提供依据，使疑难感染患者的治疗有的放矢，大大改善疑难感染患者的诊治，使患者受益。本发明具有重要的应用价值。

具体实施方式

下面结合具体实施方式对本发明进行进一步的详细描述，给出的实施例仅为了阐明本发明，而不是为了限制本发明的范围。以下提供的实施例可作为本技术领域普通技术人员进行进一步改进的指南，并不以任何方式构成对本发明的限制。

下述实施例中的实验方法，如无特殊说明，均为常规方法，按照本领域内的文献所描述的技术或条件或者按照产品说明书进行。下述实施例中所用的材料、试剂等，如无特殊说明，均可从商业途径得到。

实施例1、基于Nanopore宏基因组RNA-seq的生物信息学检测病原体的方法的建立

一、临床样本的获得

取临床样本，置于DNA/RNA Shield稳定剂(Zymo，Catalog Code：R1100-50)中保存。

注：临床样本为组织样本时，需先剪成小块再处理。临床样本为痰样本时，需先用痰消化液处理。

二、RNA提取

采用qiagenAllPrepPowerViral DNA/RNA Kit(Catalog no.28000-50)提取临床样本的RNA。具体步骤依次如下：

1、向Glass PowerBead Tube中加入200μL临床样本；

2、加入600μL PM1/β-ME；

3、置于Qiagen匀浆器，最大转速震荡10min；

4、室温、13000g离心1min，收集上清并转移至收集管中；

5、向收集管中加入150μL IRS液体涡旋混匀，4℃孵育5min；

6、13000g离心1min，收集上清并转移至收集管中；

7、加入600μL PM3液体和600μL PM4液体，涡旋混匀，得到混合液体；

8、向MB Spin Column柱子中转移625μL混合液体，13000g离心1min；

9、向MB Spin Column柱子中加入600μL PM5液体，13000g离心1min；

10、去除流出液，加入600μL PM4液体，13000g离心1min；

11、去除流出液，13000g离心2min；

12、将MB Spin Column放到一新的收集管中；

13、加入50μL RNase-free水，孵育3min；

14、13000g离心1min，收集液体即为临床样本的RNA，-70℃保存。

PM1/β-ME、IRS液体、PM4液体、PM3液体和PM5液体均为qiagenAllPrepPowerViralDNA/RNA Kit中的组件。

三、RNA测序

1、退火

(1)制备退火反应体系。退火反应体系为13μl，包括1μl浓度为10μM的引物KN8溶液、1μl dNTP(浓度为10mM)和11μl临床样本的RNA。

KN8引物：5’-GACCATCTAGCGACCTCCACNNNNNNNN-3’(SEQ ID NO:1)。(N为A、T、G和C中的任一种)

(2)取所述退火反应体系，进行退火，得到退火产物。

退火程序为：65℃5min；冰上急冷>1min；热盖105℃。

2、逆转录

(1)制备逆转录反应体系。逆转录反应体系为20μl，包括13μl退火产物、4μl superIVase(Thermo Fisher)、1μl 100mM DDT、1μl RNase 007inhibitor和1μl super IVase(Thermo Fisher)。

(2)取所述逆转录反应体系，进行逆转录，得到逆转录产物。

逆转录程序为：42℃50min；70℃10min；4℃∞。

3、cDNA的获得

向逆转录产物中加入1μl lemon酶，之后混匀，进行二链合成，得到cDNA。

反应程序为：37℃3min；75℃15min；4℃保存。

4、PCR扩增

(1)制备PCR反应体系。PCR反应体系为50μl，包括5μl cDNA、25μl 2×mix(Q5溶液)、2μl primer K和18μl无核酸水。

primer K：5’-GACCATCTAGCGACCTCCAC-3’(SEQ ID NO:2)。

(2)取所述PCR反应体系，进行PCR扩增，得到PCR扩增产物。

PCR反应程序为：98℃30s；98℃10s，55℃15s，72℃1min，35cycles；72℃10min。

5、取PCR扩增产物，采用TakaRaMiniBEST DNA Fragment Purification kit(cat#9761lot#AH70942A)进行DNA纯化(目的为去除小片段)，得到纯化产物。

6、末端修复

(1)制备反应体系。反应体系为15μl，包括Xμl纯化产物(含50ng DNA)、(12.5-X)μl无核酸水、1.75μl ultraⅡEnd-prep reaction buffer和0.75μl ultraⅡEnd-prepenzyme mix。

(2)取所述反应体系，进行末端修复，得到末段修饰后的cDNA。

末端修复的程序为：20℃3min；65℃5min。

7、加建库接头即barcode

试剂盒EXP-NBD104：barcode 1-12；试剂盒EXP-NBD114：barcode 13-24。

(1)制备反应体系。反应体系为20μl，包括4μl无核酸水、3μl末段修饰后的cDNA、2.5μl barcode、10μl Master mix和0.5μl Enhancer(目的为将barcode加在DNA末端)。

(2)取所述反应体系，反应，得到连接了barcode的样品。

反应程序为：20℃20min；65℃10min。

8、混样纯化

(1)将所有连接了barcode的样品放入一个1.5ml的LoBind tube中，然后加入磁珠(总cDNA体积：磁珠体积＝1:0.4)，充分混合；

(2)室温孵育10min后，放到磁力架上，直至澄清；

(3)去除上清(留一点)(在磁力架上操作)；

(4)加入500μl 80％乙醇水溶液洗涤，旋转EP管；

(5)去除上清，80％乙醇水溶液再洗一次；去除上清，瞬离，再吸干净，开盖放置30s，晾干；

(6)加入35μl无核酸水，室温放置5min，洗脱DNA；

(7)将LoBind tube管重新放回磁力架上，液体澄清后，将上清转移至新的LoBindtube管中(不要碰到磁珠)，得到纯化后cDNA样品(已加接头)；

(8)Qubit定量浓度(198μl qubit工作液+2μl纯化后cDNA样品(已加接头))，并记录。

9、接头连接和清洗

(1)制备反应体系。反应体系为50μl，包括Xμl纯化后cDNA样品(含30-50ng cDNA)(已加接头)、(30-X)μl无核酸水、5μl AMⅡ、10μl 5×NEB Quick Ligation buffer和5μlT

(2)取所述反应体系，21℃20min(开热盖)。

(3)纯化

(3-1)向完成步骤(2)的反应体系汇总加入20μl磁珠，室温孵育10min(10min内反复颠倒混匀)；

(3-2)移到磁力架上，直至澄清，弃上清；

(3-3)加入125μl的SFB到磁珠中，悬浮磁珠，放回磁力架，直至澄清，弃上清；

(3-4)加入125μl的SFB再洗一次，去除上清，瞬离，再吸干净；

(3-5)用15μl的EB悬浮磁珠，放回磁力架，直至澄清，将上清转移至新的LoBindtube管中(不要碰到磁珠)，得到cDNA文库；

(3-6)Qubit定量浓度(198μl qubit工作液+2μlcDNA文库)，并记录。

10、混样测序

(1)配制priming mix(测序mix)：30μl FLT直接加到1管FB中。

(2)将1000μl加样枪调到780μl，1μl慢慢的往上调排气泡，枪尖里有液体即可(大约调到800μl)。

(3)用1000μl的加样枪吸取800μl的priming mix，缓慢注入flow cell里，避免气泡，放置5min。

(4)配制文库。文库为75μl，包括37.5μl SQB、25.5μl LB(用前冲打混匀)和12μlcDNA文库。

(5)轻轻的打开SpotON sample port的盖，用1000μl的加样枪吸取200μl的priming mix，从priming port缓慢注入flow cell里，避免气泡。

(6)用200μl的加样枪吸取75μl文库一滴一滴加入到SpotON sample port里。

(7)盖上盖上机。

放10min左右loading后上机。

11、MinKnow操作

(1)Experiment：输入实验名

(2)Kit：选择对应的建库kit

(3)Basecalling：Fast basecalling(电脑性能不够好时用此模式)

Basecalling and barcoding处于“ON”状态

(4)Run options：

Time：72h Bias voltage(mV)：-180Active channel selection：“ON”

1.5h扫描1次

Start Run。

四、生物信息学分析

1、Basecalling

采用如下命令生成fastq文件：

guppy_basecaller-i./-s output/--config/raw/chb/lung_tissue_nanopore/ont-guppy-cpu/data/dna_r9.4.1_450bps_hac.cfg-r--num_callers 24--cpu_threads_per_caller 2

2、拆分

采用如下命令拆分fastq文件：

guppy_barcoder-i./-s barcoder_20200806_1--barcode_kits SQK-RPB004-t12

3、质控

采用nano_pp脚本去除低于500bp的短reads，同时统计reads数、reads长度和分布。

4、去宿主

采用nano_rm_host脚本去除人宿主序列。

5、物种比对

采用nano_classify脚本与数据库进行物种比对，得出物种名称、特异性reads数。

数据库构建的方法如下：数据库全部采用NCBI数据库收录的基因组数据，选择其中的细菌基因组库、古菌基因组库、真菌基因组库、病毒基因组库、非脊椎动物基因组库和原生生物基因组库，此外数据库还包含了人的全基因组数据。基因组数据按照分类下载后，根据物种信息，每个物种选择一个代表序列，构建代表序列数据库。然后使用mummer将同一个物种的基因组进行比对去冗余，得到这个物种的pan-genome序列。去除冗余后的所有数据，使用centrifuge-build构建比对数据库。

最终数据库中共包含细菌5208个，古菌274个，真菌322个，病毒12947个，原生生物91个，非脊椎动物216个。

6、重新比对

采用nano_realign脚本将比对上的reads重新和参考基因组进行直接的比对进一步确定鉴定物种的准确性，得出物种名称、特异性reads数、丰度。

7、基因组覆盖度

采用nano_stats脚本计算每个物种测到的reads覆盖参考基因组的覆盖度，同时以图形显示测得reads覆盖参考基因组情况。

8、报告

采用nano_remove_neg脚本将待测样本测得的reads数标化，同时和阴性对照相比，RPM比率＝RPM(样品)/RPM(阴性对照)大于10的物种为病原体。

实施例2、实施例建立的方法的有效性检测

28份临床样本均由北京大学人民医院检验科微生物实验室提供，且临床样本的提供者均知情同意。

28份临床样本的样本号、临床诊断结果和标准类型见表1中第1-3列。

1、采用常规微生物培养方法分别检测28份临床样本的病原体感染情况。

检测结果见表1中第4列。

2、采用实施例1提供的方法高通量检测28份临床样本的病原体感染情况。

检测结果见表1中第5列。

表1

结果表明，与常规微生物培养方法相比，实施例1提供的方法的敏感性为83.3％，特异性为90.9％，阳性预测值为71.4％，阴性预测值为95.2％，符合率为89.3％，诊断比值比(OR)为49.7。

以上对本发明进行了详述。对于本领域技术人员来说，在不脱离本发明的宗旨和范围，以及无需进行不必要的实验情况下，可在等同参数、浓度和条件下，在较宽范围内实施本发明。虽然本发明给出了特殊的实施例，应该理解为，可以对本发明作进一步的改进。总之，按本发明的原理，本申请欲包括任何变更、用途或对本发明的改进，包括脱离了本申请中已公开范围，而用本领域已知的常规技术进行的改变。按以下附带的权利要求的范围，可以进行一些基本特征的应用。

<110> 北京大学人民医院

<120> 一种基于Nanopore宏基因组RNA-seq的生物信息学检测病原体的方法

<160>2

<170> PatentIn version 3.5

<210>1

<211>28

<212> DNA

<213> Artificial sequence

<400>1

gaccatctag cgacctccac nnnnnnnn 28

<210>2

<211>20

<212> DNA

<213> Artificial sequence

<400>2

gaccatctag cgacctccac 20