Last updated: 2019-06-10
workflowr checks: (Click a bullet for more information) ✔ R Markdown file: up-to-date
Great! Since the R Markdown file has been committed to the Git repository, you know the exact version of the code that produced these results.
✔ Environment: empty
Great job! The global environment was empty. Objects defined in the global environment can affect the analysis in your R Markdown file in unknown ways. For reproduciblity it’s best to always run the code in an empty environment.
✔ Seed:
set.seed(20190513)
The command set.seed(20190513)
was run prior to running the code in the R Markdown file. Setting a seed ensures that any results that rely on randomness, e.g. subsampling or permutations, are reproducible.
✔ Session information: recorded
Great job! Recording the operating system, R version, and package versions is critical for reproducibility.
✔ Repository version: 028d3cc
wflow_publish
or wflow_git_commit
). workflowr only checks the R Markdown file, but you know if there are other scripts or data files that it depends on. Below is the status of the Git repository when the results were generated:
Ignored files:
Ignored: .Rhistory
Ignored: .Rproj.user/
Ignored: data/NAPA_clim_U.Rda
Ignored: data/NAPA_clim_V.Rda
Ignored: data/NAPA_clim_W.Rda
Ignored: data/NAPA_clim_emp_ice.Rda
Ignored: data/NAPA_clim_emp_oce.Rda
Ignored: data/NAPA_clim_fmmflx.Rda
Ignored: data/NAPA_clim_mldkz5.Rda
Ignored: data/NAPA_clim_mldr10_1.Rda
Ignored: data/NAPA_clim_qemp_oce.Rda
Ignored: data/NAPA_clim_qla_oce.Rda
Ignored: data/NAPA_clim_qns.Rda
Ignored: data/NAPA_clim_qsb_oce.Rda
Ignored: data/NAPA_clim_qt.Rda
Ignored: data/NAPA_clim_runoffs.Rda
Ignored: data/NAPA_clim_ssh.Rda
Ignored: data/NAPA_clim_sss.Rda
Ignored: data/NAPA_clim_sst.Rda
Ignored: data/NAPA_clim_taum.Rda
Ignored: data/NAPA_clim_vars.Rda
Ignored: data/NAPA_clim_vecs.Rda
Ignored: data/node_mean_all_anom.Rda
Ignored: data/packet_all_anom.Rda
Ignored: data/som_all_anom.Rda
Ignored: data/synoptic_states.Rda
Ignored: data/synoptic_vec_states.Rda
Untracked files:
Untracked: output/som_plot_uoce_anom.pdf
Untracked: output/som_plot_voce_anom.pdf
Unstaged changes:
Modified: analysis/MHWNWA.bib
Modified: analysis/_site.yml
Modified: code/workflow.R
Modified: output/som_plot_mldr10_1_anom.pdf
Modified: output/som_plot_qt_anom.pdf
Modified: output/som_plot_sst_anom.pdf
Modified: output/som_plot_taum_anom.pdf
Note that any generated files, e.g. HTML, png, CSS, etc., are not included in this status report because it is ok for generated content to have uncommitted changes.
In this final vignette we will go over the creation of the figures used in the publication for this research. These figures are largely adapted from the techniques seen in Oliver et al. (2018) and Schlegel et al. (2017).
# Insatll from GitHub
# .libPaths(c("~/R-packages", .libPaths()))
# devtools::install_github("fabrice-rossi/yasomi")
# Packages used in this vignette
library(jsonlite, lib.loc = "../R-packages/")
library(tidyverse) # Base suite of functions
library(ncdf4) # For opening and working with NetCDF files
library(lubridate) # For convenient date manipulation
# library(scales) # For scaling data before running SOM
library(yasomi, lib.loc = "../R-packages/") # The SOM package of choice due to PCI compliance
library(data.table) # For working with massive dataframes
# Set number of cores
doMC::registerDoMC(cores = 50)
# Disable scientific notation for numeric values
# I just find it annoying
options(scipen = 999)
# Set number of cores
doMC::registerDoMC(cores = 50)
# Disable scientific notation for numeric values
# I just find it annoying
options(scipen = 999)
# Individual regions
NWA_coords <- readRDS("data/NWA_coords_cabot.Rda")
# The NAPA variables
NAPA_vars <- readRDS("data/NAPA_vars.Rda")
# Corners of the study area
NWA_corners <- readRDS("data/NWA_corners.Rda")
# Create smaller corners to use less RAM
# This also better matches the previous South African work
# The Tasmania work had corners of roughly 2 degrees greater than the study area
NWA_corners_sub <- c(NWA_corners[1]+8, NWA_corners[2]-8, NWA_corners[3]+8, NWA_corners[4]-8)
# The base map
map_base <- ggplot2::fortify(maps::map(fill = TRUE, col = "grey80", plot = FALSE)) %>%
dplyr::rename(lon = long) %>%
mutate(group = ifelse(lon > 180, group+9999, group),
lon = ifelse(lon > 180, lon-360, lon)) %>%
select(-region, -subregion)
###########################################################################
### "3.Figures.R"
## This script creates the figures for the paper and supplemental
# 1. Load all libraries and functions used in this script
# 2. Create synoptic figure for each event
# 3. Create synoptic figure showing SOM nodes
# 4. Create lolliplots for the SOM nodes
# 5. Create dendrogram for HCA results
# 6. Create ordiplot for MDS results
# 7. Create map of study area
#############################################################################
# 1. Load all libraries and functions used in this script -----------------
source("func/synoptic.func.R")
source("func/som.func.R")
source("func/scale.bar.func.R")
library(vegan)
library(ggdendro)
library(broom)
# 2. Create synoptic figure for each event -------------------------------
# Load SACTN data
load("~/data/SACTN/AHW/SACTN_clims.Rdata")
load("data/SACTN/SACTN_events.Rdata")
load("setupParams/SACTN_site_list.Rdata")
# The files for loading
event_idx <- data.frame(event = dir("data/SOM", full.names = TRUE),
x = length(dir("data/SOM")))
# Create a synoptic atlas figure for each MHW
system.time(plyr::ddply(event_idx, c("event"), synoptic.fig, .progress = "text")) # 539 seconds
# 3. Create synoptic figure showing SOM nodes -----------------------------
load("data/node_means.Rdata")
load("data/node_all_anom.Rdata")
all.panels(node_means, node_all_anom)
# 4. Create lolliplots for the SOM nodes ----------------------------------
# Load data for figure
load("data/SACTN/SACTN_events.Rdata")
load("data/node_all_anom.Rdata")
# Merge into one dataframe
node_all <- merge(node_all_anom, SACTN_events, by = c("event", "site", "season", "event_no"))
# Calculate mean and median per node for plotting
node_h_lines <- node_all %>%
group_by(node) %>%
summarise(mean_int_cum = mean(int_cum, na.rm = T),
median_int_cum = median(int_cum, na.rm = T))
# Create the figure
ggplot(data = node_all, aes(x = date_start, y = int_cum)) +
geom_lolli() +
geom_point(aes(colour = season)) +
geom_label(aes(x = as.Date("2005-01-01"), y = 580, label = paste0("n = ", count,"/",length(node))),
size = 3, label.padding = unit(0.5, "lines")) +
geom_hline(data = node_h_lines, aes(yintercept = mean_int_cum), linetype = "dashed") +
geom_hline(data = node_h_lines, aes(yintercept = median_int_cum), linetype = "dotted") +
facet_wrap(~node) +
labs(x = "", y = "Cummulative intensity (°C·days)", colour = "Season") +
theme_grey() +
# scale_y_continuous(expand = c(0, 100)) +
theme(strip.background = element_rect(fill = NA),
panel.border = element_rect(fill = NA, colour = "black", size = 1),
axis.text = element_text(size = 12, colour = "black"),
axis.ticks = element_line(colour = "black"))
ggsave("graph/SOM_lolli.pdf", height = 9, width = 10)
# 5. Create dendrogram for HCA results ------------------------------------
# Load data
load("data/all_anom_hclust.Rdata")
load("data/all_anom_env.Rdata")
# Prep dendrogram
dhc <- as.dendrogram(all_anom_hclust)
ddata <- dendro_data(dhc, type = "rectangle")
label_data <- label(ddata)
label_data$season <- all_anom_env$season[as.numeric(as.character(label_data$label))]
label_data$type <- all_anom_env$type[as.numeric(as.character(label_data$label))]
# Plot the dendrogram
ggplot(segment(ddata)) +
geom_segment(aes(x = x, y = y, xend = xend, yend = yend)) +
geom_point(data = label_data, aes(x = x, y = y, shape = type, colour = season)) +
scale_shape_manual(values = c(1, 16)) +
theme_grey() +
theme(panel.border = element_rect(fill = NA, colour = "black", size = 1),
axis.text = element_text(size = 12, colour = "black"),
axis.ticks = element_line(colour = "black"))
ggsave("graph/HCA.pdf", height = 9, width = 9)
# 6. Create ordiplot for MDS results --------------------------------------
# Load data
load("data/all_anom_MDS.Rdata")
load("data/all_anom_env.Rdata")
load("data/SACTN/SACTN_events.Rdata")
SACTN_events$type <- NULL
# load("data/node_means.Rdata")
# load("data/node_all_anom.Rdata")
## Fit environmental variables
ord_fit <- envfit(all_anom_MDS ~ season + type, data = all_anom_env)
# ord_fit
ord_fit_df <- as.data.frame(ord_fit$factors$centroids)
ord_fit_df$factors <- c("autumn", "spring", "summer", "winter", "clim", "MHW")
# Merge event values
# SACTN_events_sub <- SACTN_events[,c(1,2)]
all_anom_env <- left_join(all_anom_env, SACTN_events, by = c("event", "season"))
# Create MDS dataframe
# mds_df <- data.frame(all_anom_MDS$points, type = all_anom_env$type,
# event = all_anom_env$event, season = all_anom_env$season)
mds_df <- data.frame(all_anom_MDS$points, all_anom_env)
mds_df$duration[is.na(mds_df$duration)] <- mean(mds_df$duration, na.rm = T)
# Plot the fits
ggplot(data = mds_df, aes(x = MDS1, y = MDS2)) +
geom_point(aes(colour = season, shape = type, size = duration)) +
geom_segment(data = ord_fit_df, aes(x = 0, y = 0, xend = NMDS1, yend = NMDS2),
arrow = arrow(angle = 40, length = unit(0.2, "cm"), type = "open"),
alpha = 1, colour = "black", size = 0.5) +
geom_text(data = ord_fit_df, aes(label = factors, x = NMDS1, y = NMDS2), size = 8) +
scale_shape_manual(name = "State", values = c(19, 15), labels = c("clim", "MHW")) +
scale_colour_discrete(name = "Season") +
scale_size_continuous(name = "Duration\n(days)", breaks = c(20, 70, 120, 170, 220)) +
guides(colour = guide_legend(order = 1),
shape = guide_legend(order = 2),
size = guide_legend(override.aes = list(shape = 15), order = 3)) +
# labs(size = "Duration") +
theme_grey() +
theme(strip.background = element_rect(fill = NA),
panel.border = element_rect(fill = NA, colour = "black", size = 1),
axis.text = element_text(size = 12, colour = "black"),
axis.ticks = element_line(colour = "black"))
ggsave("graph/MDS.pdf", height = 9, width = 12)
# Linear model: duration vs. ordination from centre point
lm_all_results <- mds_df %>%
na.omit() %>%
mutate(MDS2 = abs(MDS2)) %>%
select(-MDS1, -event, -season, -type, -coast, -site) %>%
gather(key = group,
value = measurement,
-MDS2) %>%
group_by(group) %>%
nest() %>%
mutate(model = purrr::map(data, ~lm(measurement ~ MDS2, data = .))) %>%
unnest(model %>% purrr::map(glance)) %>%
select(-data, -model)
lm_all_results[lm_all_results$adj.r.squared == max(lm_all_results$adj.r.squared, na.rm = T),]
# Visualization of the best linear model
ggplot(data = mds_df) +
geom_smooth(aes(x = duration, y = abs(MDS2)), method = "lm", se = F) +
theme_grey()
# 7. Create map of study area ---------------------------------------------
## Load data
# International borders
load("graph/africa_borders.Rdata")
# Hi-res bathy
load("graph/bathy.Rdata")
# Reanalysis data
load("data/ERA/ERA_temp_clim.Rdata")
load("data/ERA/ERA_u_clim.Rdata")
colnames(ERA_u_clim)[4] <- "val"
ERA_u_clim$variable <- "u"
load("data/ERA/ERA_v_clim.Rdata")
colnames(ERA_v_clim)[4] <- "val"
ERA_v_clim$variable <- "v"
# Remote data
load("data/OISST/OISST_temp_clim.Rdata")
load("data/AVISO/AVISO_u_clim.Rdata")
colnames(AVISO_u_clim)[4] <- "val"
AVISO_u_clim$variable <- "u"
load("data/AVISO/AVISO_v_clim.Rdata")
colnames(AVISO_v_clim)[4] <- "val"
AVISO_v_clim$variable <- "v"
# In situ time series locations
load("setupParams/SACTN_site_list.Rdata")
SACTN_site_list$order <- 1:nrow(SACTN_site_list)
# Create annual mean air-sea state
# sea_temp <- filter(OISST_temp_clim, date == "01-15") %>%
sea_temp <- data.table::data.table(OISST_temp_clim)
sea_temp <- sea_temp[, .(temp = mean(temp, na.rm = TRUE)),
by = .(x, y)] %>%
rename(lon = x, lat = y)
currents <- data.table::data.table(rbind(AVISO_u_clim, AVISO_v_clim))
currents <- currents[, .(val = mean(val, na.rm = TRUE)),
by = .(x, y, variable)] %>%
spread(key = variable, value = val) %>%
rename(lon = x, lat = y) %>%
mutate(arrow_size = ((abs(u*v)/ max(abs(u*v)))+0.2)/6)
air_temp <- data.table::data.table(ERA_temp_clim)
air_temp <- air_temp[, .(temp = mean(temp, na.rm = TRUE)),
by = .(x, y)] %>%
rename(lon = x, lat = y)
winds <- data.table::data.table(rbind(ERA_u_clim, ERA_v_clim))
winds <- winds[, .(val = mean(val, na.rm = TRUE)),
by = .(x, y, variable)] %>%
spread(key = variable, value = val) %>%
rename(lon = x, lat = y) %>%
mutate(arrow_size = ((abs(u*v)/ max(abs(u*v)))+0.3)/6)
# Reduce wind/ current vectors
lon_sub <- seq(10, 40, by = 1)
lat_sub <- seq(-40, -15, by = 1)
# currents <- currents[(currents$lon %in% lon_sub & currents$lat %in% lat_sub),]
winds <- winds[(winds$lon %in% lon_sub & winds$lat %in% lat_sub),]
# Establish the vector scalar for the currents
current_uv_scalar <- 2
# Establish the vector scalar for the wind
wind_uv_scalar <- 0.5
# Wind feature vector coordinates
cyc_atlantic <- data.frame(x = c(14.0, 16.1, 16.0), y = c(-36.0, -34.4, -32.1), xend = c(16.0, 16.1, 14.0), yend = c(-34.5, -32.2, -30.6))
cyc_indian <- data.frame(x = c(36.0, 33.9, 34.0), y = c(-31.5, -33.1, -35.4), xend = c(34.0, 33.9, 36.0), yend = c(-33.0, -35.3, -36.9))
westerlies <- data.frame(x = c(18.0, 21.1, 24.2), y = c(-38.0, -37.8, -37.8), xend = c(21.0, 24.1, 27.2), yend = c(-37.8, -37.8, -38.0))
# The top figure (sea)
fig_1_top <- ggplot(data = southern_africa_coast, aes(x = lon, y = lat)) +
# The ocean temperature
geom_raster(data = sea_temp, aes(fill = temp)) +
# The bathymetry
stat_contour(data = bathy[bathy$depth < -100 & bathy$depth > -300,],
aes(x = lon, y = lat, z = depth), alpha = 0.5,
colour = "ivory", size = 0.5, binwidth = 200, na.rm = TRUE, show.legend = FALSE) +
# The current vectors
geom_segment(data = currents, aes(xend = lon + u * current_uv_scalar, yend = lat + v * current_uv_scalar),
arrow = arrow(angle = 40, length = unit(currents$arrow_size, "cm"), type = "open"),
linejoin = "mitre", size = 0.4) +
# The land mass
geom_polygon(aes(group = group), fill = "grey70", colour = "black", size = 0.5, show.legend = FALSE) +
geom_path(data = africa_borders, aes(group = group)) +
# The legend for the vector length
geom_label(aes(x = 37.0, y = -38.0, label = "1.0 m/s\n"), size = 3, label.padding = unit(0.5, "lines")) +
geom_segment(aes(x = 36.0, y = -38.3, xend = 38.0, yend = -38.3), linejoin = "mitre",
arrow = arrow(angle = 40, length = unit(0.2, "cm"), type = "open")) +
# The in situ sites
geom_point(data = SACTN_site_list, shape = 19, size = 2.8, colour = "ivory") +
geom_text(data = SACTN_site_list[-c(3,4,7:9,18,21,23:24),], aes(label = order), size = 1.9, colour = "red") +
# Oceans
annotate("text", label = "INDIAN\nOCEAN", x = 37.00, y = -34.0, size = 4.0, angle = 0, colour = "ivory") +
annotate("text", label = "ATLANTIC\nOCEAN", x = 13.10, y = -34.0, size = 4.0, angle = 0, colour = "ivory") +
# Benguela
geom_segment(aes(x = 17.2, y = -32.6, xend = 15.2, yend = -29.5),
arrow = arrow(length = unit(0.3, "cm")), size = 0.5, colour = "ivory") +
annotate("text", label = "Benguela", x = 16.0, y = -31.8, size = 3.5, angle = 298, colour = "ivory") +
# Agulhas
geom_segment(aes(x = 33, y = -29.5, xend = 29.8, yend = -33.0),
arrow = arrow(length = unit(0.3, "cm")), size = 0.5, colour = "ivory") +
annotate("text", label = "Agulhas", x = 31.7, y = -31.7, size = 3.5, angle = 53, colour = "ivory") +
# Agulhas Bank
annotate("text", label = "Agulhas\nBank", x = 22.5, y = -35.5, size = 3.0, angle = 0, colour = "ivory") +
# Cape Peninsula
annotate("text", label = "Cape\nPeninsula", x = 17.2, y = -35, size = 3.0, angle = 0, colour = "ivory") +
# Improve on the x and y axis labels
scale_x_continuous(breaks = seq(15, 35, 5),
labels = scales::unit_format("°E", sep = ""),
position = "top") +
scale_y_continuous(breaks = seq(-35, -30, 5),
labels = c("35°S", "30°S")) +
labs(x = NULL, y = NULL) +
# Slightly shrink the plotting area
coord_cartesian(xlim = c(10.5, 39.5), ylim = c(-39.5, -25.5), expand = F) +
# Use viridis colour scheme
scale_fill_viridis(name = "Temp.\n(°C)", option = "D", breaks = c(24, 20, 16)) +
# Adjust the theme
theme_bw() +
theme(panel.border = element_rect(fill = NA, colour = "black", size = 1),
axis.text = element_text(size = 12, colour = "black"),
axis.ticks = element_line(colour = "black"))
# fig_1_top
# False Bay inset
fb_inset <- ggplot(data = sa_shore, aes(x = lon, y = lat)) +
# The land mass
geom_polygon(aes(group = PID),
fill = "grey70", colour = NA, size = 0.5, show.legend = FALSE) +
# The in situ sites
geom_point(data = SACTN_site_list, shape = 1, size = 3, colour = "black") +
geom_text(data = SACTN_site_list[-6,], aes(label = order), size = 2.0, colour = "red") +
# Text label
geom_text(aes(x = 18.65, y = -34.25, label = "False\nBay"), size = 2.7) +
# Control the x and y axes
coord_cartesian(xlim = c(18.2, 19), ylim = c(-34.5, -33.8), expand = F) +
scale_x_continuous(breaks = c(18.5), label = "18.5°E") +
scale_y_continuous(breaks = c(-34.1), label = "34.1°S") +
labs(x = NULL, y = NULL) +
# Change the theme for cleaner over-plotting
theme_bw() +
theme(plot.background = element_blank(),
axis.text = element_text(colour = "ivory"),
axis.text.y = element_text(angle = 90, hjust = 0.5),
axis.ticks = element_line(colour = "ivory"),
panel.border = element_rect(colour = "ivory"),
panel.grid = element_blank())
# fb_inset
# The bottom figure (air)
fig_1_bottom <- ggplot(data = southern_africa_coast, aes(x = lon, y = lat)) +
# The ocean temperature
geom_raster(data = air_temp, aes(fill = temp)) +
# The land mass
geom_polygon(aes(group = group), fill = NA, colour = "black", size = 0.5, show.legend = FALSE) +
geom_path(data = africa_borders, aes(group = group)) +
# The current vectors
geom_segment(data = winds, aes(xend = lon + u * wind_uv_scalar, yend = lat + v * wind_uv_scalar),
arrow = arrow(angle = 40, length = unit(winds$arrow_size, "cm"), type = "open"),
linejoin = "mitre", size = 0.4) +
# The legend for the vector length
geom_label(aes(x = 37.0, y = -38.0, label = "4.0 m/s\n"), size = 3, label.padding = unit(0.5, "lines")) +
geom_segment(aes(x = 36.0, y = -38.3, xend = 38.0, yend = -38.3), linejoin = "mitre",
arrow = arrow(angle = 40, length = unit(0.2, "cm"), type = "open")) +
# The coastal sections
geom_spoke(aes(x = 18.46520, y = -34.31050, angle = 180, radius = -2), linetype = "dotted", colour = "ivory") +
geom_spoke(aes(x = 18.46520, y = -34.31050, angle = 180, radius = 2), linetype = "dotted", colour = "ivory") +
geom_spoke(aes(x = 27.48889, y = -33.28611, angle = 40, radius = -2), linetype = "dotted", colour = "ivory") +
geom_spoke(aes(x = 27.48889, y = -33.28611, angle = 40, radius = 2), linetype = "dotted", colour = "ivory") +
annotate("text", label = "West\nCoast", x = 19.5, y = -31.2, size = 3.0, angle = 0, colour = "ivory") +
annotate("text", label = "South\nCoast", x = 23, y = -33.0, size = 3.0, angle = 0, colour = "ivory") +
annotate("text", label = "East\nCoast", x = 28, y = -31, size = 3.0, angle = 0, colour = "ivory") +
# South Atlantic Anticyclone
annotate("text", label = "SOUTH\nATLANTIC\nANTICYCLONE", x = 13.5, y = -33.5, size = 3.0, angle = 0, colour = "ivory") +
geom_curve(data = cyc_atlantic, aes(x = x, y = y, xend = xend, yend = yend), curvature = 0.2, colour = "ivory",
arrow = arrow(angle = 40, type = "open", length = unit(0.25,"cm"))) +
# South Indian Anticyclone
annotate("text", label = "SOUTH\nINDIAN\nANTICYCLONE", x = 36.5, y = -34.0, size = 3.0, angle = 0, colour = "ivory") +
geom_curve(data = cyc_indian, aes(x = x, y = y, xend = xend, yend = yend), curvature = 0.2, colour = "ivory",
arrow = arrow(angle = 40, type = "open", length = unit(0.25,"cm"))) +
# Westerlies
annotate("text", label = "WESTERLIES", x = 22.5, y = -37.0, size = 3.0, angle = 0, colour = "ivory") +
geom_curve(data = westerlies, aes(x = x, y = y, xend = xend, yend = yend), colour = "ivory",
arrow = arrow(angle = 40, type = "open", length = unit(0.25,"cm")), curvature = -0.01) +
# Improve on the x and y axis labels
scale_x_continuous(breaks = seq(15, 35, 5),
labels = scales::unit_format("°E", sep = "")) +
scale_y_continuous(breaks = seq(-35, -30, 5),
labels = c("35°S", "30°S")) +
labs(x = NULL, y = NULL) +
# Scale bar
scaleBar(lon = 22.0, lat = -29.5, distanceLon = 200, distanceLat = 50, distanceLegend = 90, dist.unit = "km",
arrow.length = 100, arrow.distance = 130, arrow.North.size = 3,
legend.colour = "ivory", arrow.colour = "ivory", N.colour = "ivory") +
# Slightly shrink the plotting area
coord_cartesian(xlim = c(10.5, 39.5), ylim = c(-39.5, -25.5), expand = F) +
# Use viridis colour scheme
scale_fill_viridis(name = "Temp.\n(°C)", option = "A", breaks = c(24, 20, 16)) +
# Adjust the theme
theme_bw() +
theme(panel.border = element_rect(fill = NA, colour = "black", size = 1),
axis.text = element_text(size = 12, colour = "black"),
axis.ticks = element_line(colour = "black"))
# fig_1_bottom
# Convert the figures to grobs
fig_1_top_grob <- ggplotGrob(fig_1_top)
fb_inset_grob <- ggplotGrob(fb_inset)
fig_1_bottom_grob <- ggplotGrob(fig_1_bottom)
# Stick them together
fig_1 <- ggplot() +
# First set the x and y axis values so we know what the ranges are
# in order to make it easier to place our facets
coord_equal(xlim = c(1, 10), ylim = c(1, 10), expand = F) +
# Then we place our facetsover one another using the coordinates we created
annotation_custom(fig_1_top_grob,
xmin = 1, xmax = 10, ymin = 5.5, ymax = 10) +
annotation_custom(fb_inset_grob,
xmin = 3.5, xmax = 5.5, ymin = 7.2, ymax = 8.8) +
annotation_custom(fig_1_bottom_grob,
xmin = 1, xmax = 10, ymin = 1, ymax = 5.5)
# save
ggsave(plot = fig_1, filename = "graph/fig_1.pdf", height = 8, width = 8)
Oliver, E. C., Lago, V., Hobday, A. J., Holbrook, N. J., Ling, S. D., and Mundy, C. N. (2018). Marine heatwaves off eastern tasmania: Trends, interannual variability, and predictability. Progress in oceanography 161, 116–130.
Schlegel, R. W., Oliver, E. C., Perkins-Kirkpatrick, S., Kruger, A., and Smit, A. J. (2017). Predominant atmospheric and oceanic patterns during coastal marine heatwaves. Frontiers in Marine Science 4, 323.
sessionInfo()
R version 3.6.0 (2019-04-26)
Platform: x86_64-pc-linux-gnu (64-bit)
Running under: Ubuntu 16.04.5 LTS
Matrix products: default
BLAS: /usr/lib/openblas-base/libblas.so.3
LAPACK: /usr/lib/libopenblasp-r0.2.18.so
locale:
[1] LC_CTYPE=en_CA.UTF-8 LC_NUMERIC=C
[3] LC_TIME=en_CA.UTF-8 LC_COLLATE=en_CA.UTF-8
[5] LC_MONETARY=en_CA.UTF-8 LC_MESSAGES=en_CA.UTF-8
[7] LC_PAPER=en_CA.UTF-8 LC_NAME=C
[9] LC_ADDRESS=C LC_TELEPHONE=C
[11] LC_MEASUREMENT=en_CA.UTF-8 LC_IDENTIFICATION=C
attached base packages:
[1] stats graphics grDevices utils datasets methods base
other attached packages:
[1] bindrcpp_0.2.2 data.table_1.11.6 yasomi_0.3
[4] proxy_0.4-22 e1071_1.7-0 lubridate_1.7.4
[7] ncdf4_1.16 forcats_0.3.0 stringr_1.3.1
[10] dplyr_0.7.6 purrr_0.2.5 readr_1.1.1
[13] tidyr_0.8.1 tibble_1.4.2 ggplot2_3.0.0
[16] tidyverse_1.2.1 jsonlite_1.6
loaded via a namespace (and not attached):
[1] tidyselect_0.2.4 haven_1.1.2 lattice_0.20-35
[4] colorspace_1.3-2 htmltools_0.3.6 yaml_2.2.0
[7] rlang_0.2.2 R.oo_1.22.0 pillar_1.3.0
[10] glue_1.3.0 withr_2.1.2 R.utils_2.7.0
[13] doMC_1.3.5 modelr_0.1.2 readxl_1.1.0
[16] foreach_1.4.4 bindr_0.1.1 plyr_1.8.4
[19] munsell_0.5.0 gtable_0.2.0 workflowr_1.1.1
[22] cellranger_1.1.0 rvest_0.3.2 R.methodsS3_1.7.1
[25] codetools_0.2-15 evaluate_0.11 knitr_1.20
[28] parallel_3.6.0 class_7.3-14 broom_0.5.0
[31] Rcpp_0.12.18 backports_1.1.2 scales_1.0.0
[34] hms_0.4.2 digest_0.6.16 stringi_1.2.4
[37] grid_3.6.0 rprojroot_1.3-2 cli_1.0.0
[40] tools_3.6.0 maps_3.3.0 magrittr_1.5
[43] lazyeval_0.2.1 crayon_1.3.4 whisker_0.3-2
[46] pkgconfig_2.0.2 xml2_1.2.0 iterators_1.0.10
[49] assertthat_0.2.0 rmarkdown_1.10 httr_1.3.1
[52] rstudioapi_0.7 R6_2.2.2 nlme_3.1-137
[55] git2r_0.23.0 compiler_3.6.0
This reproducible R Markdown analysis was created with workflowr 1.1.1