????說(shuō)到deeptools，大家很容易想到ChIP。deeptools集美貌與才華?是表觀分析的深海利器。關(guān)于deeptools的應(yīng)用，其中最為經(jīng)典的就是plotHeatmap和plotProfiles（如下圖）。幾乎每個(gè)經(jīng)典的表觀文章都能找到他們的身影。
經(jīng)典的用法主要是研究某種組蛋白或者轉(zhuǎn)錄因子是不是在TSS區(qū)域富集（如下圖）。
但殊不知deeptools不僅僅是查看組蛋白或者轉(zhuǎn)錄因子在基因區(qū)域的修飾，還能為三維基因組理論奠定基礎(chǔ)添磚加瓦，例如EP_loop（enhancer_promoter構(gòu)成的loop）的CTCF以及蛋白修飾情況（如下圖）：
例如，TAD邊界以及ABcompartment兩端的組蛋白修飾情況（如下圖）：
說(shuō)到這里，小編的手機(jī)響了，收到 一條信息～～
Hi～～fan，原來(lái)deeptools還可以用作HiC多組學(xué)分析呀，那我想拿TAD邊界的組蛋白修飾和AB compartment組蛋白修飾練練手，我具體該怎么做呢？
哈哈，別急，容我一一道來(lái)～～
首先看圖～～?
TAD邊界這幅圖呢，主要是以一個(gè)點(diǎn)為中心，而ABcompart這幅圖則是有兩個(gè)中心點(diǎn)border1，border2（箭頭所指）
這就涉及deeptools的computeMatrix兩種模式（如下圖）：
deeptools? computeMatrix模塊可以計(jì)算每個(gè)基因區(qū)域的結(jié)合得分，生成中間文件用以給plotHeatmap和plotProfiles作圖。
reference-point mode則是給定一個(gè)bed file，以某個(gè)點(diǎn)為中心開(kāi)始統(tǒng)計(jì)信號(hào)（如TSS/TES/center）。
來(lái)舉個(gè)'??'：
首先，整理一下TAD邊界的bed文件boundaries.bed??
其次，第二步通過(guò)以下命令將bam轉(zhuǎn)換成bigwig：
????bamCoverage --bam? CTCF.bam -o?CTCF.bw
? ? --binSize 10
? ? --numberOfProcessors 4?
? ?--normalizeUsingRPKM
? ? --extendReads
之后，第三步利用reference-point 計(jì)算以這個(gè)點(diǎn)為中心上下游500k范圍內(nèi)以10k為窗口計(jì)算每個(gè)窗口CHIP的信號(hào)強(qiáng)度：
computeMatrix reference-point? --referencePoint center?
? ?????-b 500000 -a 500000?
? ? ? ?-R boundaries.bed
? ? ? ? -S CTCF.bw H3K27ac.bw H3K27me3.bw H3K4me3.bw --skipZeros
? ? ? ? -o? boundaries_modification.gz?
? ? ? ?--binSize 10000
然后，第四步畫(huà)線圖：
plotProfile -m boundaries_modification.gz?
????-out boundaries_modification.pdf?
? ? ? --plotType se?
??????--plotFileFormat 'pdf'?
????--refPointLabel Boundary --numPlotsPerRow 4?
?????--perGroup?
????--plotHeight 8 --plotWidth 8 -
? ? --legendLocation lower-center?
????--yMin 0 --plotTitle ""
?????--yAxisLabel "Normalized signal"
經(jīng)過(guò)以上四步就可以收工了，快來(lái)看結(jié)果，怎么都感覺(jué)有點(diǎn)小清新呢～
scale-regions mode簡(jiǎn)單來(lái)說(shuō)會(huì)將給定bed file范圍內(nèi)的結(jié)合信號(hào)做一個(gè)統(tǒng)計(jì)（指的是一段長(zhǎng)度），并將基因長(zhǎng)度統(tǒng)一scale到設(shè)定regionBdoyLength的長(zhǎng)度，加上統(tǒng)計(jì)基因上游和下游（通過(guò)beforeRegionStartLength參數(shù)和afterRegionStartLength參數(shù)設(shè)定）n bp的信號(hào)。
同樣的針對(duì)ABcompartment先要識(shí)別連續(xù)的ABcomaprtment位點(diǎn)形成chr，start，end三列bed文件(A_com.bed ,B_com.bed)。
之后scale-regions? 計(jì)算當(dāng)前以及延伸區(qū)域的信號(hào)值～
computeMatrix scale-regions?
????--regionsFileName A_com.bed? B_com.bed
?????--scoreFileName? H3K27ac.bw? ?H3K27me3.bw? H3K36me3.bw? ?H3K4me3.bw? ????????????H3K9me3.bw??
?????--skipZeros?
????--regionBodyLength 6000000?
????--beforeRegionStartLength 2000000?
????--afterRegionStartLength 2000000?
????--startLabel border --endLabel border?
????-outFileName? AB_modification.gz??
?????--binSize 100000?
統(tǒng)計(jì)完信號(hào)之后，就可以畫(huà)圖了：
plotProfile? -m? ?AB_modification.gz??
????-out? AB_mod_lineProfile.pdf?
????--plotType? se?
????--startLabel border?
????--endLabel border?
????--plotFileFormat? pdf?
? ?--plotHeight 8 --plotWidth 12
? ? --yAxisLabel "Reads Density"
? ? ?--perGroup
咚咚咚～～，通過(guò)這樣簡(jiǎn)單的步驟，圖就畫(huà)好了，驚不驚喜？意不意外？
艾瑪，就扯到這兒吧～～，朕累了～～
