怎么使用SICER進行peak calling

怎么使用SICER進行peak calling，針對這個問題，這篇文章詳細介紹了相對應的分析和解答，希望可以幫助更多想解決這個問題的小伙伴找到更簡單易行的方法。

chip_seq數據中peak的長度范圍跨度較大，既有覆蓋幾個核小體的幾百bp的peak, 也有包含多個基因長度在上千kb的peak。比如H3K4me2和H3K4me3這兩種組蛋白修飾中peak在幾百bp左右, 而H3K27me3中則為長度在幾十到幾百kb之間。組蛋白修飾中peak長度跨度大，弱信號分散都特點，使得基于轉錄因子TF結合位點的peak calling軟件在分析這類數據時準確度較差。

SICER是一款專門針對組蛋白修飾的chip數據進行peak calling的軟件，核心思想也是基于滑動窗口和局部泊松分布的方式來識別富集區域，下圖所示為該軟件用默認參數識別到的H3K27me3的peak區域

黑色區域為ENCODE分析得到的peak區域，紅色區域為SICER分析得到的peak區域。該軟件官網如下

https://home.gwu.edu/~wpeng/Software.htm

為例方便使用，有人對該軟件進行了分裝，使用起來更加方便，源代碼托管在github上，網址如下

https://github.com/dariober/SICERpy

基本用法如下

python SICERpy \
-c input.bam \
-w 200 \
-g 3 \
-t ip.bam \
> peak.bed

-w參數表示滑動窗口的大小，默認值為200。數值越小， 識別到的peak區間長度相對越短且越分散；數值越大，會造成過渡擬合，識別到的peak區間過長，丟失掉真實的信息，示意如下

對于轉錄因子，官方推薦滑動窗口設置為50-100bp， 對于組蛋白修飾，推薦設置為200bp。

-g參數代表gap的大小，默認值為3。和windows size類似，該參數也直接影響peak區間的定義，示意如下

對于轉錄因子，官方推薦該數值和滑動窗口數值保持相同；對于組蛋白修飾，推薦值為3。

輸出文件為bed格式，共8列，每列含義如下

1. chrom

2. start

3. end

4. chip read count

5. input read count

6. pvalue

7. fold_enrichment

8. fdr

   
   

可以最后一列的fdr值，來篩選得到高可信度的peak信息，用法如下

awk '$8 < 0.01' peaks.bed > peaks.01.bed

關于怎么使用SICER進行peak calling問題的解答就分享到這里了，希望以上內容可以對大家有一定的幫助，如果你還有很多疑惑沒有解開，可以關注億速云行業資訊頻道了解更多相關知識。