【ChIP-seq 实战】五、合并bam和去除PCR重复
2022-08-11 本文已影响0人
佳奥
这里是佳奥!
获得了比对的bam文件,生成了index的bai文件,生成了比对结果stat文件。
##可以先multiqc看一下结果
multiqc ./
合并bam文件和去除PCR重复没有严格先后顺序,是否必要根据文章决定
1 合并bam文件
需要回文章看
是同一个样本的多个测序、还是多个生物学重复。(一个样本测了两次,还是两次实验同一个样本)
因为一个样品分成了多个lane进行测序,所以在进行peaks calling的时候,需要把bam进行合并。
## 如果不用循环 格式:输出 两个输入
samtools merge control.merge.bam Control_1_trimmed.bam Control_2_trimmed.bam
## 通常我们用批处理
mkdir mergeBam
cd /home/kaoku/chipseq/mouse_project/align
ls *.bam| sed 's/_[0-9]_trimmed.bam//g' |sort -u | while read id; do samtools merge ../mergeBam/$id.merge.bam $id*.bam ; done
##合并以后
(chipseq) root 16:05:07 /home/kaoku/chipseq/mouse_project/mergeBam
$ ls -lh
总用量 9.6G
-rw-r--r-- 1 root root 836M 8月 11 15:57 Control.merge.bam
-rw-r--r-- 1 root root 1.3G 8月 11 15:58 H2Aub1.merge.bam
-rw-r--r-- 1 root root 1.6G 8月 11 15:59 H3K36me3.merge.bam
-rw-r--r-- 1 root root 1.3G 8月 11 16:00 Ring1B.merge.bam
-rw-r--r-- 1 root root 1.1G 8月 11 16:01 RNAPII_8WG16.merge.bam
-rw-r--r-- 1 root root 1.5G 8月 11 16:02 RNAPII_S2P.merge.bam
-rw-r--r-- 1 root root 1.4G 8月 11 16:03 RNAPII_S5P.merge.bam
-rw-r--r-- 1 root root 215M 8月 11 16:04 RNAPII_S5PRepeat.merge.bam
-rw-r--r-- 1 root root 713M 8月 11 16:04 RNAPII_S7P.merge.bam
2 去除PCR重复
##使用软件samtools(后台运行,top查看后台)
ls *merge.bam | while read id ; do ( samtools markdup -r $id $(basename $id ".bam").rmdup.bam & );done
##把去除重复的bam文件再次建立索引,查看比对效果
ls *.rmdup.bam | xargs -i samtools index {}
ls *.rmdup.bam | while read id ; do ( samtools flagstat $id > $(basename $id ".bam").stat & );done
##结果如下,去除重复前后比较
(chipseq) root 20:58:44 /home/kaoku/chipseq/mouse_project/mergeBam
$ ls -lh
总用量 18G
-rw-r--r-- 1 root root 836M 8月 11 15:57 Control.merge.bam
-rw-r--r-- 1 root root 743M 8月 11 20:57 Control.merge.rmdup.bam
-rw-r--r-- 1 root root 1.3G 8月 11 15:58 H2Aub1.merge.bam
-rw-r--r-- 1 root root 1.1G 8月 11 20:58 H2Aub1.merge.rmdup.bam
-rw-r--r-- 1 root root 1.6G 8月 11 15:59 H3K36me3.merge.bam
-rw-r--r-- 1 root root 1.5G 8月 11 20:58 H3K36me3.merge.rmdup.bam
-rw-r--r-- 1 root root 1.3G 8月 11 16:00 Ring1B.merge.bam
-rw-r--r-- 1 root root 1006M 8月 11 20:58 Ring1B.merge.rmdup.bam
-rw-r--r-- 1 root root 1.1G 8月 11 16:01 RNAPII_8WG16.merge.bam
-rw-r--r-- 1 root root 984M 8月 11 20:58 RNAPII_8WG16.merge.rmdup.bam
-rw-r--r-- 1 root root 1.5G 8月 11 16:02 RNAPII_S2P.merge.bam
-rw-r--r-- 1 root root 1.2G 8月 11 20:58 RNAPII_S2P.merge.rmdup.bam
-rw-r--r-- 1 root root 1.4G 8月 11 16:03 RNAPII_S5P.merge.bam
-rw-r--r-- 1 root root 775M 8月 11 20:58 RNAPII_S5P.merge.rmdup.bam
-rw-r--r-- 1 root root 215M 8月 11 16:04 RNAPII_S5PRepeat.merge.bam
-rw-r--r-- 1 root root 210M 8月 11 20:57 RNAPII_S5PRepeat.merge.rmdup.bam
-rw-r--r-- 1 root root 713M 8月 11 16:04 RNAPII_S7P.merge.bam
-rw-r--r-- 1 root root 610M 8月 11 20:57 RNAPII_S7P.merge.rmdup.bam
##查看一下比对成功率
$ grep 'N/A' *.stat | grep '%'
Control.merge.rmdup.stat:12330969 + 0 mapped (85.16% : N/A)
H2Aub1.merge.rmdup.stat:17516222 + 0 mapped (96.82% : N/A)
H3K36me3.merge.rmdup.stat:22685679 + 0 mapped (98.51% : N/A)
Ring1B.merge.rmdup.stat:24901367 + 0 mapped (93.46% : N/A)
RNAPII_8WG16.merge.rmdup.stat:23397509 + 0 mapped (94.84% : N/A)
RNAPII_S2P.merge.rmdup.stat:26655659 + 0 mapped (95.36% : N/A)
RNAPII_S5P.merge.rmdup.stat:13680963 + 0 mapped (90.78% : N/A)
RNAPII_S5PRepeat.merge.rmdup.stat:3997567 + 0 mapped (82.22% : N/A)
RNAPII_S7P.merge.rmdup.stat:9759486 + 0 mapped (77.96% : N/A)
RNAPII_S7P.merge.rmdup.stat:9759486 + 0 primary mapped (77.96% : N/A)
去除PCR重复和不去除PCR重复的样本都找一次peaks看一下。
我们拿到两批文件:合并的bam文件,去除PCR重复的合并bam文件。
下一步就是使用macs2寻找peaks了!
我们下一篇再见!