############################################
library(ggplot2)
library(dplyr)
##
## Attaching package: 'dplyr'
## The following objects are masked from 'package:stats':
##
## filter, lag
## The following objects are masked from 'package:base':
##
## intersect, setdiff, setequal, union
library(aplot)
# Create simulated data
set.seed(123)
n_genes <- 1000
gene_values <- sort(rnorm(n_genes), decreasing = TRUE)
gene_names <- paste0("gene", 1:n_genes)
# Create geneList as data.frame
geneList_df <- data.frame(
name = gene_names,
value = gene_values,
stringsAsFactors = FALSE
)
head(geneList_df); dim(geneList_df)
## name value
## 1 gene1 3.241040
## 2 gene2 2.691714
## 3 gene3 2.684859
## 4 gene4 2.575450
## 5 gene5 2.571458
## 6 gene6 2.553026
## [1] 1000 2
# Create example pathway_genes with runningScore
pathway_size <- 15
selected_genes <- sample(gene_names, pathway_size)
# Simulate running score calculation (replace this with your actual runningScore)
runningScore <- cumsum(ifelse(gene_names %in% selected_genes,
1/pathway_size,
-1/(n_genes-pathway_size)))
runningScore <- runningScore/max(abs(runningScore))
# Create pathway_genes data.frame
pathway_genes <- data.frame(
name = gene_names, # All gene names
runningScore = runningScore, # Corresponding running scores
stringsAsFactors = FALSE
)
head(pathway_genes); dim(pathway_genes)
## name runningScore
## 1 gene1 -0.004484305
## 2 gene2 -0.008968610
## 3 gene3 -0.013452915
## 4 gene4 -0.017937220
## 5 gene5 -0.022421525
## 6 gene6 -0.026905830
## [1] 1000 2
unique(pathway_genes$name)
## [1] "gene1" "gene2" "gene3" "gene4" "gene5" "gene6"
## [7] "gene7" "gene8" "gene9" "gene10" "gene11" "gene12"
## [13] "gene13" "gene14" "gene15" "gene16" "gene17" "gene18"
## [19] "gene19" "gene20" "gene21" "gene22" "gene23" "gene24"
## [25] "gene25" "gene26" "gene27" "gene28" "gene29" "gene30"
## [31] "gene31" "gene32" "gene33" "gene34" "gene35" "gene36"
## [37] "gene37" "gene38" "gene39" "gene40" "gene41" "gene42"
## [43] "gene43" "gene44" "gene45" "gene46" "gene47" "gene48"
## [49] "gene49" "gene50" "gene51" "gene52" "gene53" "gene54"
## [55] "gene55" "gene56" "gene57" "gene58" "gene59" "gene60"
## [61] "gene61" "gene62" "gene63" "gene64" "gene65" "gene66"
## [67] "gene67" "gene68" "gene69" "gene70" "gene71" "gene72"
## [73] "gene73" "gene74" "gene75" "gene76" "gene77" "gene78"
## [79] "gene79" "gene80" "gene81" "gene82" "gene83" "gene84"
## [85] "gene85" "gene86" "gene87" "gene88" "gene89" "gene90"
## [91] "gene91" "gene92" "gene93" "gene94" "gene95" "gene96"
## [97] "gene97" "gene98" "gene99" "gene100" "gene101" "gene102"
## [103] "gene103" "gene104" "gene105" "gene106" "gene107" "gene108"
## [109] "gene109" "gene110" "gene111" "gene112" "gene113" "gene114"
## [115] "gene115" "gene116" "gene117" "gene118" "gene119" "gene120"
## [121] "gene121" "gene122" "gene123" "gene124" "gene125" "gene126"
## [127] "gene127" "gene128" "gene129" "gene130" "gene131" "gene132"
## [133] "gene133" "gene134" "gene135" "gene136" "gene137" "gene138"
## [139] "gene139" "gene140" "gene141" "gene142" "gene143" "gene144"
## [145] "gene145" "gene146" "gene147" "gene148" "gene149" "gene150"
## [151] "gene151" "gene152" "gene153" "gene154" "gene155" "gene156"
## [157] "gene157" "gene158" "gene159" "gene160" "gene161" "gene162"
## [163] "gene163" "gene164" "gene165" "gene166" "gene167" "gene168"
## [169] "gene169" "gene170" "gene171" "gene172" "gene173" "gene174"
## [175] "gene175" "gene176" "gene177" "gene178" "gene179" "gene180"
## [181] "gene181" "gene182" "gene183" "gene184" "gene185" "gene186"
## [187] "gene187" "gene188" "gene189" "gene190" "gene191" "gene192"
## [193] "gene193" "gene194" "gene195" "gene196" "gene197" "gene198"
## [199] "gene199" "gene200" "gene201" "gene202" "gene203" "gene204"
## [205] "gene205" "gene206" "gene207" "gene208" "gene209" "gene210"
## [211] "gene211" "gene212" "gene213" "gene214" "gene215" "gene216"
## [217] "gene217" "gene218" "gene219" "gene220" "gene221" "gene222"
## [223] "gene223" "gene224" "gene225" "gene226" "gene227" "gene228"
## [229] "gene229" "gene230" "gene231" "gene232" "gene233" "gene234"
## [235] "gene235" "gene236" "gene237" "gene238" "gene239" "gene240"
## [241] "gene241" "gene242" "gene243" "gene244" "gene245" "gene246"
## [247] "gene247" "gene248" "gene249" "gene250" "gene251" "gene252"
## [253] "gene253" "gene254" "gene255" "gene256" "gene257" "gene258"
## [259] "gene259" "gene260" "gene261" "gene262" "gene263" "gene264"
## [265] "gene265" "gene266" "gene267" "gene268" "gene269" "gene270"
## [271] "gene271" "gene272" "gene273" "gene274" "gene275" "gene276"
## [277] "gene277" "gene278" "gene279" "gene280" "gene281" "gene282"
## [283] "gene283" "gene284" "gene285" "gene286" "gene287" "gene288"
## [289] "gene289" "gene290" "gene291" "gene292" "gene293" "gene294"
## [295] "gene295" "gene296" "gene297" "gene298" "gene299" "gene300"
## [301] "gene301" "gene302" "gene303" "gene304" "gene305" "gene306"
## [307] "gene307" "gene308" "gene309" "gene310" "gene311" "gene312"
## [313] "gene313" "gene314" "gene315" "gene316" "gene317" "gene318"
## [319] "gene319" "gene320" "gene321" "gene322" "gene323" "gene324"
## [325] "gene325" "gene326" "gene327" "gene328" "gene329" "gene330"
## [331] "gene331" "gene332" "gene333" "gene334" "gene335" "gene336"
## [337] "gene337" "gene338" "gene339" "gene340" "gene341" "gene342"
## [343] "gene343" "gene344" "gene345" "gene346" "gene347" "gene348"
## [349] "gene349" "gene350" "gene351" "gene352" "gene353" "gene354"
## [355] "gene355" "gene356" "gene357" "gene358" "gene359" "gene360"
## [361] "gene361" "gene362" "gene363" "gene364" "gene365" "gene366"
## [367] "gene367" "gene368" "gene369" "gene370" "gene371" "gene372"
## [373] "gene373" "gene374" "gene375" "gene376" "gene377" "gene378"
## [379] "gene379" "gene380" "gene381" "gene382" "gene383" "gene384"
## [385] "gene385" "gene386" "gene387" "gene388" "gene389" "gene390"
## [391] "gene391" "gene392" "gene393" "gene394" "gene395" "gene396"
## [397] "gene397" "gene398" "gene399" "gene400" "gene401" "gene402"
## [403] "gene403" "gene404" "gene405" "gene406" "gene407" "gene408"
## [409] "gene409" "gene410" "gene411" "gene412" "gene413" "gene414"
## [415] "gene415" "gene416" "gene417" "gene418" "gene419" "gene420"
## [421] "gene421" "gene422" "gene423" "gene424" "gene425" "gene426"
## [427] "gene427" "gene428" "gene429" "gene430" "gene431" "gene432"
## [433] "gene433" "gene434" "gene435" "gene436" "gene437" "gene438"
## [439] "gene439" "gene440" "gene441" "gene442" "gene443" "gene444"
## [445] "gene445" "gene446" "gene447" "gene448" "gene449" "gene450"
## [451] "gene451" "gene452" "gene453" "gene454" "gene455" "gene456"
## [457] "gene457" "gene458" "gene459" "gene460" "gene461" "gene462"
## [463] "gene463" "gene464" "gene465" "gene466" "gene467" "gene468"
## [469] "gene469" "gene470" "gene471" "gene472" "gene473" "gene474"
## [475] "gene475" "gene476" "gene477" "gene478" "gene479" "gene480"
## [481] "gene481" "gene482" "gene483" "gene484" "gene485" "gene486"
## [487] "gene487" "gene488" "gene489" "gene490" "gene491" "gene492"
## [493] "gene493" "gene494" "gene495" "gene496" "gene497" "gene498"
## [499] "gene499" "gene500" "gene501" "gene502" "gene503" "gene504"
## [505] "gene505" "gene506" "gene507" "gene508" "gene509" "gene510"
## [511] "gene511" "gene512" "gene513" "gene514" "gene515" "gene516"
## [517] "gene517" "gene518" "gene519" "gene520" "gene521" "gene522"
## [523] "gene523" "gene524" "gene525" "gene526" "gene527" "gene528"
## [529] "gene529" "gene530" "gene531" "gene532" "gene533" "gene534"
## [535] "gene535" "gene536" "gene537" "gene538" "gene539" "gene540"
## [541] "gene541" "gene542" "gene543" "gene544" "gene545" "gene546"
## [547] "gene547" "gene548" "gene549" "gene550" "gene551" "gene552"
## [553] "gene553" "gene554" "gene555" "gene556" "gene557" "gene558"
## [559] "gene559" "gene560" "gene561" "gene562" "gene563" "gene564"
## [565] "gene565" "gene566" "gene567" "gene568" "gene569" "gene570"
## [571] "gene571" "gene572" "gene573" "gene574" "gene575" "gene576"
## [577] "gene577" "gene578" "gene579" "gene580" "gene581" "gene582"
## [583] "gene583" "gene584" "gene585" "gene586" "gene587" "gene588"
## [589] "gene589" "gene590" "gene591" "gene592" "gene593" "gene594"
## [595] "gene595" "gene596" "gene597" "gene598" "gene599" "gene600"
## [601] "gene601" "gene602" "gene603" "gene604" "gene605" "gene606"
## [607] "gene607" "gene608" "gene609" "gene610" "gene611" "gene612"
## [613] "gene613" "gene614" "gene615" "gene616" "gene617" "gene618"
## [619] "gene619" "gene620" "gene621" "gene622" "gene623" "gene624"
## [625] "gene625" "gene626" "gene627" "gene628" "gene629" "gene630"
## [631] "gene631" "gene632" "gene633" "gene634" "gene635" "gene636"
## [637] "gene637" "gene638" "gene639" "gene640" "gene641" "gene642"
## [643] "gene643" "gene644" "gene645" "gene646" "gene647" "gene648"
## [649] "gene649" "gene650" "gene651" "gene652" "gene653" "gene654"
## [655] "gene655" "gene656" "gene657" "gene658" "gene659" "gene660"
## [661] "gene661" "gene662" "gene663" "gene664" "gene665" "gene666"
## [667] "gene667" "gene668" "gene669" "gene670" "gene671" "gene672"
## [673] "gene673" "gene674" "gene675" "gene676" "gene677" "gene678"
## [679] "gene679" "gene680" "gene681" "gene682" "gene683" "gene684"
## [685] "gene685" "gene686" "gene687" "gene688" "gene689" "gene690"
## [691] "gene691" "gene692" "gene693" "gene694" "gene695" "gene696"
## [697] "gene697" "gene698" "gene699" "gene700" "gene701" "gene702"
## [703] "gene703" "gene704" "gene705" "gene706" "gene707" "gene708"
## [709] "gene709" "gene710" "gene711" "gene712" "gene713" "gene714"
## [715] "gene715" "gene716" "gene717" "gene718" "gene719" "gene720"
## [721] "gene721" "gene722" "gene723" "gene724" "gene725" "gene726"
## [727] "gene727" "gene728" "gene729" "gene730" "gene731" "gene732"
## [733] "gene733" "gene734" "gene735" "gene736" "gene737" "gene738"
## [739] "gene739" "gene740" "gene741" "gene742" "gene743" "gene744"
## [745] "gene745" "gene746" "gene747" "gene748" "gene749" "gene750"
## [751] "gene751" "gene752" "gene753" "gene754" "gene755" "gene756"
## [757] "gene757" "gene758" "gene759" "gene760" "gene761" "gene762"
## [763] "gene763" "gene764" "gene765" "gene766" "gene767" "gene768"
## [769] "gene769" "gene770" "gene771" "gene772" "gene773" "gene774"
## [775] "gene775" "gene776" "gene777" "gene778" "gene779" "gene780"
## [781] "gene781" "gene782" "gene783" "gene784" "gene785" "gene786"
## [787] "gene787" "gene788" "gene789" "gene790" "gene791" "gene792"
## [793] "gene793" "gene794" "gene795" "gene796" "gene797" "gene798"
## [799] "gene799" "gene800" "gene801" "gene802" "gene803" "gene804"
## [805] "gene805" "gene806" "gene807" "gene808" "gene809" "gene810"
## [811] "gene811" "gene812" "gene813" "gene814" "gene815" "gene816"
## [817] "gene817" "gene818" "gene819" "gene820" "gene821" "gene822"
## [823] "gene823" "gene824" "gene825" "gene826" "gene827" "gene828"
## [829] "gene829" "gene830" "gene831" "gene832" "gene833" "gene834"
## [835] "gene835" "gene836" "gene837" "gene838" "gene839" "gene840"
## [841] "gene841" "gene842" "gene843" "gene844" "gene845" "gene846"
## [847] "gene847" "gene848" "gene849" "gene850" "gene851" "gene852"
## [853] "gene853" "gene854" "gene855" "gene856" "gene857" "gene858"
## [859] "gene859" "gene860" "gene861" "gene862" "gene863" "gene864"
## [865] "gene865" "gene866" "gene867" "gene868" "gene869" "gene870"
## [871] "gene871" "gene872" "gene873" "gene874" "gene875" "gene876"
## [877] "gene877" "gene878" "gene879" "gene880" "gene881" "gene882"
## [883] "gene883" "gene884" "gene885" "gene886" "gene887" "gene888"
## [889] "gene889" "gene890" "gene891" "gene892" "gene893" "gene894"
## [895] "gene895" "gene896" "gene897" "gene898" "gene899" "gene900"
## [901] "gene901" "gene902" "gene903" "gene904" "gene905" "gene906"
## [907] "gene907" "gene908" "gene909" "gene910" "gene911" "gene912"
## [913] "gene913" "gene914" "gene915" "gene916" "gene917" "gene918"
## [919] "gene919" "gene920" "gene921" "gene922" "gene923" "gene924"
## [925] "gene925" "gene926" "gene927" "gene928" "gene929" "gene930"
## [931] "gene931" "gene932" "gene933" "gene934" "gene935" "gene936"
## [937] "gene937" "gene938" "gene939" "gene940" "gene941" "gene942"
## [943] "gene943" "gene944" "gene945" "gene946" "gene947" "gene948"
## [949] "gene949" "gene950" "gene951" "gene952" "gene953" "gene954"
## [955] "gene955" "gene956" "gene957" "gene958" "gene959" "gene960"
## [961] "gene961" "gene962" "gene963" "gene964" "gene965" "gene966"
## [967] "gene967" "gene968" "gene969" "gene970" "gene971" "gene972"
## [973] "gene973" "gene974" "gene975" "gene976" "gene977" "gene978"
## [979] "gene979" "gene980" "gene981" "gene982" "gene983" "gene984"
## [985] "gene985" "gene986" "gene987" "gene988" "gene989" "gene990"
## [991] "gene991" "gene992" "gene993" "gene994" "gene995" "gene996"
## [997] "gene997" "gene998" "gene999" "gene1000"
# Get number of genes
n_genes <- nrow(geneList_df)
# Create gsdata
gsdata <- data.frame(
x = 1:n_genes,
runningScore = pathway_genes$runningScore,
geneList = geneList_df$value,
position = ifelse(geneList_df$name %in% selected_genes, 1, 0),
id = "Pathway1"
)
dim(gsdata); head(gsdata)
## [1] 1000 5
## x runningScore geneList position id
## 1 1 -0.004484305 3.241040 0 Pathway1
## 2 2 -0.008968610 2.691714 0 Pathway1
## 3 3 -0.013452915 2.684859 0 Pathway1
## 4 4 -0.017937220 2.575450 0 Pathway1
## 5 5 -0.022421525 2.571458 0 Pathway1
## 6 6 -0.026905830 2.553026 0 Pathway1
# Create gsdata1
gsdata1 <- gsdata %>%
filter(position == 1) %>%
mutate(gene_name = geneList_df$name[which(gsdata$position == 1)])
# Calculate dynamic breaks
tmp_positions <- gsdata$position
v <- seq(1, sum(tmp_positions), length.out = 9)
inv <- findInterval(rev(cumsum(tmp_positions)), v)
if (min(inv) == 0) inv <- inv + 1
# Create result data frame for heatmap
result_df <- data.frame(
Interval = 1:length(unique(inv)),
xmin = which(!duplicated(inv)),
xmax = c(which(!duplicated(inv))[-1], length(inv)),
ymin = 0,
ymax = 0.3,
Mean_LogFC = quantile(geneList_df$value, probs = seq(0, 1, length.out = length(unique(inv)))),
id = "Pathway1"
)
# Plot parameters
curveCol <- c("#76BA99", "#EB4747")
htCol <- c("#08519C", "#A50F15")
lineSize <- 1.2
# 1. Enrichment Score Curve
pcurve <- ggplot(gsdata, aes(x = x, y = runningScore)) +
geom_line(aes(color = runningScore), size = lineSize) +
scale_color_gradient2(
low = curveCol[1],
mid = "yellow",
high = curveCol[2],
midpoint = mean(range(gsdata$runningScore))
) +
geom_hline(yintercept = 0, size = lineSize/2, color = "black", lty = "dashed") +
scale_x_continuous(expand = c(0, 0), limits = c(1, n_genes)) +
theme_bw(base_size = 12) +
theme(
legend.position = "none",
panel.grid = element_blank(),
axis.ticks.x = element_blank(),
axis.text.x = element_blank(),
axis.line.x = element_blank(),
axis.title.x = element_blank(),
plot.margin = unit(c(0.2, 0.2, 0, 0.2), "cm")
) +
ylab("Running Enrichment Score")
## Warning: Using `size` aesthetic for lines was deprecated in ggplot2 3.4.0.
## ℹ Please use `linewidth` instead.
## This warning is displayed once every 8 hours.
## Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
## generated.
# 2. Gene Position Ticks with Heatmap
pseg_ht <- ggplot(gsdata, aes(x = x, y = runningScore)) +
geom_segment(
data = gsdata1,
aes(x = x, xend = x, y = 0, yend = 1),
color = "black",
size = 0.5,
show.legend = FALSE
) +
geom_rect(
data = result_df,
aes(xmin = xmin, xmax = xmax, ymin = ymin, ymax = ymax, fill = Mean_LogFC),
inherit.aes = FALSE,
alpha = 0.8,
show.legend = FALSE
) +
scale_fill_gradient2(
low = htCol[1],
mid = "white",
high = htCol[2],
midpoint = 0
) +
scale_x_continuous(expand = c(0, 0), limits = c(1, n_genes)) +
scale_y_continuous(expand = c(0, 0)) +
theme_bw(base_size = 12) +
theme(
axis.ticks = element_blank(),
axis.text = element_blank(),
axis.title = element_blank(),
panel.grid = element_blank(),
plot.margin = unit(c(-.1, 0.2, 0, 0.2), "cm")
)
# 3. Ranked List Metric
prank <- ggplot(data.frame(x = 1:nrow(geneList_df), y = geneList_df$value), aes(x = x, y = y)) +
geom_col(aes(fill = y), width = 1, show.legend = FALSE) +
scale_fill_gradient2(
low = htCol[2],
mid = "white",
high = htCol[1],
midpoint = 0
) +
geom_hline(yintercept = 0, size = 0.5, color = "black", lty = "dashed") +
scale_x_continuous(expand = c(0, 0), limits = c(1, n_genes),
breaks = seq(0, n_genes, 200)) +
theme_bw(base_size = 12) +
theme(
panel.grid = element_blank(),
plot.margin = unit(c(-.1, 0.2, 0.2, 0.2), "cm")
) +
coord_cartesian(expand = 0) +
ylab("Ranked List Metric") +
xlab("Rank in Ordered Dataset")
# Combine plots
combined_plot <- aplot::plot_list(
gglist = list(pcurve, pseg_ht, prank),
ncol = 1,
heights = c(0.5, 0.2, 0.3)
)
# Display the combined plot
print(combined_plot)
## Warning: Removed 2 rows containing missing values (`geom_col()`).
