• 新冠病毒分型和突变分析(SARS-CoV2_ARTIC_Nanopore)


    新冠病毒分型和突变分析(SARS-CoV2_ARTIC_Nanopore)

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

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

    概览

    流程输入SRR14800265.fastq.gz

    测试数据下载

    SRX11133330: Nanopore sequencing of SARS-CoV-2: V-22
    1 OXFORD_NANOPORE (MinION) run: 277,605 spots, 139.3M bases, 133.2Mb downloads
    使用NCBI官方工具sra-toolkit拆分成fastq.gz文件 fastq-dump SRR14800265 --gzip
    得到SRR14800265.fastq.gz

    参考文件,默认路径/opt/ref下
    Artic-ncov2019
    artic-ncov2019 primer&参考序列

    分析流程文件(可一键导入sliverworkspace运行)及报告文件,conda环境文件下载,导入操作
    运行环境docker image based on ubuntu21.04 Conda Mamba(默认使用清华源) ssh
    获取镜像代码见下文段落
    分析软件- artic=1.2.1
    - artic-network::rampart=1.2.0
    - snakemake-minimal=5.8.1
    - pangolin=4.1.3
    输出结果按照序列一致性组装的新冠病毒序列 SRR14800265.consensus.fa
    Panglin 根据组装的序列分析得出病毒分型信息 lineage_report.csv
    根据primertrim.bam获的新冠病毒突变信息,过滤后得到 SRR14800265.pass.vcf.gz

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

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

    备注:docker运行的操作系统,推荐为Linux,windows,macOS系统下docker可能部分功能(网络)不能正常运行

    # 拉取docker镜像
    docker     pull     doujiangbaozi/sliverworkspace:latest
    
    # 查看docker 镜像
    docker     images
    
    • 1
    • 2
    • 3
    • 4
    • 5
    基础环境配置,docker-compose.yml 配置文件,可以根据需要自行修改调整
    version: "3"
    services:
      SarsCov2:
        image: doujiangbaozi/sliverworkspace:latest
        container_name: SarsCov2
        volumes:
          - /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端口
    
    
    • 1
    • 2
    • 3
    • 4
    • 5
    • 6
    • 7
    • 8
    • 9
    • 10
    • 11
    • 12
    • 13
    • 14
    • 15
    • 16
    • 17
    • 18

    基础环境运行

    # docker-compose.yml 所在目录下运行
    docker-compose up -d
    
    # 或者 
    docker-compose up -d -f /路径/docker-compose.yaml
    
    # 查看docker是否正常运行,docker-compose.yaml目录下运行
    docker-compose ps
    
    # 或者
    docker ps
    
    • 1
    • 2
    • 3
    • 4
    • 5
    • 6
    • 7
    • 8
    • 9
    • 10
    • 11

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

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

    三. 分析流程:本流程包含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
    
    • 1
    • 2
    • 3
    • 4
    • 5
    • 6
    • 7
    • 8
    • 9
    • 10
    • 11
    • 12

    2. 数据过滤

    #首次运行下载artic-ncov2019
    if  [ ! -d "/opt/ref/artic-ncov2019" ]; then
    	apt-get install -y git
        git clone https://github.com/artic-network/artic-ncov2019.git "/opt/ref/artic-ncov2019"
    fi
    
    #conda检测环境是否存在,首次运行不存在创建该环境并安装软件
    if [ ! -d "${envs}/artic-ncov2019" ]; then
      mamba env create -f /opt/ref/artic-ncov2019/environment.yml
      mamba install muscle=3.8
    fi
    
    source activate artic-ncov2019
    
    cp    -f ${data}/artic/${sn}.fastq.gz ${result}/${sn}/
    
    mkdir -p ${result}/${sn}/clean
    
    artic guppyplex --min-length 400 --max-length 700 \
    	--directory ${result}/${sn}/ \
        --output ${result}/${sn}/clean/${sn}.clean.fastq
    
    conda  deactivate
    
    • 1
    • 2
    • 3
    • 4
    • 5
    • 6
    • 7
    • 8
    • 9
    • 10
    • 11
    • 12
    • 13
    • 14
    • 15
    • 16
    • 17
    • 18
    • 19
    • 20
    • 21
    • 22
    • 23

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

    source activate artic-ncov2019
    
    if [ ! -f /opt/ref/artic-ncov2019/primer_schemes/nCoV-2019/${artic_primer_version}/nCoV-2019.reference.fasta ]; then
    	cp  -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
    fi
    
    cd ${result}/${sn}
    
    mkdir -p ${result}/${sn}/aligned
    
    minimap2 -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.bam
    
    conda  deactivate
    
    • 1
    • 2
    • 3
    • 4
    • 5
    • 6
    • 7
    • 8
    • 9
    • 10
    • 11
    • 12
    • 13
    • 14
    • 15
    • 16
    • 17
    • 18

    4. trim bam / samtools coverage

    source activate artic-ncov2019
    
    cd ${result}/${sn}
    
    if [ ! -f /opt/ref/artic-ncov2019/primer_schemes/nCoV-2019/${artic_primer_version}/nCoV-2019.scheme.bed ]; then
    	cp  -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
    fi
    
    align_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.bam
    
    align_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.bam
    
    samtools index ${result}/${sn}/aligned/${sn}.trimmed.rg.sorted.bam
    samtools index ${result}/${sn}/aligned/${sn}.primertrimmed.rg.sorted.bam
    
    samtools coverage ${result}/${sn}/aligned/${sn}.primertrimmed.rg.sorted.bam -o ${result}/${sn}/aligned/${sn}.samcov.tsv
    
    conda  deactivate
    
    • 1
    • 2
    • 3
    • 4
    • 5
    • 6
    • 7
    • 8
    • 9
    • 10
    • 11
    • 12
    • 13
    • 14
    • 15
    • 16
    • 17
    • 18
    • 19
    • 20
    • 21
    • 22
    • 23
    • 24
    • 25

    5. medaka variant

    source activate artic-ncov2019
    
    if [ ! -f /opt/ref/artic-ncov2019/primer_schemes/nCoV-2019/${artic_primer_version}/nCoV-2019.reference.fasta ]; then
    	cp  -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
    fi
    
    mkdir -p ${result}/${sn}/vcf
    
    if [ -f ${result}/${sn}/vcf/${sn}.nCoV-2019_1.hdf ];then
    	rm -f ${result}/${sn}/vcf/${sn}.nCoV-2019_1.hdf
    fi
    
    if [ -f ${result}/${sn}/vcf/${sn}.nCoV-2019_2.hdf ];then
    	rm -f ${result}/${sn}/vcf/${sn}.nCoV-2019_2.hdf
    fi
    
    medaka 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.hdf
    
    medaka 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.vcf
    
    medaka 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.hdf
    
    medaka 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.vcf
    
    artic_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.vcf
    
    bgzip -f ${result}/${sn}/vcf/${sn}.merged.vcf
    tabix -f -p vcf ${result}/${sn}/vcf/${sn}.merged.vcf.gz
    
    conda  deactivate
    
    • 1
    • 2
    • 3
    • 4
    • 5
    • 6
    • 7
    • 8
    • 9
    • 10
    • 11
    • 12
    • 13
    • 14
    • 15
    • 16
    • 17
    • 18
    • 19
    • 20
    • 21
    • 22
    • 23
    • 24
    • 25
    • 26
    • 27
    • 28
    • 29
    • 30
    • 31
    • 32
    • 33
    • 34
    • 35
    • 36
    • 37
    • 38
    • 39
    • 40
    • 41

    6. longshot vcf

    source activate artic-ncov2019
    
    longshot -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.gz
    
    conda  deactivate
    
    • 1
    • 2
    • 3
    • 4
    • 5
    • 6
    • 7
    • 8
    • 9

    7. artic_vcf_filter 过滤variant vcf

    source activate artic-ncov2019
    
    artic_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.gz
    
    conda  deactivate
    
    • 1
    • 2
    • 3
    • 4
    • 5
    • 6
    • 7
    • 8
    • 9

    8. depth mask

    source activate artic-ncov2019
    
    artic_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.txt
    
    artic_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.fasta
    
    conda  deactivate
    
    • 1
    • 2
    • 3
    • 4
    • 5
    • 6
    • 7
    • 8
    • 9
    • 10
    • 11
    • 12
    • 13

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

    source activate artic-ncov2019
    
    bcftools 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.fasta
    
    artic_fasta_header ${result}/${sn}/${sn}.consensus.fasta "${sn}"
    
    conda  deactivate
    
    • 1
    • 2
    • 3
    • 4
    • 5
    • 6
    • 7
    • 8
    • 9
    • 10

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

    #conda检测环境是否存在,首次运行不存在创建该环境并安装软件
    if [ ! -d "${envs}/pangolin" ]; then
      mamba env create -f /opt/config/pangolin.yaml
    fi
    
    source activate pangolin
    
    pangolin ${result}/${sn}/${sn}.consensus.fasta --outdir ${result}/${sn}
    
    conda  deactivate
    
    • 1
    • 2
    • 3
    • 4
    • 5
    • 6
    • 7
    • 8
    • 9
    • 10

    11. muscle (可选)

    source activate artic-ncov2019
    
    cat ${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.fasta
    
    muscle	-in ${result}/${sn}/${sn}.muscle.in.fasta -out ${result}/${sn}/${sn}.muscle.out.fasta
    
    conda  deactivate
    
    • 1
    • 2
    • 3
    • 4
    • 5
    • 6
    • 7
    • 8
    • 9

    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
    
    
    • 1
    • 2
    • 3
    • 4
    • 5
    • 6
    • 7
    • 8
    • 9
    • 10
    • 11
    • 12
  • 相关阅读:
    【保姆级】lookup-method标签实践与分析
    查找算法之二分查找
    卷麻了,面试了一个00后,绝对能称为是卷王之王....
    会员题-力扣408-有效单词缩写
    服务状态巡检:
    Python Tcp编程
    nodejs 模块
    Android---彻底掌握 Handler
    Vue组件的计算属性和普通属性区别、属性侦听器的作用
    Windows安装mysql详细步骤(通俗易懂,简单上手)
  • 原文地址:https://blog.csdn.net/weixin_39900139/article/details/128124509