forked from YuLab-SMU/biomedical-knowledge-mining-book
-
Notifications
You must be signed in to change notification settings - Fork 0
/
09_Reactome.Rmd
119 lines (74 loc) · 4.3 KB
/
09_Reactome.Rmd
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
# Reactome enrichment analysis {#reactomepa}
```{r include=FALSE}
library(knitr)
opts_chunk$set(message=FALSE, warning=FALSE, eval=TRUE, echo=TRUE, cache=TRUE)
library(clusterProfiler)
library(ReactomePA)
```
`r Biocpkg("ReactomePA")` is designed for reactome pathway based analysis [@yu_reactomepa_2016]. Reactome is an open-source, open access, manually curated and peer-reviewed pathway database.
## Supported organisms {#reactomepa-supported-organisms}
Currently `r Biocpkg("ReactomePA")` supports several model organisms, including 'celegans', 'fly', 'human', 'mouse', 'rat', 'yeast' and 'zebrafish'. The input gene ID should be Entrez gene ID. We recommend using [`clusterProfiler::bitr()`](#bitr) to convert biological IDs.
## Reactome pathway over-representation analysis {#reactomepa-ora}
Enrichment analysis is a widely used approach to identify biological
themes. `r Biocpkg("ReactomePA")` implemented `enrichPathway()` that uses [hypergeometric model](#ora-algorithm) to assess whether the number of selected genes associated with a reactome pathway is larger than expected.
```{r enrichpathway}
library(ReactomePA)
data(geneList, package="DOSE"
)
de <- names(geneList)[abs(geneList) > 1.5]
head(de)
x <- enrichPathway(gene=de, pvalueCutoff = 0.05, readable=TRUE)
head(x)
```
## Reactome pathway gene set enrichment analysis {#reactomepa-gsea}
```{r gsepathway}
y <- gsePathway(geneList,
pvalueCutoff = 0.2,
pAdjustMethod = "BH",
verbose = FALSE)
head(y)
```
## Pathway Visualization
`r Biocpkg("ReactomePA")` implemented the `viewPathway()` to visualized selected reactome pathway. More general purpose of visualization methods for ORA and GSEA results are provided in the `r Biocpkg("enrichplot")` package and are documented on [Chapter 14](#enrichplot).
(ref:viewpathwayscap) Visualize reactome pathway.
(ref:viewpathwaycap) **Visualize reactome pathway.**
```{r viewpathway, fig.height=6, fig.width=8, fig.cap="(ref:viewpathwaycap)", fig.scap="(ref:viewpathwayscap)"}
viewPathway("E2F mediated regulation of DNA replication",
readable = TRUE,
foldChange = geneList)
```
<!--
## Pathway analysis of NGS data
Pathway analysis using NGS data (eg, RNA-Seq and ChIP-Seq) can be performed by linking coding and non-coding regions to coding genes via `r Biocpkg("ChIPseeker")` package, which can annotates genomic regions to their nearest genes, host genes, and flanking genes respectivly. In addtion, it provides a function, __*seq2gene*__, that simultaneously considering host genes, promoter region and flanking gene from intergenic region that may under control via cis-regulation. This function maps genomic regions to genes in a many-to-many manner and facilitate functional analysis. For more details, please refer to `r Biocpkg("ChIPseeker")`[@yu_chipseeker_2015].
## Visualize enrichment result
We implement barplot, dotplot enrichment map and category-gene-network for visualization. It is very common to visualize the enrichment result in bar or pie chart. We believe the pie chart is misleading and only provide bar chart.
```r fig.height=6, fig.width=12}
barplot(x, showCategory=8)
```
```r fig.height=6, fig.width=12}
dotplot(x, showCategory=15)
```
Enrichment map can be viusalized by __*enrichMap*__:
```r fig.height=10, fig.width=10}
emapplot(x)
```
In order to consider the potentially biological complexities in which a gene may belong to multiple annotation categories, we developed __*cnetplot*__ function to extract the complex association between genes and diseases.
```r fig.height=8, fig.width=8}
cnetplot(x, categorySize="pvalue", foldChange=geneList)
```
## Visualize GSEA result
```r fig.height=8, fig.width=8}
emapplot(y, color="pvalue")
```
```r fig.height=7, fig.width=10}
gseaplot(y, geneSetID = "R-HSA-69242")
```
## Comparing enriched reactome pathways among gene clusters with clusterProfiler
We have developed an `R` package `r Biocpkg("clusterProfiler")`[@yu_clusterprofiler:_2012] for comparing biological themes among gene clusters. `r Biocpkg("ReactomePA")` works fine with `r Biocpkg("clusterProfiler")` and can compare biological themes at reactome pathway perspective.
```r fig.height=8, fig.width=13, eval=FALSE}
require(clusterProfiler)
data(gcSample)
res <- compareCluster(gcSample, fun="enrichPathway")
dotplot(res)
```
-->