-
Notifications
You must be signed in to change notification settings - Fork 2
/
Copy pathvisualization.Rmd
173 lines (153 loc) · 7.7 KB
/
visualization.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
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
---
title: "multi-omics_visualization"
author:
- "ddedesener"
- "DeniseSl22"
date: "21/07/22"
output:
md_document:
variant: markdown_github
always_allow_html: true
---
## Introduction
In this script, visualization of the enriched pathway which include both altered genes and metabolites is performed.
## Setup
```{r setup, warning=FALSE, message=FALSE}
# check if libraries are already installed > otherwise install it
if (!requireNamespace("BiocManager", quietly = TRUE)) install.packages("BiocManager")
if(!"RCy3" %in% installed.packages()) BiocManager::install("RCy3")
if(!"rWikiPathways" %in% installed.packages()) BiocManager::install("rWikiPathways")
if(!"RColorBrewer" %in% installed.packages()) BiocManager::install("RColorBrewer")
if(!"dplyr" %in% installed.packages()) BiocManager::install("dplyr")
#load libraries
library(RCy3)#connect cytoscape via R
library(rWikiPathways)#to get pathways from WikiPathways
library(RColorBrewer)#to manage colors with R
library(dplyr)
```
##Obtain transcriptomics and metabolomics data
```{r}
setwd('../..')
work_DIR <- getwd()
#Set location to download data for transcriptomics:
filelocation_t <- paste0(work_DIR, "/transcriptomics_analysis/3-identifier_mapping/output/")
#Obtain data from step 3
tSet_CD <- read.delim(paste0(filelocation_t, 'IDMapping_CD.tsv'), sep = "\t", na.strings=c("", "NA"))
tSet_UC <- read.delim(paste0(filelocation_t, 'IDMapping_UC.tsv'), sep = "\t", na.strings=c("", "NA"))
##Select the corresponding location to be visualized:
location_transcriptomics <- "ileum" ##Options: ileum or rectum
#filter out unused columns
if(location_transcriptomics == "ileum"){
tSet_CD <- na.omit(tSet_CD [,c(1,4:5)])
tSet_UC <- na.omit(tSet_UC [,c(1,4:5)])}else if(location_transcriptomics == "rectum"){
tSet_CD <- na.omit(tSet_CD [,c(1,6:7)])
tSet_UC <- na.omit(tSet_UC [,c(1,6:7)])}else{print("Location for transcriptomics data not recognised.")}
#Rename columns for merger later
colnames(tSet_CD) <- c('ID','log2FC','pvalues')
colnames(tSet_UC) <- c('ID','log2FC','pvalues')
#Set location to download data for metabolomics:
filelocation_m <- paste0(work_DIR, "/metabolomics_analysis/10-identifier_mapping/output/")
#Obtain data from step 10
mSet_CD <- read.csv(paste0(filelocation_m, 'mbx_mapped_data_CD.tsv'), sep = "\t", na.strings=c("", "NA"))
mSet_UC <- read.csv(paste0(filelocation_m, 'mbx_mapped_data_UC.tsv'), sep = "\t", na.strings=c("", "NA"))
#filter out unused columns
mSet_CD <- na.omit(mSet_CD [,c(2,4:5)])
mSet_UC <- na.omit(mSet_UC [,c(2,4:5)])
#Rename columns for merger later
colnames(mSet_CD) <- c('ID','log2FC','pvalues')
colnames(mSet_UC) <- c('ID','log2FC','pvalues')
```
##Combine both dataframes
```{r}
combined.data_CD <- rbind(tSet_CD, mSet_CD)
combined.data_UC <- rbind(tSet_UC, mSet_UC)
##Select disorder to visualize later on:
disorder <- "CD" ##Options are: CD,UC
if(disorder == "CD"){combined.data <- combined.data_CD}else if(disorder == "UC"){combined.data <- combined.data_UC}else{print("Disorder not recognized.")}
```
## Import pathway
```{r import, warning=FALSE, message=FALSE}
#make sure to launch cytoscape, if you get CyREST error you need to relaunch cytoscape
cytoscapePing()
#close all opened session before starting
closeSession(FALSE)
#Set up WikiPathways app in Cytoscape, v.3.3.10
if("WikiPathways" %in% commandsHelp("")) print("Success: the WikiPathways app is installed") else print("Warning: WikiPathways app is not installed. Please install the WikiPathways app before proceeding.")
if(!"WikiPathways" %in% commandsHelp("")) installApp("WikiPathways")
#pathway IDs to be visualized
pathway.id <- "WP4726"# Sphingolipid metabolism: integrated pathway is a relevant and significantly altered pathways regarding metabolomics data.
#import pathways as pathway in cytoscape
RCy3::commandsRun(paste0('wikipathways import-as-pathway id=',pathway.id ))
```
## Data upload
```{r data_upload, warning=FALSE, message=FALSE}
#get node table from imported pathway in cytoscape
ID.cols <- getTableColumns(table ="node", columns = c("XrefId","Ensembl", "ChEBI"))
#filter out rows which contain NA value for columns Ensembl and ChEBI
ID.cols <- ID.cols[!with(ID.cols, is.na(Ensembl) & is.na(ChEBI)),]
#if a row value in the Ensembl column is NA then replace it with ChEBI
ID.cols$Ensembl <- ifelse(is.na(ID.cols$Ensembl), ID.cols$ChEBI, ID.cols$Ensembl)
#use the only one column contains both Ensembl and ChEBI identifiers
ID.cols <- data.frame(ID.cols[,c(1,2)])
#change column name
colnames(ID.cols)[2] <- "omics.ID"
##Remove ':T" in mapped IDs for now:
ID.cols$omics.ID<-gsub(":T","", as.character(ID.cols$omics.ID))
#merge two data frames for adding xrefid to the combined data frame
data <- merge(combined.data, ID.cols, by.x = "ID", by.y = "omics.ID" )
#remove duplicate rows
data <- data %>% distinct(ID, .keep_all = TRUE)
colnames(data)[1] <- "omics.ID"
#load data to the imported pathway in cytoscape by key column as XrefId
loadTableData(table = "node", data = data, data.key.column = "XrefId", table.key.column = "XrefId")
```
## Visualization options
```{r visualization, warning=FALSE, message=FALSE}
#new visual style is created
RCy3::copyVisualStyle("default","pathwayStyle")
#set new style as the current style
RCy3::setVisualStyle("pathwayStyle")
#set node dimensions as fixed sizes
RCy3::lockNodeDimensions(TRUE, style.name="pathwayStyle")
#node shape mapping
RCy3::setNodeShapeMapping('Type',c('GeneProduct','Protein', 'Metabolite'),c('ELLIPSE','ELLIPSE','RECTANGLE'), style.name="pathwayStyle")
#change node height
RCy3::setNodeHeightMapping('Type',c('GeneProduct','Protein', 'Metabolite'), c(23,23,25), mapping.type = "d", style.name = "pathwayStyle")
#change node width
RCy3::setNodeWidthMapping('Type',c('GeneProduct','Protein', 'Metabolite'), c(60,60,100), mapping.type = "d", style.name = "pathwayStyle")
#set node color based on log2FC for both genes and metabolites
node.colors <- c(rev(brewer.pal(3, "RdBu")))
setNodeColorMapping("log2FC", c(-1,0,1), node.colors, default.color = "#D3D3D3", style.name = "pathwayStyle")
#Set node border width and color based on p-value
#First we need to get all p-values from node table
pvalues <- getTableColumns(table = 'node', columns = 'pvalues')
pvalues <- na.omit(pvalues)
#Create a range for all sign. p-values, and one for all not significant.
significant_pvalues <- pvalues[(pvalues < 0.05)]
not.significant_pvalues <- pvalues[(pvalues >= 0.05)]
significant_pvalues.colors <- rep("#2e9d1d", length(significant_pvalues))
not.significant_pvalues.colors <- rep("#FFFFFF", length(not.significant_pvalues))
RCy3::setNodeBorderWidthMapping('pvalues', table.column.values = NULL , c(6,6) , mapping.type = "c", style.name = "pathwayStyle")
RCy3::setNodeBorderColorMapping('pvalues', c(significant_pvalues,not.significant_pvalues), c(significant_pvalues.colors, not.significant_pvalues.colors), default.color = "#AAAAAA", mapping.type = "d", style.name = "pathwayStyle")
##Update relevant interactions to directional ones:
RCy3::setEdgeTargetArrowShapeMapping(table.column = 'EndArrow', c('mim-conversion', 'Arrow', 'mim-catalysis'), c('DELTA', 'DELTA', 'OPEN_CIRCLE'), style.name = "pathwayStyle")
work_DIR
#Save output
filename_multiomics <- paste0(work_DIR, "/visualization_multiomics/12-visualization/output/", pathway.id, "_", disorder, "_location_", location_transcriptomics,"_visualization.png")
png.file <- file.path( filename_multiomics)
exportImage(png.file, 'PNG', zoom = 500)
```
##Print session info:
```{r print_session_info}
##Print session info:
sessionInfo()
```
## Creating jupyter files
```{r writing_to_notebooks,warning=FALSE, message=FALSE }
#Jupyter Notebook file
if(!"devtools" %in% installed.packages()) BiocManager::install("devtools")
devtools::install_github("mkearney/rmd2jupyter", force=TRUE)
library(devtools)
library(rmd2jupyter)
rmd2jupyter("visualization.Rmd")
```