您的位置:新葡亰496net > 电脑系统 > 新葡亰496net:队列比对,转录组入门

新葡亰496net:队列比对,转录组入门

发布时间:2019-12-29 22:50编辑:电脑系统浏览(162)

    职分列表

    参考http://www.jianshu.com/p/681e02e7f9af
    http://www.jianshu.com/p/93f96e7538da
    任务:

    1. 比对软件超多,首先大家去搜罗一下,因为大家是带大家入门,请统风华正茂用hisat2,况且搞懂它的用法。
    2. 直白去hisat2的主页下载index文件就可以,然后把fastq格式的reads比对上去获得sam文件。
    3. 紧接着用samtools把它转为bam文件,並且排序(注意N和P二种排序不同卡塔尔(قطر‎索引好,载入IGV,再截图几个基因看看!
    4. 顺便对bam文件进行轻巧QC,参照他事他说加以考察直播自个儿的基因组种类。

    HISAT2:比对到基因组

    • 比对软件
    • hisat2的用法
    • 下载index文件
    • 比对、排序、索引
    • 质量调节
    • 载入IGV,截图几个基因

    1. 比对软件

    • HISAT2:http://ccb.jhu.edu/software/hisat2/index.shtml
      参考资料:http://blog.biochen.com/archives/337
    • STAR:https://codeload.github.com/alexdobin/STAR/zip/master
      仿照效法资料:http://www.bio-info-trainee.com/727.html
    • TopHat:http://ccb.jhu.edu/software/tophat/index.shtml
      参谋资料:http://blog.sina.com.cn/s/blog_8808cae20101amqp.html
    • RapMap:https://github.com/COMBINE-lab/RapMap
      参谋文献:https://academic.oup.com/bioinformatics/article/32/12/i192/2288985/RapMap-a-rapid-sensitive-and-accurate-tool-for
    • CIDANE:http://ccb.jhu.edu/software/cidane/
      参谋文献:https://genomebiology.biomedcentral.com/articles/10.1186/s13059-015-0865-0
    • CLASS2 :https://sourceforge.net/projects/splicebox/files/?source=navbar
      参谋文献:https://academic.oup.com/nar/article/44/10/e98/2516329/CLASS2-accurate-and-efficient-splice-variant

    新葡亰496net,TopHat第三回被刊登已然是7年前,STA奥迪Q5的比对速度是TopHat的50倍,HISAT更是STAENVISION的1.2倍。HISAT2是TopHat2/Bowti2的继任者,使用改善的BWT算法,达成了更加快的快慢和越来越少的资源占用,作者推荐TopHat2/Bowti2和HISAT的客商改变成HISAT2。
    官网:https://ccb.jhu.edu/software/hisat2/index.shtml(学习三个软件最佳的方法就是组成现有中文资料,加上阅读官方表明和HELP文书档案,平常刚开始学习的时候先使用暗许参数,不要乱调参数)

    hisat2的用法

    2. HISAT2的使用

    怎么要用index?官方网站有描述:为了用全套index代表全体基因组,HISAT2 用小的index覆盖了全副基因组,每一个index覆盖了56 Kbp的约束,覆盖整个人类基因组须求55,000 indexes,那一个index结合别的计谋能够快捷准确的比对类别。

    #index 下载
    nohup wget ftp://ftp.ccb.jhu.edu/pub/infphilo/hisat2/data/hg19.tar.gz &
    tar -zxvf *.tar.gz #解压
    # 删除压缩包
    rm -rf *.tar.gz
    
    • hisat2 -h查看援救文件:
    Usage: 
      hisat [options]* -x <bt2-idx> {-1 <m1> -2 <m2> | -U <r> | --sra-acc <SRA accession number>} [-S <sam>]
    
      <bt2-idx>  Index filename prefix (minus trailing .X.ht2).
      <m1>       Files with #1 mates, paired with files in <m2>.
                 Could be gzip'ed (extension: .gz) or bzip2'ed (extension: .bz2).
      <m2>       Files with #2 mates, paired with files in <m1>.
                 Could be gzip'ed (extension: .gz) or bzip2'ed (extension: .bz2).
      <r>        Files with unpaired reads.
                 Could be gzip'ed (extension: .gz) or bzip2'ed (extension: .bz2).
      <SRA accession number>        Comma-separated list of SRA accession numbers, e.g. --sra-acc SRR353653,SRR353654.
      <sam>      File for SAM output (default: stdout)
    
      <m1>, <m2>, <r> can be comma-separated lists (no whitespace) and can be  specified many times.  E.g. '-U file1.fq,file2.fq -U file3.fq'.
    
    • 参数表明:

    -x 指定index文件名
    -1 <m1> 双端测序第二个文件
    -2 <m2> 双端测序第二个文件
    -U 单端测序文件
    --sra-acc 登录号
    新葡亰496net:队列比对,转录组入门。-S 输出文件为sam格式
    -q 输入文件为FASTQ .fq/.fastq格式

    -比对到参照他事他说加以考查基因组,中华VNA-Seq数据是从 S福睿斯奥德赛3589957 ~ SOdyssey库罗德3589962,6个样本的PE数据,也便是有6次巡回,能够写剧本,也可以直接在命令行里运营

    for i in `seq 56 62`
    do
        hisat2 -t -x /opt/NfsDir/UserDir/qin/qin/Data/annotation/hg19/genome -1 /opt/NfsDir/UserDir/qin/qin/Data/RNAseq/SRR35899${i}.sra_1.fastq.gz -2 /opt/NfsDir/UserDir/qin/qin/Data/RNAseq/SRR35899${i}.sra_2.fastq.gz -S /opt/NfsDir/UserDir/qin/qin/Data/RNAseq/align/SRR35899${i}.sam &
    done
    

    运行时刻如下:

    新葡亰496net 1

    Paste_Image.png

    • 很耗CPU啊!用的服务器!

    新葡亰496net 2

    Paste_Image.png

    • 结果:

    新葡亰496net 3

    Paste_Image.png

    SAM(sequence Alignment/mapping卡塔尔数据格式是日前高通量测序中寄放比对数据的正经八百格式,当然她能够用来存放未比对的数额。如今管理SAM格式的工具关键是SAMTools。samtools功能多多,在本次作业中,大家最首要学会将sam文件转变为bam文件,并对bam文件实行sorted(在那之中有三种排序方式N和P),最终创设目录。

    Program: samtools (Tools for alignments in the SAM format)
    Version: 1.3.1-58-gcbee45e (using htslib 1.3.2-228-g0c32631)
    
    Usage:   samtools <command> [options]
    
    Commands:
      -- Indexing
         dict           create a sequence dictionary file
         faidx          index/extract FASTA
         index          index alignment
    
      -- Editing
         calmd          recalculate MD/NM tags and '=' bases
         fixmate        fix mate information
         reheader       replace BAM header
         rmdup          remove PCR duplicates #移除PCR产生的重复序列
         targetcut      cut fosmid regions (for fosmid pool only)
         addreplacerg   adds or replaces RG tags
    
      -- File operations
         collate        shuffle and group alignments by name
         cat            concatenate BAMs
         merge          merge sorted alignments
         mpileup        multi-way pileup
         sort           sort alignment file
         split          splits a file by read group
         quickcheck     quickly check if SAM/BAM/CRAM file appears intact
         fastq          converts a BAM to a FASTQ #格式转换
         fasta          converts a BAM to a FASTA
    
      -- Statistics
         bedcov         read depth per BED region #bed文件的测序深度
         depth          compute the depth 
         flagstat       simple stats
         idxstats       BAM index stats
         phase          phase heterozygotes
         stats          generate stats (former bamcheck)
    
      -- Viewing
         flags          explain BAM flags
         tview          text alignment viewer
         view           SAM<->BAM<->CRAM conversion
         depad          convert padded BAM to unpadded BAM
    

    最首要功效:
    view: BAM-SAM/SAM-BAM 调换和领取部分比对
    sort: 比对排序
    merge: 聚合多少个排序比对
    index: 索引排序比对
    faidx: 创设FASTA索引,提取部分体系
    tview: 文本格式查看系列
    pileup: 产生基于地点的结果和 consensus/indel calling

    # 首先将比对后的sam文件转换成bam文件
    # 利用的是samtools的view选项,参数-S 输入sam文件;参数-b 指定输出的文件为bam;最后重定向写入bam文件
    $ for ((i=56;i<=62;i  ));do samtools view -S /opt/NfsDir/UserDir/qin/qin/Data/RNAseq/align/SRR35899${i}.sam -b > /opt/NfsDir/UserDir/qin/qin/Data/RNAseq/align/SRR35899${i}.bam;done
    # 将所有的bam文件进行排序
    $ nohup for (( i=58;i<=62;i   )); do samtools sort /opt/NfsDir/UserDir/qin/qin/Data/RNAseq/align/SRR35899${i}.bam -o /opt/NfsDir/UserDir/qin/qin/Data/RNAseq/align/SRR35899${i}.sorted.bam;done
    # 将所有的排序文件建立索引,索引文件.bai后缀
    $ for ((i=56;i<=62;i  ));do samtools index /opt/NfsDir/UserDir/qin/qin/Data/RNAseq/align/SRR35899${i}.sorted.bam;done
    
    合在一块跑:
    for i in `seq 56 58`
    do
        samtools view -S SRR35899${i}.sam -b > SRR35899${i}.bam
        samtools sort SRR35899${i}.bam -o SRR35899${i}_sorted.bam
        samtools index SRR35899${i}_sorted.bam
    done
    

    新葡亰496net 4

    Paste_Image.png

    直接去hisat2的主页下载index文件就可以,然后把fastq格式的reads比对上去取得sam文件。 接着用samtools把它转为bam文件,而且排序(注意N和P三种排序差别卡塔尔索引好,载入IGV,再截图多少个基因看看! 顺便对bam文件进行简短QC

    本作业是比对到基因组,所以使用gapped or splices mapper,此流程已经更新。TopHat第三遍被刊登已然是7年前,STAPAJERO的比对速度是TopHat的50倍,HISAT更是STA昂Cora的1.2倍。HISAT2是TopHat2/Bowti2的继任者,使用改正的BWT算法,完成了越来越快的进程和越来越少的能源占用,小编推荐TopHat2/Bowti2和HISAT的客商更动来HISAT2。

    比对质量控制(QC)

    常用工具备

    • Picard https://broadinstitute.github.io/picard/
    • RSeQC http://rseqc.sourceforge.net/
    • Qualimap http://qualimap.bioinfo.cipf.es/
    java -jar /opt/NfsDir/BioDir/picard-tools-1.119/picard.jar CollectMultipleMetrics 
          I=/opt/NfsDir/UserDir/qin/qin/Data/RNAseq/align/SRR3589956.sorted.bam 
          O=/opt/NfsDir/UserDir/qin/qin/Data/RNAseq/align/multiple_metrics 
          R=/opt/NfsDir/PublicDir/reference/ucsc.hg19.fasta 
    

    统计bam文件:

    /opt/NfsDir/BioDir/RSeQC/RSeQC-2.6.4/scripts/bam_stat.py -i /opt/NfsDir/UserDir/qin/qin/Data/RNAseq/align/SRR3589956.sorted.bam
    

    升迁出错:

    新葡亰496net 5

    Paste_Image.png

    反省开掘是路线和称号写错了,修改后就足以了
    对bam文件进行计算解析

    新葡亰496net 6

    Paste_Image.png

    #下载hg19_RefSeq.bed文件
    https://sourceforge.net/projects/rseqc/files/BED/Human_Homo_sapiens/hg19_RefSeq.bed.gz/download
    #查看基因组覆盖率
    /opt/NfsDir/BioDir/RSeQC/RSeQC-2.6.4/scripts/read_distribution.py -i /opt/NfsDir/UserDir/qin/qin/Data/RNAseq/align/SRR3589956.sorted.bam -r /opt/NfsDir/UserDir/qin/qin/Data/annotation/hg19/hg19_RefSeq.bed
    

    新葡亰496net 7

    Paste_Image.png

    1. 下载index文件(小鼠)
      axel ftp://ftp.ccb.jhu.edu/pub/infphilo/hisat2/data/mm10.tar.gz
      tar -xvzf mm10.tar.gz

      下载注释文件
      axel ftp://ftp.sanger.ac.uk/pub/gencode/Gencode_mouse/release_M10/gencode.vM10.annotation.gtf.gz

    2. 队列比对
      a. 配置 1核8G内存
      hisat2 -t -x ref/mm10/genome -1 rawdata/SRR3589960.sra_1.fastq.gz -2 rawdata/SRR3589960.sra_2.fastq.gz -S align/SRR3589960.sam

    官网:

    ![](https://upload-images.jianshu.io/upload_images/10354448-7e1d4e5f0a826b5d.jpg)
    
    111111111111.jpg
    
    
    
    运行时间38分钟,源文件1.6G
    

    下载index文件

    更改配置4核32G
    for ((i=60;i<=62;i=i ));do hisat2 -t -x ref/mm10/genome -1 rawdata/SRR35899${i}.sra_1.fastq.gz -2 rawdata/SRR35899${i}.sra_2.fastq.gz -S align/SRR35899${i}.sam;done

    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
    

    疏解: -t 记录时间 -x hg19(indexState of Qatar文件路径 -1 -2 测序的五个fastq文件 -S 比对结果输出路线 -U 单端测序文件
    reference=~/wikiwei/human/ref/index/hg19/genome
    hisat2 -t -x $reference -U SRR957679.fastq -S siSUZ12_1.sam 2>siSUZ12_1.log

    -c:断点续传

    参照文章

    比对、排序、索引

    1. 转录组入门(mac版本)
    2. https://www.cnblogs.com/freescience/p/7342895.html

    把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,HavalSeQC,Qualimap还会有一大堆

    本文由新葡亰496net发布于电脑系统,转载请注明出处:新葡亰496net:队列比对,转录组入门

    关键词: