chipseq生信

【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了!

我们下一篇再见!

上一篇下一篇

猜你喜欢

热点阅读