這篇文章主要講解了“如何使用ChIPseeker進(jìn)行peak注釋”,文中的講解內(nèi)容簡單清晰,易于學(xué)習(xí)與理解,下面請大家跟著小編的思路慢慢深入,一起來研究和學(xué)習(xí)“如何使用ChIPseeker進(jìn)行peak注釋”吧!
創(chuàng)新互聯(lián)建站長期為上千客戶提供的網(wǎng)站建設(shè)服務(wù),團(tuán)隊(duì)從業(yè)經(jīng)驗(yàn)10年,關(guān)注不同地域、不同群體,并針對不同對象提供差異化的產(chǎn)品和服務(wù);打造開放共贏平臺,與合作伙伴共同營造健康的互聯(lián)網(wǎng)生態(tài)環(huán)境。為宜良企業(yè)提供專業(yè)的成都網(wǎng)站建設(shè)、網(wǎng)站設(shè)計(jì),宜良網(wǎng)站改版等技術(shù)服務(wù)。擁有10年豐富建站經(jīng)驗(yàn)和眾多成功案例,為您定制開發(fā)。
ChIPseeker是使用的最廣泛的peak注釋軟件之一,提供了以下多種功能
peak在染色體和TSS位點(diǎn)附近分布情況可視化
peak關(guān)聯(lián)基因注釋以及在基因組各種元件上的分布
獲取GEO數(shù)據(jù)庫中peak的bed文件
多個peak文件的比較和overlap分析
首先我們需要輸入peak文件,支持兩種格式,第一種是BED格式,最少只需要3列內(nèi)容記錄peak的染色體位置就可以了,示意如下
當(dāng)然也可以有多余的列,只需要符合BED格式的標(biāo)準(zhǔn)即可;另外一種和MACS的peak calling輸出結(jié)果類似,第一行為表頭,示意如下
通過函數(shù)readPeaks
讀取peak文件,用法如下
peak <- readPeakFile("peak.bed")
函數(shù)根據(jù)文件名稱的后綴來判斷是否為bed格式,建議BED格式的輸入文件后綴統(tǒng)一成.bed
, 當(dāng)然壓縮文件也是支持的,比如.bed.gz
;如果不是BED格式的輸入,文件名稱則不能使用BED格式對應(yīng)的后綴。
下面來詳細(xì)看下幾個主要功能的代碼和結(jié)果展示
用法如下
covplot(peak, chr = c("chr1", "chr2"))
輸出結(jié)果示意如下
用法如下
library(TxDb.Hsapiens.UCSC.hg19.knownGene)
txdb <- TxDb.Hsapiens.UCSC.hg19.knownGene
# 定義TSS上下游的距離
promoter <- getPromoters(TxDb=txdb, upstream=3000, downstream=3000)
tagMatrix <- getTagMatrix(peak, windows=promoter)
tagHeatmap(tagMatrix, xlim=c(-3000, 3000), color="red")
熱圖每一行代表一個基因,展示的是所有基因TSS兩側(cè)的分布,除了熱圖外,還可以對所有基因取均值,用折線圖來展示TSS兩側(cè)分布情況,用法如下
plotAvgProf(
tagMatrix,
xlim=c(-3000, 3000),
xlab="Genomic Region (5'->3')",
ylab = "Read Count Frequency")
輸出結(jié)果示意如下
用法如下
peakAnno <- annotatePeak(
peak,
tssRegion = c(-3000, 3000),
TxDb = txdb,
annoDb = "org.Hs.eg.db")
write.table(
as.data.frame(peakAnno),
"peak.annotation.tsv",
sep="\t",
row.names = F,
quote = F)
注釋文件內(nèi)容如下
給出了關(guān)聯(lián)的基因以及對應(yīng)的基因組區(qū)域的類別,根據(jù)這個結(jié)果,可以提取關(guān)聯(lián)基因進(jìn)行下游的功能富集分析,比如提取geneid
這一列,用clusterProfiler
進(jìn)行GO/KEGG等功能富集分析。
注釋的結(jié)果還提供了多種可視化方式,其中餅圖最為常見,用法如下
plotAnnoPie(peakAnno)
輸出結(jié)果示意如下
以hg19為例,首先查詢對應(yīng)的GEO編號信息,用法如下
> hg19 <- getGEOInfo(genome="hg19", simplify=TRUE)
> head(hg19)
series_id gsm organism
111 GSE16256 GSM521889 Homo sapiens
112 GSE16256 GSM521887 Homo sapiens
113 GSE16256 GSM521883 Homo sapiens
114 GSE16256 GSM1010966 Homo sapiens
115 GSE16256 GSM896166 Homo sapiens
116 GSE16256 GSM910577 Homo sapiens
由于列數(shù)太多,上述結(jié)果只展示了部分信息,對于每個bed文件,會列出對應(yīng)的描述信息,方便篩選感興趣的peak進(jìn)行下載,可以根據(jù)GSM
編號進(jìn)行下載,用法如下
downloadGSMbedFiles("GSM521889", destDir="hg19")
也可以根據(jù)下載一個基因組對應(yīng)的所有peak文件,用法如下
downloadGEObedFiles(genome="hg19", destDir="hg19")
peak的overlap分析不僅可以探究生物學(xué)重復(fù)樣本間的一致性,還可以進(jìn)一步識別多種蛋白或者轉(zhuǎn)錄因子在調(diào)控網(wǎng)絡(luò)中的作用,如果兩個蛋白的chip結(jié)果overlap顯著,很可能這兩個蛋白構(gòu)成了復(fù)合體,或者兩種蛋白具有相互作用,這對于探究其調(diào)控機(jī)制有相當(dāng)大的幫助。用法如下
enrichPeakOverlap(
queryPeak = peak_setA,
targetPeak = c(peak_setB, peak_setC),
TxDb = txdb,
pAdjustMethod = "BH",
nShuffle = 1000,
chainFile = NULL,
verbose = FALSE)
依次將query的peak與target中的每一個peak文件進(jìn)行overlap分析,計(jì)算出一個p值代表兩個peak之間overlap的程度,p值越小,overlap的程度越高。
ChIPseeker除了peak基因注釋的基本功能外,整合了GEO的下載功能與peak的overlap分析,可以方便的將自己的chip_seq數(shù)據(jù)與GEO的公共數(shù)據(jù)集進(jìn)行比較分析。
感謝各位的閱讀,以上就是“如何使用ChIPseeker進(jìn)行peak注釋”的內(nèi)容了,經(jīng)過本文的學(xué)習(xí)后,相信大家對如何使用ChIPseeker進(jìn)行peak注釋這一問題有了更深刻的體會,具體使用情況還需要大家實(shí)踐驗(yàn)證。這里是創(chuàng)新互聯(lián),小編將為大家推送更多相關(guān)知識點(diǎn)的文章,歡迎關(guān)注!
網(wǎng)站題目:如何使用ChIPseeker進(jìn)行peak注釋
標(biāo)題URL:http://chinadenli.net/article20/iegoco.html
成都網(wǎng)站建設(shè)公司_創(chuàng)新互聯(lián),為您提供動態(tài)網(wǎng)站、靜態(tài)網(wǎng)站、企業(yè)建站、品牌網(wǎng)站建設(shè)、網(wǎng)站導(dǎo)航、網(wǎng)站內(nèi)鏈
聲明:本網(wǎng)站發(fā)布的內(nèi)容(圖片、視頻和文字)以用戶投稿、用戶轉(zhuǎn)載內(nèi)容為主,如果涉及侵權(quán)請盡快告知,我們將會在第一時間刪除。文章觀點(diǎn)不代表本網(wǎng)站立場,如需處理請聯(lián)系客服。電話:028-86922220;郵箱:631063699@qq.com。內(nèi)容未經(jīng)允許不得轉(zhuǎn)載,或轉(zhuǎn)載時需注明來源: 創(chuàng)新互聯(lián)