GGPLOT guide for R
By Etienne Aumont
- Packages used
- Data preparation
- Scatterplot
- Boxplots
- Barplot
- Multiple ROC plot
- Heatmaps
- Grid arrange and saving your file
library(ggplot2)
library(gplots)
library(grid) #to arrange several graphs into one grid
library(ggsignif) #for group comparison bars within plots
library(pROC) #for ROC
library(pheatmap) #for heatmaps
library(paletteer) #fun optional color palettes
library(RColorBrewer) #optional color package
library(ComplexHeatmap) #other heatmap package to install this one, you need to go through bioconductor:
if (!require("BiocManager", quietly = TRUE))
install.packages("BiocManager")
BiocManager::install("ComplexHeatmap")
I typically correct my variables for covariates used in the analyses that my graphs depict so that they are closer to the data in my analyses. There are 2 steps to this: 1) residualize for covariates and 2) Add mean of the original variable to residuals to obtain corrected values
Step 1 of covariate correction example
data_tPa$AB_Jack_resid2 <- residuals(lm(AB_Jack~ Age, data=data_tPa))
data_tPa$AB_Jack_resid2 <- residuals(lm(AB_Jack_resid2~ sex, data=data_tPa))
data_tPa$AB_Jack_resid2 <- residuals(lm(AB_Jack_resid2~ ApoE4, data=data_tPa))
If there are missing data, you may add ", na.action=na.exclude" at the end of each line. However, missing data has to be removed before running step 2
Step 2 of covariate correction example
M <- mean(data_tPa$AB_Jack)
data_tPa$AB_Jack_resid2 <- data_tPa$AB_Jack_resid2+M
p1 <- ggplot(data_tPa, aes(x=AB_Jack_resid2, y=Braak1_DP_resid2)) +
geom_point()+ #Here, points may be customized (shapes, color, etc), but variables for different shape/colors must be identified in the first line
geom_smooth(method='lm', color="black")+
#To add the regression line and the confidence interval
theme_classic() + #A minimalist theme, there are many others
ylab("Tau-PET rate of change") +xlab("Baseline global amyloid-PET SUVR")+ggtitle("Braak I")+
#To name the axies and the graph
theme(legend.position = "none")+
annotate(geom="text", x=1.4, y=0.4, label= "β = 0.096, p = 0.383", size = 5)+
#Adds text on the graph. The x and y coordinates for this need to be adjusted as needed
theme(plot.title = element_text(hjust = 0.5)) +
theme(plot.title = element_text(size = 20, face = "bold"), legend.title=element_text(size = 15, face="bold")
, legend.text=element_text(size = 15, face="bold")) +
theme(axis.text.x = element_text(color = "grey20", size = 15, angle = 0, hjust = .5, vjust = .5, face = "bold"),
axis.text.y = element_text(color = "grey20", size = 15, angle = 0, hjust = 1, vjust = 0, face = "bold"),
axis.title.x = element_text(color = "grey20", size = 16, angle = 0, hjust = .5, vjust = 0, face = "bold"),
axis.title.y = element_text(color = "grey20", size = 16, angle = 90, hjust = .5, vjust = 1, face = "bold"),
plot.tag = element_text(color = "grey20", size = 15, angle = 0, hjust = 1, vjust = 0, face = "bold"))+
#These are the characteristics of your text throughout the figure. You can adjust their position, color, etc.
coord_cartesian(ylim = c(-0.3, 0.4)) #if you need to specify the limits of your graph
Let's start by obtaining p-values for group comparisons
result <- signif(pvalue_AmyTau_fdr[12], digits = 2) #See the regression guide to store p-values in a matrix
df_p_val <- data.frame(
group1 = "0",
group2 = "1",
label = result
)
Then, creating the graph itself
p3 <- ggplot(data_no_NA_left, aes(x=Tau_status, y=L_CA2CA3_resid2)) +
geom_boxplot(outlier.shape = NA)+
#It is very important to set outliers as NA if you want to include individual datapoints! Otherwise, outliers will look like datapoints.
geom_jitter(size = 0.8)+
#Adding individual datapoints. The jitter is so they can be aligned with the boxes within the right group.
scale_x_discrete(labels = c("T-", "T+"))+
#Renaming the data from "0"s and "1"s into something more specific and informative
add_pvalue(df_p_val,
xmin = "group1",
xmax = "group2",
label = "label",
y.position = 172, label.size = 3.8)+
#Here, we use the data specified above
theme_classic() +ylab("Volume (adjusted voxel count)") +xlab("")+ggtitle("Left CA2/CA3")+
theme(legend.position = "none",
plot.tag = element_text(color = "grey20", size = 10, angle = 0, hjust = 1, vjust = 0, face = "bold"))+
theme(plot.title = element_text(hjust = 0.5)) +
theme(plot.title = element_text(size = 12, face = "bold"), legend.title=element_text(size = 10, face="bold")
, legend.text=element_text(size = 10, face="bold")) +
theme(axis.text.x = element_text(color = "grey20", size = 10, angle = 0, hjust = 0.5, vjust = 0.5, face = "bold"),
axis.text.y = element_text(color = "grey20", size = 10, angle = 0, hjust = 1, vjust = 0, face = "bold"),
axis.title.x = element_text(color = "grey20", size = 10, angle = 0, hjust = 0.5, vjust = 0, face = "bold"),
axis.title.y = element_text(color = "grey20", size = 10, angle = 90, hjust = 0.5, vjust = 2, face = "bold")) +
coord_cartesian(ylim = c(75, 175))
Here, we compare different groups to one another
p2<-ggplot(data, aes(x= reorder(CogABC, CogABCNum), y=iFilA.wt_res, fill=CogABCNum2)) +
#reorder is to bring back CogABC (grouping) in the right order so that "AD dementia" is not first. Fill is for the color of the bar
stat_summary(fun.data=mean_sdl, geom="bar") +
labs(tag = "B")+ #Add a label to the upper left corner of the plot
scale_fill_paletteer_d("beyonce::X54")+ #This is the color palette I chose to identify the bars (fill).
#There are thousands you can choose from in the paletteer package. See https://github.com/EmilHvitfeldt/paletteer for more details
#You can also explore color palettes and obtain codes for plots here: https://r-graph-gallery.com/color-palette-finder
scale_x_discrete(labels = c("Non-AD", "Preclinical AD", "Prodromal AD", "AD dementia"))+ #Renaming the 4 CogABC categories
geom_signif(comparisons = list(c("NonAD", "AD+")),annotation = c("*"),y_position = 3400, size = 1, textsize = 7)+
geom_signif(comparisons = list(c("NonAD", "MCI+")),map_signif_level=TRUE,y_position = 3000, size = 1, textsize = 7)+
geom_signif(comparisons = list(c("NonAD", "CN+")),map_signif_level=TRUE,y_position = 2600, size = 1, textsize = 7)+
#Adding lines between 2 bars to specify if they are significantly different to one another. The names in "" are the groups in your dataframe
geom_jitter(shape=18, size=2)+
#To overlay individual datapoints over the bars. You may resize them or change the shape.
stat_summary(fun.data=mean_cl_boot, geom="errorbar", width=0.3, size = .8)+
#To add errorbars
theme_classic() +ylab("Predicted iFLNA relative optical density\n") +xlab("clinicopathologic stages of AD")+
ggtitle("Insoluble FLNA by clinicopathologic stages of AD")+
annotate(geom="text", x=.8, y=3900, label= "ρ = .386*", size = 7)+
theme(legend.position = "none",
plot.tag = element_text(color = "grey20", size = 20, angle = 0, hjust = 1, vjust = 0, face = "bold"))+
theme(plot.title = element_text(hjust = 0.5)) +
theme(plot.title = element_text(size = 20, face = "bold"), legend.title=element_text(size = 15, face="bold")
, legend.text=element_text(size = 15, face="bold")) +
theme(axis.text.x = element_text(color = "grey20", size = 15, angle = 0, hjust = .5, vjust = .5, face = "bold"),
axis.text.y = element_text(color = "grey20", size = 15, angle = 0, hjust = 1, vjust = 0, face = "bold"),
axis.title.x = element_text(color = "grey20", size = 20, angle = 0, hjust = .5, vjust = 0, face = "bold"),
axis.title.y = element_text(color = "grey20", size = 18, angle = 90, hjust = .5, vjust = 0, face = "bold"))
p2<-ggplot(data, aes(x= reorder(CogABC, CogABCNum), y=iFilA.wt_res, fill=ApoE4)) +
stat_summary(fun.data=mean_sdl, geom="bar", position=position_dodge()) +
#position_dodge is to split the CogABC groups into the 2 fill categories (ApoE4 carriers or not)
labs(tag = "B")+
scale_fill_paletteer_d("beyonce::X54",name = "APOE ε4" ,breaks=c("0", "1"),
labels=c("Noncarrier", "Carrier"))+
#to define the color of the fill (ApoE) categories and switch the names from 0 or 1 to a proper name for the legend
scale_x_discrete(labels = c("Non-AD", "Preclinical AD", "Prodromal AD", "AD dementia"))+ #Renaming the 4 CogABC categories
geom_point(shape=18, size=2, position=position_jitterdodge(jitter.width = .6, dodge.width = .9))+
#position_jitterdodge is to split the individual data points of each fill categories into different columns
stat_summary(fun.data=mean_cl_boot, geom="errorbar", position=position_dodge(.9), width=0.3, size = .8)+
theme_classic() +ylab("Predicted iFLNA relative optical density\n") +xlab("clinicopathologic stages of AD")+
ggtitle("Insoluble FLNA by clinicopathologic stages & APOE")+
theme(legend.position = "right",
plot.tag = element_text(color = "grey20", size = 20, angle = 0, hjust = 1, vjust = 0, face = "bold"))+
theme(plot.title = element_text(hjust = 0.5)) +
theme(plot.title = element_text(size = 20, face = "bold"), legend.title=element_text(size = 15, face="bold")
, legend.text=element_text(size = 15, face="bold")) +
theme(axis.text.x = element_text(color = "grey20", size = 15, angle = 0, hjust = .5, vjust = .5, face = "bold"),
axis.text.y = element_text(color = "grey20", size = 15, angle = 0, hjust = 1, vjust = 0, face = "bold"),
axis.title.x = element_text(color = "grey20", size = 20, angle = 0, hjust = .5, vjust = 0, face = "bold"),
axis.title.y = element_text(color = "grey20", size = 18, angle = 90, hjust = .5, vjust = 0, face = "bold"))
Default ROC plots with pROC can be inputted into ggroc, which generates objects of the saame nature as ggplot outputs. Step 1: Create the ROC objects and put them into a list
roci<-list()
roci[["MCI"]] <- roc(data$MCIAUC,data$iFilA.wt_res)
roci[["NCI"]] <- roc(data$CNAUC,data$iFilA.wt_res)
roci[["All"]] <- roc(data$ABCAUC,data$iFilA.wt_res)
Step 2: Generate the plot using the list of ROCs
p1<-ggroc(roci, size=1)+
geom_abline(intercept=1,slope=1)+
labs(tag = "A")+
theme_classic() +ggtitle("AD detection by insoluble Filamin A")+
theme(plot.tag = element_text(color = "grey20", size = 20, angle = 0, hjust = 1, vjust = 0, face = "bold"))+
theme(plot.title = element_text(hjust = 0.5)) +
annotate(geom="text", x=0.122, y=.54, label= "AUC = .818*", size = 7)+
annotate(geom="text", x=0.133, y=.495, label= "AUC = .556", size = 7)+
annotate(geom="text", x=0.12, y=.45, label= "AUC = .727*", size = 7)+ #coordinates will need to be adjusted depending on the figure size
theme(plot.title = element_text(size = 20, face = "bold"), legend.title=element_blank(), legend.text=element_text(size = 20)) +
theme(axis.text.x = element_text(color = "grey20", size = 15, angle = 0, hjust = .5, vjust = .5, face = "bold"),
axis.text.y = element_text(color = "grey20", size = 15, angle = 0, hjust = 1, vjust = 0, face = "bold"),
axis.title.x = element_text(color = "grey20", size = 20, angle = 0, hjust = .5, vjust = 0, face = "bold"),
axis.title.y = element_text(color = "grey20", size = 20, angle = 90, hjust = .5, vjust = 0, face = "bold"))
This plot is different because it plots values obtained from the data. They are not found in the data itself. Step 1 is to created a CSV file of regression results that looks like this:
Amyloid-PET Braak 1 tau Braak 2 tau Braak 3 tau
Left CA1 -0.03904966 -0.06382974 -0.03046591 0.05466574
Right CA1 -0.23329494 -0.21828464 -0.16341319 -0.09826552
Left CA2/CA3 -0.0940731 -0.09515193 -0.07504647 -0.16213403
Right CA2/CA3 -0.1780232 -0.1194428 -0.1133248 -0.18880813
Left DG -0.21374795 -0.21778841 -0.23151909 -0.21443269
Right DG -0.3112549 -0.29603521 -0.28034264 -0.29740854
Left sub -0.06225126 -0.13241333 -0.1837268 -0.15344105
Right sub -0.32845427 -0.32064072 -0.32701141 -0.26208623
Left SRLM -0.15826004 -0.20863016 -0.17782254 -0.19753768
Right SRLM -0.33523554 -0.28892321 -0.22910347 -0.24390613
Step 2 is to load and format the data so that the column names are properly recognized as such
RegMat1<- read.csv('/Users/eaumo/Desktop/Labo_PRN/Article 2/Figures, tables & supplementary material/Hidden/Reg_CS.csv')
RegMat1<- column_to_rownames(RegMat1, var="X")
RegMat1 <- as.matrix(RegMat1)
colnames(RegMat1) <- gsub("\\.", " ", colnames(RegMat1))
The next step is to create a color scale. I had to do a lot of trial and error here, with a lot of weird things happening to the color intervals, so this section might not be ideal.
paletteLength <- 8
myColor <- colorRampPalette(c("red", "white", "white", "blue"))(paletteLength)
myBreaks <- c(seq(-0.5, -0.1, length.out=ceiling(5)), seq(-0.099, 0.099, length.out=1),
seq(0.1, 0.5, length.out=floor(5)))
The last step is to put all of the elements together. Here I used ComplexHeatmaps because of the extra functionality, and I used it using the syntax of the pheatmap package because it offers additional options. Cell size was fixed to better adjust the figure size.
p1<-ComplexHeatmap::pheatmap(RegMat1, cluster_rows = FALSE, cluster_cols = FALSE, color = myColor,
breaks = myBreaks, display_numbers = TRUE, fontsize_number = 12,
column_names_side = c("top"), row_names_side = c("left"),
name = "STD Beta", fontsize = 12, border_color = NA,
cellwidth = 50, cellheight = 40,
main = " \na) Baseline PET SUVR with baseline \n hippocampal subfield volumes ")
Here, numbers 1 to 4 in the layout matrix are associated with the rank of the plots loaded (p1 = 1, p2 = 4, etc.)
p0 <- grid.arrange(p1, p4, p3, p2,
widths = c(1,1.2), #to customize the width of the plots. Here, column 2 will be 20% wider than column 1
layout_matrix = rbind(c(1, 2),
c(3, 4)),
top = text_grob("Regression of tau-PET rates of change with baseline global amyloid-PET ", size = 25)
)
ggsave('/path/to/the/figure/folder/Fig2.tiff', plot = p0, scale = 2, width = 1800, height = 2100, units = c("px"), dpi = 300)
#This ggsave is particularly useful when you need a high DPI (for publication for example)
Heatmap objects must first be transformed into grob objects to fit into a grid.arrange
grob1 = grid.grabExpr(draw(p1))
grob2 = grid.grabExpr(draw(p2))
grob3 = grid.grabExpr(draw(p3))
p0 <- grid.arrange(grob1, grob2, grob3,
layout_matrix = rbind(c(1, 2, 3)),
top = text_grob("Regression heatmap of tau and amyloid-PET with hippocampal subfields", size = 20)
)
