samtools 命令详解
samtools 常用命令总结
samtools
view
重要参数释义：
-b：输出bam格式，用于后续分析
-C：输出CRAM文件
-1：快速压缩（需要-b）
-u：输出未压缩的bam文件，节约时间，占据较多磁盘空间（需要-b）
-h：默认输出sam文件不带表头，该参数设定后输出带表头信息
-H：仅仅输出表头信息
-c：仅打印匹配数
-o：输出文件（stdout标准输出）
-U：输出没有经过过滤选择的reads
-t：制表分隔符文件（需要提供额外的参考数据，比如参考基因组的索引）
-L：仅包括和bed文件有重叠的reads
-r：仅输出在STR读段组中的reads
-R：仅输出特定reads
-q：定位的质量大于INT[默认0]
-l：仅输出在STR 库中的reads
-F：获得比对上（mapped）的过滤设置[默认0]
-f：获得未必对上（unmapped）的过滤设置[默认0]
-T：使用fasta格式的参考序列
例子如下
- bam文件转换为sam文件
samtools view -h smallNA06985 > test.sam
- sam文件转换为bam文件
# converting sam to bam and sort as well
ls *.sam | while read i
do (samtools sort -O bam -@ 10 -o $(basename ${i} ".sam").bam ${i})
done
# building index
ls *.bam | xargs -i samtools index {}
ls *.bam | while read i; do samtools flagstat -@ 10 $i > $(basename ${i} ".bam").flagstat; done
- 提取比对到参考基因组上的数据
samtools view -bF 4 test.bam > test.F.bam
- 提取没有比对到参考基因组上的数据
samtools view -bf 4 test.bam > test.f.bam
- 双端reads都比对到参考基因组上的数据
samtools view -bF 12 test.bam > test.12.bam
## 提取paired reads中两条reads都比对到参考序列上的比对结果，只需要把两个4+8的值12作为过滤参数即可，但其实我们双端测序的数据是需要取PE map到同一个fragment的reads
samtools view -bf 2 test.bam > test.f2.bam
- 单端测序筛掉低质量的数据
samtools view -bF 4 abc.bam > abc.F.bam
- 单端reads1比对到参考基因组上的数据
samtools view -bF 4 -f 8 test .bam > test1.bam
- 单端reads2比对到参考基因组上的数据
samtools view -bF 8 -f 4 test.bam > test2.bam
merge
bam文件的header section用不同的tag表示不同的信息，主要有@HD，说明符合标准的版本、对比序列的排列顺序；@SQ，参考序列说明；@RG，比对上的序列（read）说明；@PG，使用的程序说明；@CO，任意的说明信息
merge命令将2个或2个以上的已经sort了的bam文件融合成一个bam文件。融合后的文件不需要则是已经sort过了的。
Usage: samtools merge [-nr] [-h inh.sam] <out.bam> <in1.bam> <in2.bam>[...]
Options:
-n sort by read names
-r attach RG tag (inferred from file names) #@RG，比对上的序列（read）说明
-u uncompressed BAM output
-f overwrite the output BAM if exist
-1 compress level 1
-R STR merge file in the specified region STR [all]
-h FILE copy the header in FILE to <out.bam> [in1.bam] #指定FILE内的’@’头复制到输出bam文件中并替换输出文件的文件头。否则，输出文件的文件头将从第一个输入文件复制过来；
Note: Samtools' merge does not reconstruct the @RG dictionary in the header. Users must provide the correct header with -h, or uses Picard which properly maintains the header dictionary in merging.
ls *.rmdup.bam| cut -d "_" -f 1-2 | sort -u | while read i; do (samtools merge $i.merge.bam $i*.rmdup.bam); done
待合并的bam文件，必须有与其对应的index文件。
samtools merge [-nur1f] [-h inh.sam] [-R reg] [-b <list>] <out.bam> <in1.bam> [<in2.bam><in3.bam>…<inN.bam]
-c 当多个输入文件包含相同的@RG头ID时，只保留第一个到合并后输出的文件。当合并多个相同样本的不同文件时，非常有用。
-p 与-c参数类似，对于要合并的每一个文件中的@PG ID只保留第一个文件中的@PG。
合并test_L1.bam和test_L2.bam文件
samtools merge -h test.sam
test_L1_L2.bam
test_L1.sorted.bam
test_L2.sroted.bam
合并test_L1.bam和test_L2.bam文件中的指定区域chr7
samtools merge -h test.sam
-R chr7
test_L1_L2.bam
test_L1.sorted.bam
test_L2.sroted.bam
在实际操作中，我是使用yz教我的办法，在刚开始的时候就用cat命令，将两个测序文件merge到一起
ls
HiC_rep1_run1_1.fastq.gz HiC_rep1_run3_1.fastq.gz HiC_rep2_run1_1.fastq.gz HiC_rep2_run3_1.fastq.gz
HiC_rep1_run1_2.fastq.gz HiC_rep1_run3_2.fastq.gz HiC_rep2_run1_2.fastq.gz HiC_rep2_run3_2.fastq.gz
HiC_rep1_run2_1.fastq.gz HiC_rep1_run4_1.fastq.gz HiC_rep2_run2_1.fastq.gz
HiC_rep1_run2_2.fastq.gz HiC_rep1_run4_2.fastq.gz HiC_rep2_run2_2.fastq.g
ls HiC*gz | cut -d "_" -f 1-2 | sort -u | while read i; do zcat ${i}_*_1.fastq.gz | gzip -c > ${i}_1.fastq.gz; zcat ${i}_*_2.fastq.gz | gzip -c > ${i}_2.fastq.gz; done
