## ----setup, include = TRUE---------------------------------------------------- knitr::opts_chunk$set(collapse = TRUE, comment = "#>") ext_file <- function(...) { path <- system.file("extdata", ..., package = "gbif.range") if (nzchar(path)) { return(path) } normalizePath(file.path("..", "inst", "extdata", ...), mustWork = TRUE) } library(gbif.range) ## ----data-retrieval, eval = FALSE--------------------------------------------- # # Study extent, world boundaries, and packaged terrestrial ecoregions # shp.lonlat <- terra::vect(ext_file("shp_lonlat.shp")) # countries <- terra::vect(ext_file("world_countries.shp")) # eco.terra <- read_ecoreg( # ecoreg_name = "eco_terra", save_dir = tempdir() # ) # # # Custom, high-resolution ecoregions for the Alps study extent # # ('rst' --> 5x5-km resolution) # rst <- terra::rast(ext_file("rst.tif")) # my.eco <- make_ecoreg(env = rst, nclass = 200, format = "SpatVector") # # # Continent extent to keep only terrestrial my.eco # contExtL.shp <- terra::aggregate( # terra::crop(countries, terra::ext(rst)) # ) # contExtS.shp <- terra::crop( # contExtL.shp, # terra::ext(shp.lonlat) # ) # my.ecoS <- intersect(my.eco, contExtS.shp) # # # GBIF occurrences for the two focal species used throughout # obs.arcto <- get_gbif( # sp_name = "Arctostaphylos alpinus", # geo = shp.lonlat, # grain = 1 # ) # obs.pt <- get_gbif(sp_name = "Panthera tigris") # # # Diversity rasters used in Figure 3 # iucn.robin <- terra::rast(ext_file("iucn_div_robin.tif")) # gf.robin <- terra::rast(ext_file("gf_div_robin.tif")) ## ----fig1a-left, eval = FALSE------------------------------------------------- # # Regional range at ~5x5-km resolution ('my.eco' resolution) # range.arcto <- get_range( # occ_coord = obs.arcto, # ecoreg = my.eco, # ecoreg_name = "EcoRegion", # format = "SpatVector" # ) # # # Assign colors to ecoregions # range.arctoS <- terra::crop(range.arcto$rangeOutput, contExtS.shp) # col.palette <- grDevices::colorRampPalette(c( # "#a6cee3", "#1f78b4", "#b2df8a", "#33a02c", # "#fb9a99", "#e31a1c", "#fdbf6f", "#ff7f00", # "#cab2d6", "#6a3d9a", "#ffff99", "#b15928" # )) # colcol <- col.palette(length(my.ecoS)) # set.seed(7) # my.ecoS$color <- sample( # paste0(colcol, ""), # length(my.ecoS), # replace = FALSE # ) # # pt.col <- terra::extract( # x = my.ecoS, # y = as.data.frame( # obs.arcto[, c( # "decimalLongitude", "decimalLatitude" # )] # ) # ) # pt.plot <- obs.arcto[ # !is.na(pt.col$color), # c("decimalLongitude", "decimalLatitude") # ] # pt.col2 <- pt.col[!is.na(pt.col$color), "color"] # pt.col3 <- grDevices::adjustcolor( # pt.col2, red.f = 0.6, green.f = 0.6, # blue.f = 0.6 # ) # # # Plot # terra::plot( # my.ecoS, # col = paste0(my.ecoS$color, "99"), # border = NA, # axes = FALSE # ) # terra::plot( # merge_range(range.arctoS), # border = "black", # lwd = 1, # col = "#00000099", # add = TRUE # ) # graphics::points( # pt.plot, col = pt.col2, pch = 16, cex = 1 # ) # graphics::points( # pt.plot, col = pt.col3, pch = 16, cex = 0.6 # ) ## ----fig-1a-left, echo = FALSE, out.width = "50%"----------------------------- knitr::include_graphics("../man/figures/Fig1a_left.png") ## ----fig1a-right, eval = FALSE------------------------------------------------ # # Global range # range.tiger <- get_range( # occ_coord = obs.pt, # ecoreg = eco.terra, # ecoreg_name = "ECO_NAME", # format = "SpatVector" # ) # # # Assign colors to ecoregions # ext.tiger.eco <- terra::ext(range.tiger$rangeOutput) # ext.tiger.eco <- c( # ext.tiger.eco[1] - 2, ext.tiger.eco[2] + 2, # ext.tiger.eco[3] - 2, ext.tiger.eco[4] + 2 # ) # eco.local <- terra::crop(eco.terra, ext.tiger.eco) # col.palette <- grDevices::colorRampPalette(c( # "#a6cee3", "#1f78b4", "#b2df8a", "#33a02c", # "#fb9a99", "#e31a1c", "#fdbf6f", "#ff7f00", # "#cab2d6", "#6a3d9a", "#ffff99", "#b15928" # )) # colcol <- col.palette(length(eco.local)) # set.seed(3) # eco.local$color <- sample( # paste0(colcol, ""), # length(eco.local), # replace = FALSE # ) # pt.coords <- as.data.frame( # obs.pt[, c("decimalLongitude", "decimalLatitude")] # ) # pt.col <- terra::extract(eco.local, pt.coords) # pt.plot <- obs.pt[ # !is.na(pt.col$color), # c("decimalLongitude", "decimalLatitude") # ] # pt.col2 <- pt.col[!is.na(pt.col$color), "color"] # pt.col3 <- grDevices::adjustcolor( # pt.col2, red.f = 0.6, green.f = 0.6, # blue.f = 0.6 # ) # out.plot <- terra::extract(range.tiger$rangeOutput, pt.coords) # op.na <- is.na(out.plot[, 2]) # out.plot <- obs.pt[ # op.na, c("decimalLongitude", "decimalLatitude") # ] # # # Plot # terra::plot( # eco.local, # col = eco.local$color, # border = NA, # axes = FALSE # ) # terra::plot( # merge_range(range.tiger), # border = "black", # lwd = 2, # col = "#63636399", # add = TRUE # ) # graphics::points( # pt.plot, col = pt.col2, pch = 16, cex = 1.5 # ) # graphics::points( # pt.plot, col = pt.col3, pch = 16, cex = 0.8 # ) # graphics::points( # out.plot, col = "black", pch = 4, cex = 1.5, # lwd = 2 # ) ## ----fig-1a-right, echo = FALSE, out.width = "50%"---------------------------- knitr::include_graphics("../man/figures/Fig1a_right.png") ## ----fig2a, eval = FALSE------------------------------------------------------ # # ------------------------------------------ # # Arctostaphylos alpinus (custom ecoregion) # # ------------------------------------------ # # # Create pseudo-absences # # First remove observations considered outliers in get_range() # xy.df <- range.arcto$init.args$occ_coord # r.ext <- terra::ext(range.arcto$rangeOutput) # Xrm.cond <- xy.df$decimalLongitude >= r.ext[1] & # xy.df$decimalLongitude <= r.ext[2] # Yrm.cond <- xy.df$decimalLatitude >= r.ext[3] & # xy.df$decimalLatitude <= r.ext[4] # xy.df <- xy.df[Xrm.cond & Yrm.cond, ] # # Sample n regular background points over the range extent # x.interv <- (r.ext[2] - r.ext[1]) / (sqrt(1e4) - 1) # y.interv <- (r.ext[4] - r.ext[3]) / (sqrt(1e4) - 1) # lx <- seq(r.ext[1], r.ext[2], x.interv) # ly <- seq(r.ext[3], r.ext[4], y.interv) # bp.xy <- expand.grid( # decimalLongitude = lx, # decimalLatitude = ly # ) # # # Combine observations with background # obs.xy <- xy.df[, c("decimalLongitude", "decimalLatitude")] # all.xy <- rbind(obs.xy, bp.xy) # all.xy$Pres <- 0 # all.xy[1:nrow(obs.xy), "Pres"] <- 1 # # # Run block-cv # xy.pres <- all.xy$Pres # cv.strat <- make_blocks( # nfolds = 5, # df = all.xy[, c("decimalLongitude", "decimalLatitude")], # nblocks = 5 * 2, # pres = xy.pres # ) # all.xy$bcv <- cv.strat # all.xy[all.xy$bcv %in% 1, "col"] <- "#e41a1c" # all.xy[all.xy$bcv %in% 2, "col"] <- "#377eb8" # all.xy[all.xy$bcv %in% 3, "col"] <- "#4daf4a" # all.xy[all.xy$bcv %in% 4, "col"] <- "#984ea3" # all.xy[all.xy$bcv %in% 5, "col"] <- "#ff7f00" # all.xy[all.xy$Pres %in% 1, "col"] <- substr( # grDevices::adjustcolor( # all.xy[all.xy$Pres %in% 1, "col"], # red.f = 0.5, green.f = 0.5, blue.f = 0.5 # ), # 1, 7 # ) # # # Evaluate # ar.test <- cv_range( # range_object = range.arcto, # cv = "block-cv", # nfolds = 5, # nblocks = 2 # ) # # # Use the broader (L) extent and country borders so # # the evaluation text fits on the plot # ext.L <- terra::ext(contExtL.shp) # world.local.ar <- terra::crop(countries, ext.L) # world.local.ar <- terra::aggregate(world.local.ar) # pres <- all.xy[all.xy$Pres %in% 1, ] # abs.pts <- all.xy[all.xy$Pres %in% 0, ] # pres_coords <- as.data.frame( # pres[, c("decimalLongitude", "decimalLatitude")] # ) # id.in <- terra::extract( # range.arcto$rangeOutput, pres_coords # ) # pres <- pres[!is.na(id.in[, 2]), ] # # # Plot # terra::plot( # world.local.ar, col = "#bcd1bc", # axes = FALSE, lwd = 1 # ) # graphics::points( # abs.pts, col = paste0(abs.pts$col, "50"), # pch = 16, cex = 0.5 # ) # terra::plot( # merge_range(range.arctoS), # border = "black", lwd = 1.7, # col = "#63636370", add = TRUE # ) # graphics::points( # pres, col = paste0(pres$col, "90"), # pch = 16, cex = 1.3 # ) # # # Evaluation text anchored to the L extent's bottom-left corner # txt.x <- ext.L[1] + 0.02 * (ext.L[2] - ext.L[1]) # txt.y1 <- ext.L[3] + 0.09 * (ext.L[4] - ext.L[3]) # txt.y2 <- ext.L[3] + 0.04 * (ext.L[4] - ext.L[3]) # # graphics::text( # txt.x, txt.y1, # paste( # "Mean TSS =", # round(tail(ar.test[, "TSS"], 1), 2) # ), # cex = 1.2, font = 2, adj = 0 # ) # graphics::text( # txt.x, txt.y2, # paste( # "Mean Precision =", # round(tail(ar.test[, "Precision"], 1), 2) # ), # cex = 1.2, font = 2, adj = 0 # ) ## ----fig-2a, echo = FALSE, out.width = "50%"---------------------------------- knitr::include_graphics("../man/figures/Fig2a_arcto_cv.png") ## ----fig2b, eval = FALSE------------------------------------------------------ # # ------------------------------------------ # # Panthera tigris (packaged ecoregion) # # ------------------------------------------ # # ext.temp.t2 <- terra::ext(range.tiger$rangeOutput) # ext.temp.t2 <- c( # ext.temp.t2[1] - 0.2, ext.temp.t2[2] + 0.2, # ext.temp.t2[3] - 0.2, ext.temp.t2[4] + 0.2 # ) # # # Create pseudo-absences # # First remove observations considered outliers in get_range() # xy.df <- range.tiger$init.args$occ_coord # r.ext <- terra::ext(range.tiger$rangeOutput) # Xrm.cond <- xy.df$decimalLongitude >= r.ext[1] & # xy.df$decimalLongitude <= r.ext[2] # Yrm.cond <- xy.df$decimalLatitude >= r.ext[3] & # xy.df$decimalLatitude <= r.ext[4] # xy.df <- xy.df[Xrm.cond & Yrm.cond, ] # # Sample n regular background points over the range extent # x.interv <- (r.ext[2] - r.ext[1]) / (sqrt(1e4) - 1) # y.interv <- (r.ext[4] - r.ext[3]) / (sqrt(1e4) - 1) # lx <- seq(r.ext[1], r.ext[2], x.interv) # ly <- seq(r.ext[3], r.ext[4], y.interv) # bp.xy <- expand.grid( # decimalLongitude = lx, # decimalLatitude = ly # ) # # # Combine observations with background # obs.xy <- xy.df[, c("decimalLongitude", "decimalLatitude")] # all.xy <- rbind(obs.xy, bp.xy) # all.xy$Pres <- 0 # all.xy[1:nrow(obs.xy), "Pres"] <- 1 # # # Run block-cv # xy.pres <- all.xy$Pres # cv.strat <- make_blocks( # nfolds = 5, # df = all.xy[, c("decimalLongitude", "decimalLatitude")], # nblocks = 5 * 2, # pres = xy.pres # ) # all.xy$bcv <- cv.strat # all.xy[all.xy$bcv %in% 1, "col"] <- "#e41a1c" # all.xy[all.xy$bcv %in% 2, "col"] <- "#377eb8" # all.xy[all.xy$bcv %in% 3, "col"] <- "#4daf4a" # all.xy[all.xy$bcv %in% 4, "col"] <- "#984ea3" # all.xy[all.xy$bcv %in% 5, "col"] <- "#ff7f00" # all.xy[all.xy$Pres %in% 1, "col"] <- substr( # grDevices::adjustcolor( # all.xy[all.xy$Pres %in% 1, "col"], # red.f = 0.5, green.f = 0.5, blue.f = 0.5 # ), # 1, 7 # ) # # # Evaluate # pt.test <- cv_range( # range_object = range.tiger, # cv = "block-cv", # nfolds = 5, # nblocks = 2 # ) # # # Use extent and country borders # world.local.ti <- terra::crop(countries, ext.temp.t2) # world.local.ti <- terra::aggregate(world.local.ti) # pres <- all.xy[all.xy$Pres %in% 1, ] # abs.pts <- all.xy[all.xy$Pres %in% 0, ] # pres_coords <- as.data.frame( # pres[, c("decimalLongitude", "decimalLatitude")] # ) # id.in <- terra::extract(range.tiger$rangeOutput, pres_coords) # pres <- pres[!is.na(id.in[, 2]), ] # # # Plot # terra::plot( # world.local.ti, col = "#bcd1bc", # axes = FALSE, lwd = 1 # ) # graphics::points( # abs.pts, col = paste0(abs.pts$col, "50"), # pch = 16, cex = 0.6 # ) # terra::plot( # merge_range(range.tiger), # border = "black", lwd = 2, # col = "#63636370", add = TRUE # ) # graphics::points( # pres, col = paste0(pres$col, "80"), # pch = 16, cex = 1.6 # ) # # # Evaluation text anchored to the tiger extent's top-left corner # txt.x.t <- ext.temp.t2[1] + 0.02 * (ext.temp.t2[2] - ext.temp.t2[1]) # txt.y1.t <- ext.temp.t2[4] - 0.05 * (ext.temp.t2[4] - ext.temp.t2[3]) # txt.y2.t <- ext.temp.t2[4] - 0.10 * (ext.temp.t2[4] - ext.temp.t2[3]) # # graphics::text( # txt.x.t, txt.y1.t, # paste( # "Mean TSS =", # round(tail(pt.test[, "TSS"], 1), 2) # ), # cex = 1.5, font = 2, adj = 0 # ) # graphics::text( # txt.x.t, txt.y2.t, # paste( # "Mean Precision =", # round(tail(pt.test[, "Precision"], 1), 2) # ), # cex = 1.5, font = 2, adj = 0 # ) ## ----fig-2b, echo = FALSE, out.width = "50%"---------------------------------- knitr::include_graphics("../man/figures/Fig2b_tiger_cv.png") ## ----fig2c, eval = FALSE------------------------------------------------------ # # ------------------------------------- # # Arctostaphylos alpinus (custom ecoregion) # # ------------------------------------- # # root.dir <- ext_file() # # res5km <- evaluate_range( # root_dir = root.dir, # valData_dir = "SDM", # ecoRM_dir = "EcoRM", # verbose = TRUE, # print_map = FALSE, # valData_type = "TIFF", # mask = NULL, # res_fact = NULL # ) # # # Plot plant # terra::plot( # contExtL.shp, col = "#dce8dc", # axes = FALSE, lwd = 1 # ) # terra::plot( # merge_range(range.arctoS), # lwd = 0.1, col = "#63636350", add = TRUE # ) # terra::plot( # res5km$overlay_list[[1]], # col = c( # "#d1c845", "#e39c59", "#8d60cc", "#6cba6c" # ), # breaks = c(-0.5, 0.5, 1.5, 2.5, 3.5), # axes = FALSE, legend = FALSE, # las = 1, add = TRUE # ) # # # Evaluation text anchored to the L extent's bottom-left corner # ext.L <- terra::ext(contExtL.shp) # txt.x <- ext.L[1] + 0.02 * (ext.L[2] - ext.L[1]) # txt.y1 <- ext.L[3] + 0.09 * (ext.L[4] - ext.L[3]) # txt.y2 <- ext.L[3] + 0.04 * (ext.L[4] - ext.L[3]) # # graphics::text( # txt.x, txt.y1, # paste( # "Mean TSS =", # round(res5km$df_eval[1, "TSS_ecoRM"], 2) # ), # cex = 1, font = 2, adj = 0 # ) # graphics::text( # txt.x, txt.y2, # paste( # "Mean Precision =", # round(res5km$df_eval[1, "Prec_ecoRM"], 2) # ), # cex = 1, font = 2, adj = 0 # ) # terra::plot( # contExtL.shp, # border = "#383d38", axes = FALSE, # lwd = 1, add = TRUE # ) ## ----fig-2c, echo = FALSE, out.width = "50%"---------------------------------- knitr::include_graphics("../man/figures/Fig2c_arcto_eval.png") ## ----fig2d, eval = FALSE------------------------------------------------------ # # ------------------------------------- # # Panthera tigris (packaged ecoregion) # # ------------------------------------- # # # Plot tiger # terra::plot( # world.local.ti, col = "#dce8dc", # axes = FALSE, lwd = 1 # ) # toPlot <- terra::mask( # res5km$overlay_list[[6]], world.local.ti # ) # terra::plot( # toPlot, # col = c( # "#d1c845", "#e39c59", "#8d60cc", "#6cba6c" # ), # breaks = c(-0.5, 0.5, 1.5, 2.5, 3.5), # axes = FALSE, legend = FALSE, # las = 1, add = TRUE # ) # # # Evaluation text anchored to the tiger extent's top-left corner # txt.x.t <- ext.temp.t2[1] + 0.02 * (ext.temp.t2[2] - ext.temp.t2[1]) # txt.y1.t <- ext.temp.t2[4] - 0.05 * (ext.temp.t2[4] - ext.temp.t2[3]) # txt.y2.t <- ext.temp.t2[4] - 0.10 * (ext.temp.t2[4] - ext.temp.t2[3]) # # graphics::text( # txt.x.t, txt.y1.t, # paste( # "Mean TSS =", # round(res5km$df_eval[6, "TSS_ecoRM"], 2) # ), # cex = 1.5, font = 2, adj = 0 # ) # graphics::text( # txt.x.t, txt.y2.t, # paste( # "Mean Precision =", # round(res5km$df_eval[6, "Prec_ecoRM"], 2) # ), # cex = 1.5, font = 2, adj = 0 # ) # terra::plot( # world.local.ti, # border = "#383d38", axes = FALSE, # lwd = 1, add = TRUE # ) # leg.x <- ext.temp.t2[1] + 0.05 * (ext.temp.t2[2] - ext.temp.t2[1]) # leg.y <- ext.temp.t2[3] + 0.20 * (ext.temp.t2[4] - ext.temp.t2[3]) # graphics::legend( # leg.x, leg.y, # legend = c( # "True Presences", "True Absences", # "False Presences", "False Absences" # ), # fill = c( # "#6cba6c", "#d1c845", "#e39c59", "#8d60cc" # ), # bg = NA, box.col = NA, # cex = 1.1, # x.intersp = 0.2 # ) ## ----fig-2d, echo = FALSE, out.width = "50%"---------------------------------- knitr::include_graphics("../man/figures/Fig2d_tiger_eval.png") ## ----fig3a, eval = FALSE------------------------------------------------------ # # ------------------------------------- # # World maps: IUCN vs gbif.range richness # # ------------------------------------- # # # CRS # robin <- paste( # "+proj=robin +lon_0=0 +x_0=0 +y_0=0", # "+datum=WGS84 +units=m +no_defs +type=crs" # ) # # # Reproject countries # countries.robin <- terra::project(countries, robin) # # # Boundary box # bb <- terra::as.polygons( # terra::ext(-180, 180, -90, 90), # crs = "EPSG:4326" # ) # bb <- terra::densify(bb, interval = 1) # bb.robin <- terra::project(bb, robin) # # # Shared diverging color palette for both maps # col <- grDevices::colorRampPalette( # c("#67a9cf", "#f7f7f7", "#ef8a62") # ) # # # Plot IUCN map # max.iucn <- round(terra::minmax(iucn.robin)[2]) # terra::plot( # iucn.robin, # axes = FALSE, legend = FALSE, # col = col(10), smooth = TRUE, # mar = c(1, 1, 1, 5) # ) # terra::plot(countries.robin, add = TRUE, lwd = 2) # terra::plot(bb.robin, add = TRUE, lwd = 3) # oldpar <- graphics::par(xpd = NA, lwd = 3) # cscl( # colors = col(10), # crds = c(-9501111, 9734033, -11800000, -13000000), # zrng = c(0, max.iucn), # tickle = -0.3, cx = 1.1, lablag = -1.3, # tria = "b", at = seq(0, max.iucn, 10), # horiz = TRUE, title = "IUCN richness", # labs = seq(0, max.iucn, 10), titlag = 2 # ) # graphics::par(oldpar) # # # Plot GBIF.RANGE map # max.gf <- round(terra::minmax(gf.robin)[2]) # terra::plot( # gf.robin, # axes = FALSE, legend = FALSE, # col = col(10), smooth = TRUE, # mar = c(1, 1, 1, 5) # ) # terra::plot(countries.robin, add = TRUE, lwd = 2) # terra::plot(bb.robin, add = TRUE, lwd = 3) # oldpar <- graphics::par(xpd = NA, lwd = 3) # cscl( # colors = col(10), # crds = c(-9501111, 9734033, -11800000, -13000000), # zrng = c(0, max.gf), # tickle = -0.3, cx = 1.1, lablag = -1.3, # tria = "b", at = seq(0, max.gf, 10), # horiz = TRUE, title = "gbif.range richness", # labs = seq(0, max.gf, 10), titlag = 2 # ) # graphics::par(oldpar) ## ----fig-3-maps, echo = FALSE, out.width = "50%", fig.show = "hold"----------- knitr::include_graphics(c( "../man/figures/Fig3_iucn.png", "../man/figures/Fig3_gbifrange.png" )) ## ----fig3b, eval = FALSE------------------------------------------------------ # # ------------------------------------- # # Scatter plots: richness and area agreement # # ------------------------------------- # # # Load area table # data(area_data) # # # Extract data to plot # cor.ras <- terra::rast(list(gf.robin, iucn.robin)) # names(cor.ras) <- c("RANGE", "IUCN") # set.seed(1) # samp.div <- terra::spatSample( # cor.ras, 5000, replace = FALSE, na.rm = TRUE # ) # dat.plot <- list(samp.div, area_data[, -1]) # strings <- c("richness", "areas") # # # Plot diversities relationship, side by side # oldpar <- graphics::par(mfrow = c(1, 2), mar = c(5, 5.5, 4, 2)) # lapply(seq_along(dat.plot), function(x) { # # # Extract the data # xy <- dat.plot[[x]] # sex <- 2.5 # col <- "#d95f0240" # add <- "" # # # Log or not depending on the output # if (strings[x] == "areas") { # xy <- log(xy + 1) # sex <- 2.5 # col <- "#d95f0240" # add <- "(log)" # } # # # Plot points # graphics::plot( # xy[, 2], xy[, 1], # cex.axis = 1.7, cex = sex, col = col, # xlim = c(min(xy[, 2]), max(xy[, 2])), # ylim = c(min(xy[, 1]), max(xy[, 1])), # pch = 16, # xlab = paste("IUCN", strings[x], add), # ylab = paste("gbif.range", strings[x], add), # cex.lab = 1.7, font.lab = 1 # ) # # # Run Linear Models and spearman's correlation # lm.div <- stats::lm(xy[, 1] ~ xy[, 2]) # # # Lin's concordance correlation coefficient # # (population covariance/variance) # ccc_lin <- function(x, y) { # mx <- mean(x) # my <- mean(y) # vx <- mean((x - mx)^2) # vy <- mean((y - my)^2) # sxy <- mean((x - mx) * (y - my)) # 2 * sxy / (vx + vy + (mx - my)^2) # } # ccc.div <- ccc_lin(xy[, 2], xy[, 1]) # cor.div <- stats::cor(xy[, 2], xy[, 1]) # adj.r2 <- summary(lm.div)[[9]] # # # Plot text # text_cor1 <- bquote("ccc" == .(round(ccc.div, 2))) # text_cor2 <- bquote("r" == .(round(cor.div, 2))) # fig_label( # text_cor1, # region = "plot", pos = "topleft", # bty = "n", font = 2, col = "#121212", # cex = 2, margin = 0.02 # ) # fig_label( # text_cor2, # region = "plot", pos = "bottomright", # bty = "n", font = 2, col = "#6f69c2", # cex = 2, margin = 0.02 # ) # # # Plot relationship # graphics::lines( # xy[, 2], lm.div$fit, lwd = 7, col = "#7570b3" # ) # graphics::abline( # a = 0, b = 1, col = "#252525", # lwd = 5, lty = 2 # ) # }) # graphics::par(oldpar) ## ----fig-3-scatter, echo = FALSE, out.width = "80%", fig.show = "hold"-------- knitr::include_graphics(c( "../man/figures/Fig3_scatterplots.png" ))