bwa samtools bowtie2大杂烩

作者: 11的雾 | 来源:发表于2018-09-27 09:32 被阅读53次

    比对,先要建index

    bwa index /cygene/database/TCR/TRAVBV.seq.fa
    一般bwa比对完后需要用samtools转换成bam或者排序什么的,都需要用samtools首先对ref进行建index
    samtools faidx TRAVBV.seq.fa
    基础参数介绍

    -R STR        read group header line such as '@RG\tID:foo\tSM:bar'  [null]
    

    如: "@RG\tID:0\tPL:BGISEQ-500\tPU:0\tLB:UN\tSM:0\tCN:BGI"

    bwa mem 搭配samtools比对,及sort排序:

    PE(paired end)

    bwa mem -t 18 /cygene/database/TCR/TCR.all.vdjab.fa /cygene/data/ANCCR180494_PM-CR180494-04_AHNTL2CCXY_2018-10-01/Cleandata/RP01G9E1L3/RP01G9E1L3_R1.fq.gz /cygene/data/ANCCR180494_PM-CR180494-04_AHNTL2CCXY_2018-10-01/Cleandata/RP01G9E1L3/RP01G9E1L3_R2.fq.gz |samtools view -buhS -t /cygene/database/TCR/TCR.all.vdjab.fa.fai - |samtools sort -n -o Cleandata_sample1_Rd12.vdjab.bam -

    SE (single end)

    bwa mem -t 18 /cygene/database/TCR/TCR.all.vdjab.fa /cygene/data/ANCCR180494_PM-CR180494-03_000000000-C2R58_2018-09-14/Rawdata/RP01G9E1L1/RP01G9E1L1_R1.fq.gz |samtools view -buhS -t /cygene/database/TCR/TCR.all.vdjab.fa.fai - |samtools sort -n -o sample1_R1.vdjab.bam -

    bwa aln

    bwa aln用于短read,
    PE:比对两次,再合并
    SE:

    比对结果处理

    比对结果中没有比对上的flag为4
    查看bam/sam文件用 samtools view xxx.bam |le
    过滤没有比对上的结果 用-F:

    -F INT   only include reads with none of the FLAGS in INT present [0]
    

    例子1:从bam中过滤掉没有比对上的信息,并将比对上的部分保存到sam中,
    samtools view -F 4 xxx.bam >xxx.mapped.sam
    例子2:从bam中过滤没有比对上的信息,并保存到bam中:(要加-h,表示输出header)
    samtools view -F 4 -h xxx.bam |samtools view -h -o xxx.mapped.bam -
    注意:虽然有些没有比对上的结果 的flag值是141,77等,但是都可以用4过滤掉。

    image.png

    ========================================
    提取为fastq:主要参数为-f,-F,-G三个。
    提取没有比对上的信息到fastq用:
    samtools fastq -f 4 xxx.bam > xxx.unmapped.fastq
    提取比对上的reads到fastq用
    samtools fastq -F 4 xxx.bam > xxx.mapped.fastq
    查询flag值的意义:https://broadinstitute.github.io/picard/explain-flags.html

    samtools sort 报错案例:

    使用samtools 为bam建立index之前,要对bam进行sort,而且sort时候不能加-n参数(sort by read name)。否则samtools index会报错:

    [E::hts_idx_push] NO_COOR reads not in a single block at the end 280 -1
    samtools index: failed to create index for "sample1_R1.vdjab.bam"
    

    bam与sam格式相互转换

    BAM 转SAM :samtools view -h -o out.sam out.bam
    SAM转BAM :samtools view -bS out.sam >out.bam
    -b 意思是输出为BAM format
    -S 意思是输入为SAM,如果@SQ 缺剩, 要写-t,所以如果没有@SQ
    samtools faidx ref.fa
    samtools view -bt ref.fa.fai out.sam > out.bam

    bowtie2

    bowtie2建index用:
    bowtie2-build -f TRAVBV.seq.fa TRAVBV.seq前缀就可以
    全局局部比对要加参数--local
    只比对正链或只比对负链用--norc,--nofw
    single end数据用-U比对实例:
    bowtie2 -p 20 -x /cygene/database/TCR/TRAVBV.seq -U /cygene/data/ANCCR180494_PM-CR180494-03_000000000-C2R58_2018-09-14/Rawdata/RP01G9E1L1/RP01G9E1L1_R1.fq.gz -S sample1_r1_bowtie.fw.sam --nofw --local

    相关文章

      网友评论

        本文标题:bwa samtools bowtie2大杂烩

        本文链接:https://www.haomeiwen.com/subject/zpisoftx.html