snipt

Ctrl+h for KB shortcuts

Text only

Annotate BED file Chipseeker

library(ChIPseeker)
library(TxDb.Hsapiens.UCSC.hg38.knownGene)
txdb <- TxDb.Hsapiens.UCSC.hg38.knownGene
library(clusterProfiler)
library(org.Hs.eg.db)


peak <- readPeakFile("Mec1_ATM_shRNA_c_deseq2_sig_up_only.bed")

peakAnno <- annotatePeak(peak, tssRegion=c(-2000, 2000),
                         TxDb=txdb, annoDb="org.Hs.eg.db")

promoter <- getPromoters(TxDb=txdb, upstream=3000, downstream=3000)
tagMatrix <- getTagMatrix(peak, windows=promoter)
tagHeatmap(tagMatrix, xlim=c(-3000, 3000), color="red")

plotAvgProf(tagMatrix, xlim=c(-3000, 3000),
            xlab="Genomic Region (5'->3')", ylab = "Read Count Frequency")

plotAnnoPie(peakAnno)

plotAnnoBar(peakAnno)

vennpie(peakAnno)

upsetplot(peakAnno)

upsetplot(peakAnno, vennpie=TRUE)

plotDistToTSS(peakAnno,
              title="Distribution of transcription factor-binding loci\nrelative to TSS")



genes = (as.data.frame(peakAnno)$geneId)
#genes = na.omit(genes)
#names(genes) = sub("_", "\n", names(genes))


kk <- enrichKEGG(genes, organism = "hsa", pvalueCutoff=0.05, pAdjustMethod="none", 
                 qvalueCutoff=0.2)

head(summary(kk))
dotplot(kk, showCategory=25)


ego <- clusterProfiler::enrichGO(gene          = genes,
                                 OrgDb         = org.Hs.eg.db,
                                 ont           = "BP",
                                 pAdjustMethod = "BH",
                                 pvalueCutoff  = 0.05,
                                 qvalueCutoff  = 0.2, 
                                 readable      = TRUE)
head(summary(ego)[,-8])
dotplot(ego, showCategory=20)
clusterProfiler::dotplot(ego, showCategory=15)
cnetplot(ego, categorySize="geneNum", foldChange=geneList)
https://snipt.net/embed/0f28ed72a7c2b9c38c96b428780f9ecb/
/raw/0f28ed72a7c2b9c38c96b428780f9ecb/
0f28ed72a7c2b9c38c96b428780f9ecb
text
Text only
58
2019-08-21T05:01:52
True
False
False
May 14, 2019 at 11:06 AM
/api/public/snipt/148346/
annotate-bed-file-chipseeker
<table class="highlighttable"><tr><td class="linenos"><div class="linenodiv"><pre><a href="#L-1"> 1</a> <a href="#L-2"> 2</a> <a href="#L-3"> 3</a> <a href="#L-4"> 4</a> <a href="#L-5"> 5</a> <a href="#L-6"> 6</a> <a href="#L-7"> 7</a> <a href="#L-8"> 8</a> <a href="#L-9"> 9</a> <a href="#L-10">10</a> <a href="#L-11">11</a> <a href="#L-12">12</a> <a href="#L-13">13</a> <a href="#L-14">14</a> <a href="#L-15">15</a> <a href="#L-16">16</a> <a href="#L-17">17</a> <a href="#L-18">18</a> <a href="#L-19">19</a> <a href="#L-20">20</a> <a href="#L-21">21</a> <a href="#L-22">22</a> <a href="#L-23">23</a> <a href="#L-24">24</a> <a href="#L-25">25</a> <a href="#L-26">26</a> <a href="#L-27">27</a> <a href="#L-28">28</a> <a href="#L-29">29</a> <a href="#L-30">30</a> <a href="#L-31">31</a> <a href="#L-32">32</a> <a href="#L-33">33</a> <a href="#L-34">34</a> <a href="#L-35">35</a> <a href="#L-36">36</a> <a href="#L-37">37</a> <a href="#L-38">38</a> <a href="#L-39">39</a> <a href="#L-40">40</a> <a href="#L-41">41</a> <a href="#L-42">42</a> <a href="#L-43">43</a> <a href="#L-44">44</a> <a href="#L-45">45</a> <a href="#L-46">46</a> <a href="#L-47">47</a> <a href="#L-48">48</a> <a href="#L-49">49</a> <a href="#L-50">50</a> <a href="#L-51">51</a> <a href="#L-52">52</a> <a href="#L-53">53</a> <a href="#L-54">54</a> <a href="#L-55">55</a> <a href="#L-56">56</a> <a href="#L-57">57</a></pre></div></td><td class="code"><div class="highlight"><pre><span></span><span id="L-1"><a name="L-1"></a>library(ChIPseeker) </span><span id="L-2"><a name="L-2"></a>library(TxDb.Hsapiens.UCSC.hg38.knownGene) </span><span id="L-3"><a name="L-3"></a>txdb &lt;- TxDb.Hsapiens.UCSC.hg38.knownGene </span><span id="L-4"><a name="L-4"></a>library(clusterProfiler) </span><span id="L-5"><a name="L-5"></a>library(org.Hs.eg.db) </span><span id="L-6"><a name="L-6"></a> </span><span id="L-7"><a name="L-7"></a> </span><span id="L-8"><a name="L-8"></a>peak &lt;- readPeakFile(&quot;Mec1_ATM_shRNA_c_deseq2_sig_up_only.bed&quot;) </span><span id="L-9"><a name="L-9"></a> </span><span id="L-10"><a name="L-10"></a>peakAnno &lt;- annotatePeak(peak, tssRegion=c(-2000, 2000), </span><span id="L-11"><a name="L-11"></a> TxDb=txdb, annoDb=&quot;org.Hs.eg.db&quot;) </span><span id="L-12"><a name="L-12"></a> </span><span id="L-13"><a name="L-13"></a>promoter &lt;- getPromoters(TxDb=txdb, upstream=3000, downstream=3000) </span><span id="L-14"><a name="L-14"></a>tagMatrix &lt;- getTagMatrix(peak, windows=promoter) </span><span id="L-15"><a name="L-15"></a>tagHeatmap(tagMatrix, xlim=c(-3000, 3000), color=&quot;red&quot;) </span><span id="L-16"><a name="L-16"></a> </span><span id="L-17"><a name="L-17"></a>plotAvgProf(tagMatrix, xlim=c(-3000, 3000), </span><span id="L-18"><a name="L-18"></a> xlab=&quot;Genomic Region (5&#39;-&gt;3&#39;)&quot;, ylab = &quot;Read Count Frequency&quot;) </span><span id="L-19"><a name="L-19"></a> </span><span id="L-20"><a name="L-20"></a>plotAnnoPie(peakAnno) </span><span id="L-21"><a name="L-21"></a> </span><span id="L-22"><a name="L-22"></a>plotAnnoBar(peakAnno) </span><span id="L-23"><a name="L-23"></a> </span><span id="L-24"><a name="L-24"></a>vennpie(peakAnno) </span><span id="L-25"><a name="L-25"></a> </span><span id="L-26"><a name="L-26"></a>upsetplot(peakAnno) </span><span id="L-27"><a name="L-27"></a> </span><span id="L-28"><a name="L-28"></a>upsetplot(peakAnno, vennpie=TRUE) </span><span id="L-29"><a name="L-29"></a> </span><span id="L-30"><a name="L-30"></a>plotDistToTSS(peakAnno, </span><span id="L-31"><a name="L-31"></a> title=&quot;Distribution of transcription factor-binding loci\nrelative to TSS&quot;) </span><span id="L-32"><a name="L-32"></a> </span><span id="L-33"><a name="L-33"></a> </span><span id="L-34"><a name="L-34"></a> </span><span id="L-35"><a name="L-35"></a>genes = (as.data.frame(peakAnno)$geneId) </span><span id="L-36"><a name="L-36"></a>#genes = na.omit(genes) </span><span id="L-37"><a name="L-37"></a>#names(genes) = sub(&quot;_&quot;, &quot;\n&quot;, names(genes)) </span><span id="L-38"><a name="L-38"></a> </span><span id="L-39"><a name="L-39"></a> </span><span id="L-40"><a name="L-40"></a>kk &lt;- enrichKEGG(genes, organism = &quot;hsa&quot;, pvalueCutoff=0.05, pAdjustMethod=&quot;none&quot;, </span><span id="L-41"><a name="L-41"></a> qvalueCutoff=0.2) </span><span id="L-42"><a name="L-42"></a> </span><span id="L-43"><a name="L-43"></a>head(summary(kk)) </span><span id="L-44"><a name="L-44"></a>dotplot(kk, showCategory=25) </span><span id="L-45"><a name="L-45"></a> </span><span id="L-46"><a name="L-46"></a> </span><span id="L-47"><a name="L-47"></a>ego &lt;- clusterProfiler::enrichGO(gene = genes, </span><span id="L-48"><a name="L-48"></a> OrgDb = org.Hs.eg.db, </span><span id="L-49"><a name="L-49"></a> ont = &quot;BP&quot;, </span><span id="L-50"><a name="L-50"></a> pAdjustMethod = &quot;BH&quot;, </span><span id="L-51"><a name="L-51"></a> pvalueCutoff = 0.05, </span><span id="L-52"><a name="L-52"></a> qvalueCutoff = 0.2, </span><span id="L-53"><a name="L-53"></a> readable = TRUE) </span><span id="L-54"><a name="L-54"></a>head(summary(ego)[,-8]) </span><span id="L-55"><a name="L-55"></a>dotplot(ego, showCategory=20) </span><span id="L-56"><a name="L-56"></a>clusterProfiler::dotplot(ego, showCategory=15) </span><span id="L-57"><a name="L-57"></a>cnetplot(ego, categorySize=&quot;geneNum&quot;, foldChange=geneList) </span></pre></div> </td></tr></table>
--- 
+++ 
@@ -0,0 +1,57 @@
+library(ChIPseeker)
+library(TxDb.Hsapiens.UCSC.hg38.knownGene)
+txdb <- TxDb.Hsapiens.UCSC.hg38.knownGene
+library(clusterProfiler)
+library(org.Hs.eg.db)
+
+
+peak <- readPeakFile("Mec1_ATM_shRNA_c_deseq2_sig_up_only.bed")
+
+peakAnno <- annotatePeak(peak, tssRegion=c(-2000, 2000),
+                         TxDb=txdb, annoDb="org.Hs.eg.db")
+
+promoter <- getPromoters(TxDb=txdb, upstream=3000, downstream=3000)
+tagMatrix <- getTagMatrix(peak, windows=promoter)
+tagHeatmap(tagMatrix, xlim=c(-3000, 3000), color="red")
+
+plotAvgProf(tagMatrix, xlim=c(-3000, 3000),
+            xlab="Genomic Region (5'->3')", ylab = "Read Count Frequency")
+
+plotAnnoPie(peakAnno)
+
+plotAnnoBar(peakAnno)
+
+vennpie(peakAnno)
+
+upsetplot(peakAnno)
+
+upsetplot(peakAnno, vennpie=TRUE)
+
+plotDistToTSS(peakAnno,
+              title="Distribution of transcription factor-binding loci\nrelative to TSS")
+
+
+
+genes = (as.data.frame(peakAnno)$geneId)
+#genes = na.omit(genes)
+#names(genes) = sub("_", "\n", names(genes))
+
+
+kk <- enrichKEGG(genes, organism = "hsa", pvalueCutoff=0.05, pAdjustMethod="none", 
+                 qvalueCutoff=0.2)
+
+head(summary(kk))
+dotplot(kk, showCategory=25)
+
+
+ego <- clusterProfiler::enrichGO(gene          = genes,
+                                 OrgDb         = org.Hs.eg.db,
+                                 ont           = "BP",
+                                 pAdjustMethod = "BH",
+                                 pvalueCutoff  = 0.05,
+                                 qvalueCutoff  = 0.2, 
+                                 readable      = TRUE)
+head(summary(ego)[,-8])
+dotplot(ego, showCategory=20)
+clusterProfiler::dotplot(ego, showCategory=15)
+cnetplot(ego, categorySize="geneNum", foldChange=geneList)
  • SaelM
  • 3 months, 1 week ago