<- read.table(file = "../data/counts_raw.txt", sep = "\t", header = T)
gc_raw <- read.table(file = "../data/counts_filtered.txt", sep = "\t", header = T)
gc_filt <- read.table(file = "../data/counts_vst.txt", sep = "\t", header = T)
gc_vst <- read.table(file = "../data/counts_deseq2.txt", sep = "\t", header = T)
gc_deseq <- read.table("../data/metadata.csv", header = T, sep = ";")
md rownames(md) <- md$Sample_ID
<- function(x) sqrt(var(x)/length(x))
se
<-
gene_counts_all %>% gather(Sample_ID, Raw, -Gene) %>%
gc_raw full_join(gc_filt %>% gather(Sample_ID, Filtered, -Gene), by = c("Gene", "Sample_ID")) %>%
full_join(gc_vst %>% gather(Sample_ID, VST, -Gene), by = c("Gene", "Sample_ID")) %>%
full_join(gc_deseq %>% gather(Sample_ID, DESeq2, -Gene), by = c("Gene", "Sample_ID")) %>%
gather(Method, count, Raw:DESeq2) %>%
filter(!is.na(count)) %>%
full_join(md, by = "Sample_ID")
$Time <- factor(gene_counts_all$Time, levels = c("t0","t2","t6","t24"))
gene_counts_all$Replicate <- factor(gene_counts_all$Replicate, levels = c("A","B","C"))
gene_counts_all$Method <- factor(gene_counts_all$Method, levels = c("Raw","Filtered","DESeq2","VST"))
gene_counts_all
%>%
gene_counts_all group_by(Time, Replicate, Method) %>%
summarise(mean=mean(log10(count +1)),se=se(log10(count +1))) %>%
ggplot(aes(x= Time, y= mean, fill = Replicate)) +
geom_bar(position = position_dodge2(), stat = "identity") +
geom_errorbar(aes(ymin=mean-se, ymax=mean+se), position = position_dodge2(.9, padding = .6)) +
facet_wrap(~Method, scales = "free")
The solutions given below are just one way of obtaining the desired plots! There are probably several different ways you could code to get the same plots.
ggplot basics
Exercise I and II
For these exercises the solutions are already part of the material, all you need to do is to replace the filename counts_raw.txt
with each of the different files to look at the differences between different normalization methods.
Exercise III
Task 3.1
Task 3.2
%>%
gene_counts_all group_by(Time, Replicate, Method) %>%
ggplot() +
geom_boxplot(mapping = aes(x = Sample_Name, y = log10(count + 1), fill = Time)) +
facet_wrap(~Method*Replicate, ncol = 3, scales = "free")
Advanced ggplot
Exercise I
<- gc_raw %>%
gc_long gather(Sample_ID, count, -Gene) %>%
full_join(md, by = "Sample_ID") %>%
select(Sample_ID, everything()) %>%
select(-c(Gene,count), c(Gene,count))
$Sample_Name <- factor(gc_long$Sample_Name, levels = c("t0_A","t0_B","t0_C","t2_A","t2_B","t2_C","t6_A","t6_B","t6_C","t24_A","t24_B","t24_C"))
gc_long$Time <- factor(gc_long$Time, levels = c("t0","t2","t6","t24"))
gc_long$Replicate <- factor(gc_long$Replicate, levels = c("A","B","C"))
gc_long
%>%
gc_long group_by(Time, Replicate) %>%
summarise(mean=mean(log10(count +1)),se=se(log10(count +1))) %>%
ggplot(aes(x=Time, y=mean, color = Replicate)) +
facet_wrap(~Replicate) +
geom_line(aes(group=1), stat= "identity", size = 2) +
scale_x_discrete(limits= c("t0", "t2", "t24")) +
scale_y_continuous(limits = c(0.4,0.8), breaks = seq(0.4,0.8,0.05)) +
guides(color="none") +
ylab(label = "mean(log10(count + 1))") +
theme_light() +
theme(axis.text = element_text(face="bold", size=12),
axis.title = element_text(face="bold", color = "#C84DF9", size=14),
axis.ticks = element_blank())
Exercise II
library(ggpubr)
<- ggplot(data=iris,mapping=aes(x=Sepal.Length, y = Sepal.Width, color = Species)) +
p4 geom_point(size = 3, alpha = 0.6) +
theme_classic(base_size = 12) +
border()
<- ggplot(data=iris,mapping=aes(Sepal.Length, fill = Species)) +
d1 geom_density(alpha = 0.6) +
theme_classic() +
clean_theme() +
theme(legend.position = "none")
<- ggplot(data=iris,mapping=aes(Sepal.Width, fill = Species)) +
d2 geom_density(alpha = 0.6) +
theme_classic() +
clean_theme() +
theme(legend.position = "none") +
rotate()
ggarrange(d1, NULL, p4, d2,
ncol = 2, nrow = 2, align = "hv",
widths = c(3, 1), heights = c(1, 3),
common.legend = TRUE)
PCA and Gene Expression
Exercise I
row.names(gc_vst) <- gc_vst$Gene
$Gene <- NULL
gc_vst<- dist(t(gc_vst))
gc_dist <- cmdscale(gc_dist,eig=TRUE, k=2)
gc_mds <- gc_mds$eig
Eigenvalues <- Eigenvalues / sum(Eigenvalues)
Variance <- 100 * signif(Variance[1], 3)
Variance1 <- 100 * signif(Variance[2], 3)
Variance2
<- gc_mds$points %>%
gc_mds_long as.data.frame() %>%
rownames_to_column("Sample_ID") %>%
full_join(md, by = "Sample_ID")
$Sample_Name <- factor(gc_mds_long$Sample_Name, levels = c("t0_A","t0_B","t0_C","t2_A","t2_B","t2_C","t6_A","t6_B","t6_C","t24_A","t24_B","t24_C"))
gc_mds_long$Time <- factor(gc_mds_long$Time, levels = c("t0","t2","t6","t24"))
gc_mds_long$Replicate <- factor(gc_mds_long$Replicate, levels = c("A","B","C"))
gc_mds_long
ggplot(gc_mds_long, aes(x=V1, y=V2, color = Time)) +
geom_point(size = 3, aes(shape = Replicate)) +
xlab(paste("PCO1: ", Variance1, "%")) +
ylab(paste("PCO2: ", Variance2, "%")) +
geom_vline(xintercept = 0, linetype=2) +
geom_hline(yintercept = 0, linetype=2) +
theme_bw() +
theme(panel.grid.major = element_blank(),
panel.grid.minor = element_blank())
Exercise II
Task 2.1
<- read.table("../data/Time_t2_vs_t0.txt", sep = "\t", header = TRUE, row.names = 1)
t2_vs_t0
%>%
t2_vs_t0 ggplot() +
geom_point(aes(x = baseMean, y = log2FoldChange), color = "grey70") +
geom_point(data=filter(t2_vs_t0, padj < 0.05), aes(x = baseMean, y = log2FoldChange), color = "blue") +
geom_hline(yintercept = 0) +
scale_x_log10("Mean of normalized counts") +
ylab("log fold change") +
theme_bw(base_size = 14)
Task 2.2
library(ggrepel)
<- read.table("../data/human_biomaRt_annotation.csv", sep = ";", header = TRUE, row.names = 1)
gene_info names(gene_info) <- c('gene_id', 'gene_name')
<- t2_vs_t0 %>%
de_w_names rownames_to_column('gene_id') %>%
left_join(gene_info, by = "gene_id")
ggplot(de_w_names, aes(x = log2FoldChange, y = -log10(padj))) +
geom_point(color = "grey70") +
geom_text_repel(data=filter(de_w_names, padj < 0.05 & abs(log2FoldChange) > 1.5), aes(x = log2FoldChange, y = -log10(padj), label=gene_name)) +
geom_hline(yintercept = 1.3, linetype = 2) +
geom_point(data=filter(t2_vs_t0, padj < 0.05), color = "red") +
ylab("Adjusted P-values in -log10") +
theme_bw(base_size = 14)
Dynamic Plotting
Exercise 2.1
library(palmerpenguins)
library(leaflet)
leaflet() %>%
addTiles() %>%
#addProviderTiles(providers$Esri.WorldImagery) %>%
setView(lng = -65.5, lat = -65.5, zoom = 6) %>%
addCircleMarkers(data = island_coordinates, popup = c("Biscoe", "Dream", "Torgersen"), color = "green")
Exercise 3.1
library(ggplot2)
library(plotly)
plot_ly(msleep, x = ~sleep_rem, y = ~sleep_total,
color = ~vore, colors = "Set1", size=2, type = "scatter", mode = "markers")
plot_ly(msleep, x = ~sleep_rem, y = ~sleep_total,
color = ~vore, colors = "Set1", size=2, type = "scatter", mode = "lines", line = list(width = 4))
library(ggplot2)
library(plotly)
<- plot_ly(msleep, x = ~vore, color = ~vore, colors = c("#88CCEE", "#CC6677", "#DDCC77", "#117733")) %>%
p1 add_histogram()
<- plot_ly(msleep, x = ~vore, y = ~sleep_total, color = ~vore, type = "box", size=2, colors = c("#88CCEE", "#CC6677", "#DDCC77", "#117733"), showlegend = F)
p2
<- subplot(p1, p2)
fig
%>%
fig layout(title = "Mammals Sleep Data",
xaxis = list(title = "Diet"),
yaxis = list(title = "Number of mammals"),
xaxis2 = list(title = "Diet"),
yaxis2 = list(title = "Total sleep (hours)"),
legend=list(title=list(text="Diet")))
Exercise 4.1
library(ggiraph)
library(ggplot)
<- ggplot(msleep,aes(x = vore, y = sleep_total, color = vore, tooltip = vore)) +
PlotData geom_boxplot_interactive() +
theme_minimal() +
scale_color_manual(values = c("#88CCEE", "#CC6677", "#DDCC77", "#117733", "#999999")) +
labs(title = "Mammals Sleep Data", x = "Diet", y = "Total Sleep (hours)", color = "Diet")
girafe(ggobj = PlotData)
Exercise 5.1
library(DT)
<- read.table("../../data/counts_raw.txt", header = T, row.names = 1, sep = "\t")
gc <- read.table("../../data/metadata.csv", header = T, sep = ";")
md
datatable(gc, options = list(pageLength = 10), caption = "Table 1: Raw gene counts") %>%
formatStyle(
'Sample_1',
background = styleColorBar(gc$Sample_1, 'steelblue'),
backgroundSize = '100% 90%',
backgroundRepeat = 'no-repeat',
backgroundPosition = 'center') %>%
formatStyle('Sample_3', color = styleInterval(c(8500, 15000), c('white', 'blue', 'green')),
backgroundColor = styleInterval(11000, c('gray', 'orange')))
datatable(md, options = list(pageLength = 10, order = list(list(4, 'desc'))), caption = "Table 1: Metadata")
Exercise 6.1
library(plotly)
library(ggplot2)
library(crosstalk)
<- SharedData$new(msleep)
shared_msleep
<- ggplot(data = shared_msleep, aes(x = sleep_rem, y = sleep_total, color = vore)) +
gg_plot1 geom_point() + theme_bw()
<- ggplotly(gg_plot1)
plotly_plot1
bscols(widths = c(3, NA),
list(
filter_checkbox("vore", "Diet", shared_msleep, ~vore, inline = TRUE),
filter_slider("sleep_total", "Total Sleep (hours)", shared_msleep, ~sleep_total, width = "100%"),
filter_select("conservation", "Conservation Status", shared_msleep, ~conservation)
), plotly_plot1)
Exercise 7.1
---
title: "OJS example"
author: "John Doe"
date: last-modified
date-format: "YYYY-MMM-DD"
format: html
---
```{r}
library(ggplot2)
```
```{r}
dt_sleep <- msleep
ojs_define(ojsd = dt_sleep)
```
```{ojs}
ojsdata = transpose(ojsd)
```
```{ojs}
viewof raw_table = Inputs.table(ojsdata)
```
```{ojs}
viewof diets = Inputs.checkbox(
["herbi", "carni", "omni", "insecti"],
{ value: ["herbi", "carni"],
label: "Diet:"
}
)
```
```{ojs}
filtered = ojsdata.filter(function(mammal) {
return diets.includes(mammal.vore);
})
```
```{ojs}
Plot.plot({
marks: [
Plot.dot(filtered, {
x:"sleep_rem",
y: "sleep_total",
fill: "vore"
})
],
grid: true
})