• 转录组入门(5):序列比对


    任务列表
    • 比对软件
    • hisat2的用法
    • 下载index文件
    • 比对、排序、索引
    • 质量控制
    • 载入IGV,截图几个基因
    hisat2的用法
    本作业是比对到基因组,所以使用gapped or splices mapper,此流程已经更新。TopHat首次被发表已经是7年前,STAR的比对速度是TopHat的50倍,HISAT更是STAR的1.2倍。HISAT2是TopHat2/Bowti2的继任者,使用改进的BWT算法,实现了更快的速度和更少的资源占用,作者推荐TopHat2/Bowti2和HISAT的用户转换到HISAT2。
    官网:https://ccb.jhu.edu/software/hisat2/index.shtml(学习一个软件最好的方法就是结合现有中文资料,加上阅读官方说明书和HELP文档,一般刚开始学习的时候先使用默认参数,不要乱调参数)
    下载index文件
    cd ~/reference
    mkdir -p index/hisat && cd index/hisat
    wget -c ftp://ftp.ccb.jhu.edu/pub/infphilo/hisat2/data/hg19.tar.gz
    wget -c ftp://ftp.ccb.jhu.edu/pub/infphilo/hisat2/data/mm10.tar.gz
    tar zxvf hg19.tar.gz
    tar xvzf mm10.tar.gz
    -c:断点续传
    比对、排序、索引
    把fastq格式的reads比对上去得到sam文件,接着用samtools把它转为bam文件,并且排序(注意N和P两种排序区别)索引好(可以使用管道实现,省去中间SAM保存的过程,直接输出BAM文件)
    编写bash脚本:map.sh
    #! usr/bin/bash
    set -u
    set -e
    set -o pipefail
    hg19_ref=/mnt/hgfs/2017/reference/index/hisat/hg19/genome
    mm10_ref=/mnt/hgfs/2017/reference/index/hisat/mm10/genome
    data_path=/mnt/hgfs/2017/rna_seq/data
    NUM_THREADS=25
    ls --color=never Homo*1.fastq.gz | while read id;do(~/biosoft/hisat2-2.1.0/hisat2 -t -p $NUM_THREADS -x $hg19_ref -1 $data_path/${id%_*}_1.fastq.gz -2 $data_path/${id%_*}_2.fastq.gz 2 > ${id%_*}_map.log | samtools view -Sb  - > ${id%_*}.bam);done
    ls --color=never Mus*1.fastq.gz | while read id;do(~/biosoft/hisat2-2.1.0/hisat2 -t -p $NUM_THREADS -x $mm10_ref -1 $data_path/${id%_*}_1.fastq.gz -2 $data_path/${id%_*}_2.fastq.gz 2 > ${id%_*}_map.log | samtools view -Sb  - > ${id%_*}.bam);done
    ls --color=never *.bam | while read id;do(samtools sort --threads $NUM_THREADS $id -o ${id%.*}_sorted.bam);done
    ls --color=never *_sorted.bam | while read id;do(samtools index $id);done
    运行脚本: 
    bash map.sh
    质量控制
    对bam文件进行简单QC
    Reads比对后的质量控制(评估比对质量的指标)
    比对上的reads占总reads的百分比;
    Reads比对到外显子和参考链上的覆盖度是否一致;
    比对到基因组序列,多重比对reads;
    相关质控软件除了Picard,RSeQC,Qualimap还有一大堆
  • 相关阅读:
    quartz 中JobExecutionContext的使用
    Memcached 集群架构方面的问题
    Spring+Quartz 集群
    Spring Quartz 持久化解决方案
    不同版本(2.3,2.4,2.5) web.xml 的web-app头信息
    Spring Framework 4.0.0发布,首次支持Java 8
    Serializable java序列化
    Quartz 有状态的JobDataMap
    JobDataMap 不能被序列化如何解决研究中
    Spring-3.2.5 + Quartz-2.2.1 集群实例(Tomcat+Memcached+Quartz集群session共享)
  • 原文地址:https://www.cnblogs.com/freescience/p/7342895.html
Copyright © 2020-2023  润新知