deeptools使用过程中出现问题集

1想要画一个chip-seq的热图，首先要通过bamcompare获得针对IP和input样品的处理后的bigwig文件，

通过运行 命令

bamcompare -b1 IP bam -b2 inputbam -o log2ratio_IP1bw，

结果报错

“IP bam does not appear to have an index  You must index the file first!”

在网上找了下解决方法是运行命令：

samtools index IP bam

后，重新运行上述bamcompare命令

结果可行。

注：input bam 也需要进行同样的处理。其实这个index文件也就是测序公司给的bambai文件。由于我这个分析是在服务器上运行的，上传的时候没有把bai文件也传上去，所以出现了这个问题。

2前面没有考虑到生物学重复的问题，如果ChIP有两个生物学重复，是不是要合并两个生物学重复后作图？这里就涉及了将两个生物学重复的BAM文件进行merge，需要使用samtools的merge命令：

Usage: samtools merge [-nr] [-h inhsam] <outbam> <in1bam> <in2bam>[]

Options: -n      sort by read names

        -r      attach RG tag (inferred from file names)

        -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 <outbam> [in1bam]

1、sort
sort对bam文件进行排序。
Usage: samtools sort [-n] [-m <maxMem>] <inbam> <outprefix>
-m 参数默认下是 500,000,000 即500M（不支持K，M，G等缩写）。对于处理大数据时，如果内存够用，则设置大点的值，以节约时间。
-n 设定排序方式按short reads的ID排序。默认下是按序列在fasta文件中的顺序（即header）和序列从左往右的位点排序。
例子：
$ samtools sort abcbam abcsort
$ samtools view abcsortbam | less -S
2、
samtools还有个非常重要的命令mpileup，以前为pileup。该命令用于生成bcf文件，再使用bcftools进行SNP和Indel的分析。bcftools是samtool中附带的软件，在samtools的安装文件夹中可以找到。
最常用的参数有2： -f 来输入有索引文件的fasta参考序列； -g 输出到bcf格式。用法和最简单的例子如下
Usage: samtools mpileup [-EBug] [-C capQcoef] [-r reg] [-f infa] [-l list] [-M capMapQ] [-Q minBaseQ] [-q minMapQ] inbam [in2bam []]
$ samtools mpileup -f genomefasta abcbam > abctxt
$ samtools mpileup -gSDf genomefasta abcbam > abcbcf
$ samtools mpileup -guSDf genomefasta abcbam | \
bcftools view -cvNg - > abcvcf
mpileup不使用-u或-g参数时，则不生成二进制的bcf文件，而生成一个文本文件(输出到标准输出)。该文本文件统计了参考序列中每个碱基位点的比对情况；该文件每一行代表了参考序列中某一个碱基位点的比对结果。比如：
scaffold_1 2841 A 11 ,,,, BHIGDGIJFF
scaffold_1 2842 C 12 ,$,,,^I CFGEGEGGCFF+
scaffold_1 2843 G 11 ,,, FDDDDCDDD+
scaffold_1 2844 G 11 ,,, FAAAAA<AA+
scaffold_1 2845 G 11 ,,, F656666166
scaffold_1 2846 A 11 ,,, (11111)11
scaffold_1 2847 A 11 ,,+9acggtgaag+9ACGGTGAAT+9ACGGTGAAG+9ACGGTGAAG,+9acggtgaag+9ACGGTGAAG+9ACGGTGAAG+9ACGGTGAAG+9ACGGTGAAG+9ACGGTGAAG %+-)
scaffold_1 2848 N 11 agGGGgGGGGG !!$!!!!!!!!
scaffold_1 2849 A 11 c$,, !0000000000
scaffold_1 2850 A 10 ,, 353333333
mpileup生成的结果包含6行：参考序列名；位置；参考碱基；比对上的reads数；比对情况；比对上的碱基的质量。其中第5列比较复杂,解释如下：
1 ‘’代表与参考序列正链匹配。
2 ‘,’代表与参考序列负链匹配。
3 ‘ATCGN’代表在正链上的不匹配。
4 ‘atcgn’代表在负链上的不匹配。
5 ‘’代表模糊碱基
6 ‘^’代表匹配的碱基是一个read的开始；’^'后面紧跟的ascii码减去33代表比对质量；这两个符号修饰的是后面的碱基，其后紧跟的碱基(,ATCGatcgNn)代表该read的第一个碱基。
7 ‘$’代表一个read的结束，该符号修饰的是其前面的碱基。
8 正则式’\+[0-9]+[ACGTNacgtn]+’代表在该位点后插入的碱基；比如上例中在scaffold_1的2847后插入了9个长度的碱基acggtgaag。表明此处极可能是indel。
9 正则式’-[0-9]+[ACGTNacgtn]+’代表在该位点后缺失的碱基；

---