[samtools] sam格式与bam格式互换,提取未匹配reads,转为fastq

Converting a SAM file to a BAM file

First, if you use the Unix command

head test.sam

The first 10 lines on your terminal after typing "head test.sam", should be lines starting with an "@" sign, which is an indicator for a header line. If you don't see lines starting with the "@" sign, the header information is most likely missing.

If the header information is absent不存在 from the SAM file use the command below, where reference.fa is the reference fasta file used to map the reads:

samtools view -bT reference.fa test.sam > test.bam

If the header information is available:

samtools view -bS test.sam > test.bam

Sorting a BAM file

samtools sort test.bam -o test_sorted

Creating a BAM index file

samtools index test_sorted.bam test_sorted.bai

Converting a BAM file to a SAM file

Note: remember to use -h to ensure the converted SAM file contains the header information. Generally, I suggest storing only sorted BAM files as they use much less disk space and are faster to process.

samtools view -h NA06984.chrom16.ILLUMINA.bwa.CEU.low_coverage.20100517.bam > NA06984.chrom16.ILLUMINA.bwa.CEU.low_coverage.20100517.sam

Simple stats

samtools flagstat NA06984.chrom16.ILLUMINA.bwa.CEU.low_coverage.20100517.bam

10182494 in total
0 QC failure
223627 duplicates
9861117 mapped (96.84%)
10095646 paired in sequencing
5049066 read1
5046580 read2
8174084 properly paired (80.97%)
9452892 with itself and mate mapped
321377 singletons (3.18%)
215316 with mate mapped to a different chr
126768 with mate mapped to a different chr (mapQ>=5)
For more statistics of SAM or BAM files have a look at the SAMStat program.

Interpreting the BAM flags

(也可以利用FLAG值含义解释工具:https://www.plob.org/2012/02/04/1697.html

Here are some common BAM flags:

163: 10100011 in binary

147: 10010011 in binary

99: 1100011 in binary

83: 1010011 in binary

Interpretation解释 of 10100011 (reading the binary from left to right):

1the read is paired in sequencing, no matter whether it is mapped in a pair
1the read is mapped in a proper pair (depends on the protocol, normally inferred during alignment)
0the query sequence itself is unmapped
0the mate is unmapped
0strand of the query (0 for forward; 1 for reverse strand)
1strand of the mate
0the read is the first read in a pair
1the read is the second read in a pair

163 second read of a pair on the positive strand with negative strand mate

147 second read of a pair on the negative strand with positive strand mate

99 first read of a pair on the forward strand with negative strand mate

83 first read of a pair on the reverse strand with positive strand mate

Extracting only the first read from paired end BAM files

samtools view -h -f 0x0040 test.bam > test_first_pair.sam

0x0040 is hexadecimal十六进制 for 64 (i.e. 16 * 4), which is binary for 1000000, corresponding to the read in the first read pair.

Filtering out unmapped reads in BAM files

samtools view -h -F 4 blah.bam > blah_only_mapped.sam

Creating FASTQ files from a BAM file

I found this great tool at http://www.hudsonalpha.org/gsl/software/bam2fastq.php

For example to extract ONLY unaligned from a bam file:

bam2fastq -o blah_unaligned.fastq --no-aligned blah.bam

To extract ONLY aligned reads from a bam file:

bam2fastq -o blah_aligned.fastq --no-unaligned blah.bam

转载于:https://www.cnblogs.com/xiaofeiIDO/p/6424649.html

评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值