新冠病毒分型和突变分析（SARS-CoV2_ARTIC_Nanopore）

一. 本文使用Artic官方提供环境对Nanopore minion SARS-Cov-2测序数据，对新冠病毒突变及分型鉴定

二. 概览：按照惯例，先上一张概览图，浏览下分析流程步骤

环境搭建: 为了快速完成环境搭建，节省95%以上时间。

本文使用docker + conda (mamba) 作为基础分析环境，镜像获取：docker/docker-compoes 的安装及镜像构建见《基于docker的生信基础环境镜像构建》,docker镜像基于ubuntu21.04构建，并安装有conda/mamba，ssh服务。并尝试初次运行时初始化安装所需软件下载所需文件（作为代价首次运行时间会较长，切需网络通畅），即实现自动初始化的分析流程。

# 拉取docker镜像
docker     pull     doujiangbaozi/sliverworkspace:latest# 查看docker 镜像
docker     images

基础环境配置，docker-compose.yml 配置文件，可以根据需要自行修改调整

version: "3"
services:SarsCov2:image: doujiangbaozi/sliverworkspace:latestcontainer_name: SarsCov2volumes:- /media/sliver/Data/data:/opt/data:rw                               #挂载input数据，artic目录下- /media/sliver/Manufacture/SC2/envs:/root/mambaforge-pypy3/envs:rw  #挂载envs conda环境目录- /media/sliver/Manufacture/SC2/config:/opt/config:rw                #挂载config,conda配置文件目录- /media/sliver/Manufacture/SC2/ref:/opt/ref:rw                      #挂载reference目录- /media/sliver/Manufacture/SC2/result:/opt/result:rw                #挂载中间文件和输出结果目录ports:- "9024:9024"                                                        #ssh连接端口可以按需修改environment:- TZ=Asia/Shanghai                                                   #设置时区- PS=20191124                                                        #修改默认ssh密码- PT=9024                                                            #修改默认ssh端口

基础环境运行

# docker-compose.yml 所在目录下运行
docker-compose up -d# 或者 
docker-compose up -d -f /路径/docker-compose.yaml# 查看docker是否正常运行,docker-compose.yaml目录下运行
docker-compose ps# 或者
docker ps

docker 容器使用，类似于登录远程服务器

# 登录docker，使用的是ssh服务，可以本地或者远程部署使用
ssh root@192.168.6.6 -p9024# 看到如下，显示如下提示即正常登录
(base) root@SliverWorkstation:~# 

三. 分析流程：本流程包含artic流程详细步骤，不感兴趣可以直接跳至文章末尾

1. 变量设置

#样本编号
export sn=SRR14800265
#数据输入目录
export data=/opt/data
#数据输出、中间文件目录
export result=/opt/result
#conda安装的环境目录
export envs=/root/mambaforge-pypy3/envs	
#artic primer 版本V1,V2,V3,V4,V4.1
export artic_primer_version=4.1
#设置可用线程数
export threads=8

2. 数据过滤

#首次运行下载artic-ncov2019
if  [ ! -d "/opt/ref/artic-ncov2019" ]; thenapt-get install -y gitgit clone https://github.com/artic-network/artic-ncov2019.git "/opt/ref/artic-ncov2019"
fi#conda检测环境是否存在，首次运行不存在创建该环境并安装软件
if [ ! -d "${envs}/artic-ncov2019" ]; thenmamba env create -f /opt/ref/artic-ncov2019/environment.ymlmamba install muscle=3.8
fisource activate artic-ncov2019cp    -f ${data}/artic/${sn}.fastq.gz ${result}/${sn}/mkdir -p ${result}/${sn}/cleanartic guppyplex --min-length 400 --max-length 700 \--directory ${result}/${sn}/ \--output ${result}/${sn}/clean/${sn}.clean.fastqconda  deactivate

3. 比对到参考基因组上，得到bam文件并排序

source activate artic-ncov2019if [ ! -f /opt/ref/artic-ncov2019/primer_schemes/nCoV-2019/${artic_primer_version}/nCoV-2019.reference.fasta ]; thencp  -f /opt/ref/artic-ncov2019/primer_schemes/nCoV-2019/${artic_primer_version}/SARS-CoV-2.reference.fasta \/opt/ref/artic-ncov2019/primer_schemes/nCoV-2019/${artic_primer_version}/nCoV-2019.reference.fasta
ficd ${result}/${sn}mkdir -p ${result}/${sn}/alignedminimap2 -a -x map-ont -t ${threads} \/opt/ref/artic-ncov2019/primer_schemes/nCoV-2019/${artic_primer_version}/nCoV-2019.reference.fasta \${result}/${sn}/clean/${sn}.clean.fastq  | samtools view -bS -F 4 - | samtools sort -o ${result}/${sn}/aligned/${sn}.sorted.bam -samtools index ${result}/${sn}/aligned/${sn}.sorted.bamconda  deactivate

4. trim bam / samtools coverage

source activate artic-ncov2019cd ${result}/${sn}if [ ! -f /opt/ref/artic-ncov2019/primer_schemes/nCoV-2019/${artic_primer_version}/nCoV-2019.scheme.bed ]; thencp  -r /opt/ref/artic-ncov2019/primer_schemes/nCoV-2019/${artic_primer_version}/SARS-CoV-2.scheme.bed \/opt/ref/artic-ncov2019/primer_schemes/nCoV-2019/${artic_primer_version}/nCoV-2019.scheme.bed
fialign_trim --normalise 200 /opt/ref/artic-ncov2019/primer_schemes/nCoV-2019/${artic_primer_version}/nCoV-2019.scheme.bed \--start --remove-incorrect-pairs --report ${result}/${sn}/aligned/${sn}.alignreport.txt < ${result}/${sn}/aligned/${sn}.sorted.bam \2> ${result}/${sn}/aligned/${sn}.alignreport.er | \samtools sort -T ${result}/${sn}/aligned/temp - -o ${result}/${sn}/aligned/${sn}.trimmed.rg.sorted.bamalign_trim --normalise 200 /opt/ref/artic-ncov2019/primer_schemes/nCoV-2019/${artic_primer_version}/nCoV-2019.scheme.bed \--remove-incorrect-pairs --report ${result}/${sn}/aligned/${sn}.alignreport.txt < ${result}/${sn}/aligned/${sn}.sorted.bam \2> ${result}/${sn}/aligned/${sn}.alignreport.er| \samtools sort -T ${result}/${sn}/aligned/temp - -o ${result}/${sn}/aligned/${sn}.primertrimmed.rg.sorted.bamsamtools index ${result}/${sn}/aligned/${sn}.trimmed.rg.sorted.bam
samtools index ${result}/${sn}/aligned/${sn}.primertrimmed.rg.sorted.bamsamtools coverage ${result}/${sn}/aligned/${sn}.primertrimmed.rg.sorted.bam -o ${result}/${sn}/aligned/${sn}.samcov.tsvconda  deactivate

5. medaka variant

source activate artic-ncov2019if [ ! -f /opt/ref/artic-ncov2019/primer_schemes/nCoV-2019/${artic_primer_version}/nCoV-2019.reference.fasta ]; thencp  -f /opt/ref/artic-ncov2019/primer_schemes/nCoV-2019/${artic_primer_version}/SARS-CoV-2.reference.fasta \/opt/ref/artic-ncov2019/primer_schemes/nCoV-2019/${artic_primer_version}/nCoV-2019.reference.fasta
fimkdir -p ${result}/${sn}/vcfif [ -f ${result}/${sn}/vcf/${sn}.nCoV-2019_1.hdf ];thenrm -f ${result}/${sn}/vcf/${sn}.nCoV-2019_1.hdf
fiif [ -f ${result}/${sn}/vcf/${sn}.nCoV-2019_2.hdf ];thenrm -f ${result}/${sn}/vcf/${sn}.nCoV-2019_2.hdf
fimedaka consensus --model r941_min_high_g351 \--threads ${threads} --chunk_len 800 --chunk_ovlp 400 \--RG 1 ${result}/${sn}/aligned/${sn}.trimmed.rg.sorted.bam \${result}/${sn}/vcf/${sn}.nCoV-2019_1.hdfmedaka variant /opt/ref/artic-ncov2019/primer_schemes/nCoV-2019/${artic_primer_version}/nCoV-2019.reference.fasta \${result}/${sn}/vcf/${sn}.nCoV-2019_1.hdf ${result}/${sn}/vcf/${sn}.nCoV-2019_1.vcfmedaka consensus --model r941_min_high_g351 \--threads ${threads} --chunk_len 800 --chunk_ovlp 400 \--RG 2 ${result}/${sn}/aligned/${sn}.trimmed.rg.sorted.bam \${result}/${sn}/vcf/${sn}.nCoV-2019_2.hdfmedaka variant /opt/ref/artic-ncov2019/primer_schemes/nCoV-2019/${artic_primer_version}/nCoV-2019.reference.fasta \${result}/${sn}/vcf/${sn}.nCoV-2019_2.hdf ${result}/${sn}/vcf/${sn}.nCoV-2019_2.vcfartic_vcf_merge ${result}/${sn}/vcf/${sn} /opt/ref/artic-ncov2019/primer_schemes/nCoV-2019/${artic_primer_version}/nCoV-2019.scheme.bed \2> ${result}/${sn}/vcf/${sn}.primersitereport.txt \nCoV-2019_1:${result}/${sn}/vcf/${sn}.nCoV-2019_1.vcf nCoV-2019_2:${result}/${sn}/vcf/${sn}.nCoV-2019_2.vcfbgzip -f ${result}/${sn}/vcf/${sn}.merged.vcf
tabix -f -p vcf ${result}/${sn}/vcf/${sn}.merged.vcf.gzconda  deactivate

6. longshot vcf

source activate artic-ncov2019longshot -P 0 -F -A --no_haps \--bam ${result}/${sn}/aligned/${sn}.primertrimmed.rg.sorted.bam \--ref /opt/ref/artic-ncov2019/primer_schemes/nCoV-2019/${artic_primer_version}/nCoV-2019.reference.fasta \--out ${result}/${sn}/vcf/${sn}.longshoted.vcf \--potential_variants ${result}/${sn}/vcf/${sn}.merged.vcf.gzconda  deactivate

7. artic_vcf_filter 过滤variant vcf

source activate artic-ncov2019artic_vcf_filter --medaka ${result}/${sn}/vcf/${sn}.longshoted.vcf \${result}/${sn}/vcf/${sn}.pass.vcf \${result}/${sn}/vcf/${sn}.fail.vcf
bgzip -f ${result}/${sn}/vcf/${sn}.pass.vcf
tabix -p vcf ${result}/${sn}/vcf/${sn}.pass.vcf.gzconda  deactivate

8. depth mask

source activate artic-ncov2019artic_make_depth_mask --store-rg-depths \/opt/ref/artic-ncov2019/primer_schemes/nCoV-2019/${artic_primer_version}/nCoV-2019.reference.fasta \${result}/${sn}/aligned/${sn}.trimmed.rg.sorted.bam \${result}/${sn}/${sn}.coverage_mask.txtartic_mask /opt/ref/artic-ncov2019/primer_schemes/nCoV-2019/${artic_primer_version}/nCoV-2019.reference.fasta \${result}/${sn}/${sn}.coverage_mask.txt \${result}/${sn}/vcf/${sn}.fail.vcf \${result}/${sn}/${sn}.preconsensus.fastaconda  deactivate

9. 使用 bcftools 获取一致性序列

source activate artic-ncov2019bcftools consensus \-f ${result}/${sn}/${sn}.preconsensus.fasta ${result}/${sn}/vcf/${sn}.pass.vcf.gz \-m ${result}/${sn}/${sn}.coverage_mask.txt \-o ${result}/${sn}/${sn}.consensus.fastaartic_fasta_header ${result}/${sn}/${sn}.consensus.fasta "${sn}"conda  deactivate

10. 使用Pangolin获取序列分型信息

#conda检测环境是否存在，首次运行不存在创建该环境并安装软件
if [ ! -d "${envs}/pangolin" ]; thenmamba env create -f /opt/config/pangolin.yaml
fisource activate pangolinpangolin ${result}/${sn}/${sn}.consensus.fasta --outdir ${result}/${sn}conda  deactivate

11. muscle （可选）

source activate artic-ncov2019cat ${result}/${sn}/${sn}.consensus.fasta \/opt/ref/artic-ncov2019/primer_schemes/nCoV-2019/${artic_primer_version}/nCoV-2019.reference.fasta \> ${result}/${sn}/${sn}.muscle.in.fastamuscle	-in ${result}/${sn}/${sn}.muscle.in.fasta -out ${result}/${sn}/${sn}.muscle.out.fastaconda  deactivate

12. 报告(可选)

11. 使用IGV Browser查看突变信息（可选）

四. 参考链接

https://github.com/artic-network/artic-ncov2019

#官方建议分析过程#数据过滤
artic guppyplex --min-length 400 --max-length 700 \--directory /opt/result/SRR14800265/ \--output /opt/result/SRR14800265/clean/SRR14800265.clean.fastq
#获取一致性序列和突变数据
artic minion --normalise 200 --threads 32 \--medaka --medaka-model r941_min_high_g351 \--scheme-directory /opt/ref/artic-ncov2019/primer_schemes \--read-file /opt/result/SRR14800265/clean/SRR14800265.clean.fastq nCoV-2019/V4.1 SRR14800265