#' Assessing whether phenotypic changes in response to climate are adaptive
#'
#' This R package aims at providing the data and documenting the analysis behind
#' the results from the paper entitled 'Adaptive responses of animals to climate
#' change are most likely insufficient' by Radchuk et al. Nature Communications
#' (2019).
#'
#' This package has not been conceived for general use!
#'
#' All main functions of this package contain a short documentation and examples
#' which could be useful for those who try to understand our code. Type
#' \code{ls("package:adRes")} for a list of all exported functions.
#'
#' We recommend you to follow the examples below to reproduce the results of our
#' analysis and understand the structure of our workflow.
#'
#' You may also directly explore the files contained in the package after
#' uncompressing the content of the *.tar.gz file (link available on the GitHub
#' page \url{https://github.com/radchukv/adRes}). You can use
#' the R function \code{\link{untar}} to extract the content of the tarball.
#'
#' The package also contains several fitted models directly, to spare you the
#' time required to fit them.
#'
#' In the examples below, we provide the workflow leading the results presented
#' in the paper.
#'
#' @name adRes-package
#' @aliases adRes-package adRes
#' @docType package
#'
#' @importFrom dplyr %>%
#'
#' @references
#' Radchuk et al. (2019) Adaptive responses of animals to climate change
#' are most likely insufficient". Nature Communications.
#'
#' @keywords package
#'
#' @examples
#' #~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~#
#' #### Settings for the workflow ####
#' #~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~#
#' nb_cores <- 2L ## increase the number for using more cores
#' digit <- 3L
#' par_ini <- par(no.readonly = TRUE) ## default parameters for display
#'
#'
#'
#' #~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~#
#' #### Temperature: Condition 1 ####
#' #~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~#
#' ## for PRCS
#' mod_T_prcs <- fit_all(data = dat_Clim, temperature = TRUE,
#' precipitation = FALSE, phenology = TRUE,
#' morphology = TRUE, condition = '1',
#' nb_cores = nb_cores, rand_trait = FALSE,
#' fixed = NULL, digit = digit)
#'
#' ## for PRC
#' \dontrun{
#' mod_T_prc <- fit_all(data = dat_Clim_prc, temperature = TRUE,
#' precipitation = FALSE, phenology = TRUE,
#' morphology = TRUE, condition = '1',
#' nb_cores = nb_cores, rand_trait = FALSE,
#' fixed = NULL, digit = digit)
#'
#' ## forest plot as in Fig. 2
#' par(par_ini)
#' plot_T_cond1 <- plot_forest(meta_obj1 = mod_T_prcs$meta_res,
#' meta_obj2 = NULL,
#' list_extra_meta_obj =
#' list(mod_T_prcs$meta_res, mod_T_prc$meta_res),
#' sort = c("Coord"),
#' increasing = FALSE,
#' labels = c(traits = FALSE,
#' fitness = FALSE,
#' country = TRUE,
#' authors = FALSE))
#' mtext('Across', side = 2, line = 3,
#' at = -1, las = 2, cex = 0.9, col = 'black')
#' mtext('studies', side = 2, line = 3,
#' at = -2, las = 2, cex = 0.9, col = 'black')
#' }
#'
#' ## plots of raw data and extracted slopes per study, as in Suppl. Fig. S11
#' par(par_ini); par(oma = c(3, 5, 0, 0), mar = c(1, 2, 3, 1))
#' plot_raw(data = dat_Clim, temperature = TRUE,
#' precipitation = FALSE, phenology = TRUE, morphology = TRUE,
#' condition = '1', id_to_do = c(1:20))
#'
#' ## in the full dataset there are 41 ids (21_41). In the publicly shared 36
#' par(par_ini); par(oma = c(3, 5, 0, 0), mar = c(1, 2, 3, 1))
#' plot_raw(data = dat_Clim, temperature = TRUE,
#' precipitation = FALSE, phenology = TRUE,
#' morphology = TRUE, condition = '1',
#' id_to_do = c(21:36))
#'
#'
#'
#' #~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~#
#' #### Temperature: Condition 2, phenology ####
#' #~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~#
#' ## for PRCS
#' mod_phen_T_prcs <- fit_all(data = dat_Trait, temperature = TRUE,
#' precipitation = FALSE, phenology = TRUE,
#' morphology = FALSE, condition = '2',
#' nb_cores = nb_cores, rand_trait = FALSE,
#' fixed = NULL, digit = digit)
#'
#' ## fixed effects
#' mod_phen_T_prcs_Type <- fit_all(data = dat_Trait,
#' temperature = TRUE, precipitation = FALSE,
#' phenology = TRUE, morphology = FALSE,
#' condition = '2', nb_cores = nb_cores,
#' rand_trait = FALSE, fixed = 'Trait_Cat',
#' digit = digit)
#'
#' ## for PRC
#' \dontrun{
#' mod_phen_T_prc <- fit_all(data = dat_Trait_prc,
#' temperature = TRUE, precipitation = FALSE,
#' phenology = TRUE, morphology = FALSE,
#' condition = '2', nb_cores = nb_cores,
#' rand_trait = FALSE, fixed = NULL,
#' digit = digit)
#'
#' ## fixed effects
#' mod_phen_T_prc_Type <- fit_all(data = dat_Trait_prc,
#' temperature = TRUE, precipitation = FALSE,
#' phenology = TRUE, morphology = FALSE,
#' condition = '2', nb_cores = nb_cores,
#' rand_trait = FALSE, fixed = 'Trait_Cat',
#' digit = digit)
#'
#' mod_phen_T_prc_Taxon <- fit_all(data = dat_Trait_prc,
#' temperature = TRUE, precipitation = FALSE,
#' phenology = TRUE, morphology = FALSE,
#' condition = '2', nb_cores = nb_cores,
#' rand_trait = FALSE, fixed = 'Taxon',
#' digit = digit)
#' }
#'
#'
#'
#' #~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~#
#' #### Temperature: Condition 2, morphology ####
#' #~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~#
#' ## for PRCS
#' mod_morph_T_prcs <- fit_all(data = dat_Trait, temperature = TRUE,
#' precipitation = FALSE, phenology = FALSE,
#' morphology = TRUE, condition = '2',
#' nb_cores = nb_cores, rand_trait = FALSE,
#' fixed = NULL, digit = digit)
#'
#' ## fixed effects
#' mod_morph_T_prcs_Type <- fit_all(data = dat_Trait,
#' temperature = TRUE, precipitation = FALSE,
#' phenology = FALSE, morphology = TRUE,
#' condition = '2', nb_cores = nb_cores,
#' rand_trait = FALSE, fixed = 'Morph_type',
#' digit = digit)
#'
#' ## for PRC
#' \dontrun{
#' mod_morph_T_prc <- fit_all(data = dat_Trait_prc,
#' temperature = TRUE, precipitation = FALSE,
#' phenology = FALSE, morphology = TRUE,
#' condition = '2', nb_cores = nb_cores,
#' rand_trait = FALSE, fixed = NULL,
#' digit = digit)
#'
#' ## fixed effects
#' mod_morph_T_prc_Type <- fit_all(data = dat_Trait_prc,
#' temperature = TRUE, precipitation = FALSE,
#' phenology = FALSE, morphology = TRUE,
#' condition = '2', nb_cores = nb_cores,
#' rand_trait = FALSE, fixed = 'Morph_type',
#' digit = digit)
#'
#' mod_morph_T_prc_Taxon <- fit_all(data = dat_Trait_prc,
#' temperature = TRUE, precipitation = FALSE,
#' phenology = FALSE, morphology = TRUE,
#' condition = '2', nb_cores = nb_cores,
#' rand_trait = FALSE, fixed = 'Taxon',
#' digit = digit)
#'
#' mod_morph_T_prc_Therm <- fit_all(data = dat_Trait_prc,
#' temperature = TRUE, precipitation = FALSE,
#' phenology = FALSE, morphology = TRUE,
#' condition = '2', nb_cores = nb_cores,
#' rand_trait = FALSE, fixed = 'Blood',
#' digit = digit)
#'
#'
#' ## forest plot as in Fig. 3
#' par(par_ini)
#' plot_T_cond2 <- plot_forest(meta_obj1 = mod_phen_T_prcs$meta_res,
#' meta_obj2 = mod_morph_T_prcs$meta_res,
#' list_extra_meta_obj =
#' list(mod_phen_T_prcs$meta_res,
#' mod_morph_T_prcs$meta_res,
#' mod_phen_T_prc$meta_res,
#' mod_phen_T_prc_Taxon$meta_res,
#' mod_morph_T_prc$meta_res,
#' mod_morph_T_prc_Taxon$meta_res),
#' sort = c('Species', 'Trait_Categ_det',
#' 'PaperID'),
#' increasing = TRUE,
#' labels = c(traits = TRUE,
#' fitness = FALSE,
#' country = FALSE,
#' authors = FALSE))
#' mtext('Across studies, ', at = -0.6, side = 2,
#' line = 4, las = 2, cex = 0.9, col = 'black')
#' mtext('PRCS dataset', at = -2.2, side = 2,
#' line = 4.3, las = 2, cex = 0.9, col = 'black')
#' mtext('Across studies, ', at = -7.4, side = 2,
#' line = 4, las = 2, cex = 0.9, col = 'black')
#' mtext('PRC dataset ', at = -9, side = 2,
#' line = 4.3, las = 2, cex = 0.9, col = 'black')
#' }
#'
#' ## plots of raw data and extracted slopes per study, as in Suppl. Fig. S12-S13
#' par(par_ini); par(oma = c(2, 2, 0, 0), mar = c(1, 2, 3, 1))
#' plot_raw(data = dat_Trait, temperature = TRUE,
#' precipitation = FALSE, phenology = TRUE,
#' morphology = FALSE, condition = '2',
#' id_to_do = c(1:20))
#'
#' par(par_ini); par(oma = c(2, 2, 0, 0), mar = c(1, 2, 3, 1))
#' plot_raw(data = dat_Trait, temperature = TRUE,
#' precipitation = FALSE, phenology = TRUE,
#' morphology = FALSE, condition = '2',
#' id_to_do = c(21:38))
#'
#' par(par_ini); par(oma = c(2, 2, 0, 0), mar = c(1, 2, 3, 1))
#' plot_raw(data = dat_Trait, temperature = TRUE,
#' precipitation = FALSE, phenology = FALSE,
#' morphology = TRUE, condition = '2',
#' id_to_do = c(1:4))
#'
#'
#'
#' #~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~#
#' #### Temperature: Condition 3 ####
#' #~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~#
#' ## phenology
#' mod_Sel_T_phen <- fit_all(data = dat_Sel, temperature = TRUE,
#' precipitation = FALSE, phenology = TRUE,
#' morphology = FALSE, condition = '3',
#' nb_cores = nb_cores, rand_trait = FALSE,
#' fixed = NULL, digit = digit)
#'
#' ## fixed effects
#' mod_Sel_T_phen_Fitn <- fit_all(data = dat_Sel, temperature = TRUE,
#' precipitation = FALSE, phenology = TRUE,
#' morphology = FALSE, condition = '3',
#' nb_cores = nb_cores, rand_trait = FALSE,
#' fixed = 'Fitness_Categ', digit = digit)
#'
#' mod_Sel_T_phen_Gener <- fit_all(data = dat_Sel, temperature = TRUE,
#' precipitation = FALSE, phenology = TRUE,
#' morphology = FALSE, condition = '3',
#' nb_cores = nb_cores, rand_trait = FALSE,
#' fixed = 'GenerationLength_yr',
#' digit = digit)
#' ## morphology
#' mod_Sel_T_morph <- fit_all(data = dat_Sel, temperature = TRUE,
#' precipitation = FALSE, phenology = FALSE,
#' morphology = TRUE, condition = '3',
#' nb_cores = nb_cores, rand_trait = FALSE,
#' fixed = NULL, digit = digit)
#' ## fixed effects
#' mod_Sel_T_morph_Fitn <- fit_all(data = dat_Sel, temperature = TRUE,
#' precipitation = FALSE, phenology = FALSE,
#' morphology = TRUE, condition = '3',
#' nb_cores = nb_cores, rand_trait = FALSE,
#' fixed = 'Fitness_Categ', digit = digit)
#'
#' mod_Sel_T_morph_Gener <- fit_all(data = dat_Sel, temperature = TRUE,
#' precipitation = FALSE, phenology = FALSE,
#' morphology = TRUE, condition = '3',
#' nb_cores = nb_cores, rand_trait = FALSE,
#' fixed = 'GenerationLength_yr',
#' digit = digit)
#' ## forest plot as in Fig. 4
#' par(par_ini)
#' plot_Sel_T_cond3 <- plot_forest(meta_obj1 = mod_Sel_T_phen$meta_res,
#' meta_obj2 = mod_Sel_T_morph$meta_res,
#' list_extra_meta_obj =
#' list(mod_Sel_T_phen$meta_res,
#' mod_Sel_T_morph$meta_res,
#' mod_Sel_T_phen_Fitn$meta_res),
#' sort = c('Species', 'Fitness_Categ', 'PaperID'),
#' increasing = TRUE,
#' labels = c(traits = TRUE,
#' fitness = TRUE,
#' country = FALSE,
#' authors = FALSE))
#' mtext('Across', side = 2, line = 4.5, at = -2.2,
#' las = 2, cex = 0.85, col = 'black')
#' mtext('studies', side = 2, line = 4.5, at = -3.4,
#' las = 2, cex = 0.85, col = 'black')
#'
#' ## plot of raw data for supplementary, as in S14-S15
#' par(par_ini); par(oma = c(2, 2, 0, 0), mar = c(1, 2, 3, 1))
#' plot_raw(data = dat_Sel, temperature = TRUE,
#' precipitation = FALSE, phenology = TRUE, morphology = FALSE,
#' condition = '3', id_to_do = c(1:20))
#'
#' par(par_ini); par(oma = c(2, 2, 0, 0), mar = c(1, 2, 3, 1))
#' plot_raw(data = dat_Sel, temperature = TRUE,
#' precipitation = FALSE, phenology = TRUE, morphology = FALSE,
#' condition = '3', id_to_do = c(21:36))
#'
#' par(par_ini); par(oma = c(2, 2, 0, 0), mar = c(1, 2, 3, 1))
#' plot_raw(data = dat_Sel, temperature = TRUE,
#' precipitation = FALSE, phenology = FALSE, morphology = TRUE,
#' condition = '3', id_to_do = c(1:9))
#'
#' ## change of selection over years, phenology
#' mod_Sel_T_phen_year <- fit_all(data = dat_Sel, temperature = TRUE,
#' precipitation = FALSE, phenology = TRUE,
#' morphology = FALSE, condition = '3b',
#' nb_cores = nb_cores, rand_trait = FALSE,
#' fixed = NULL, digit = digit)
#' ## change of selection over years, morphology
#' mod_Sel_T_morph_year <- fit_all(data = dat_Sel, temperature = TRUE,
#' precipitation = FALSE, phenology = FALSE,
#' morphology = TRUE, condition = '3b',
#' nb_cores = nb_cores, rand_trait = FALSE,
#' fixed = NULL, digit = digit)
#'
#' par(par_ini)
#' plot_Sel_T_overYrs <- plot_forest(meta_obj1 = mod_Sel_T_phen_year$meta_res,
#' meta_obj2 = mod_Sel_T_morph_year$meta_res,
#' list_extra_meta_obj =
#' list(mod_Sel_T_phen_year$meta_res,
#' mod_Sel_T_morph_year$meta_res),
#' sort = c('Species', 'Study_Authors',
#' 'Fitness_Categ'),
#' increasing = TRUE,
#' labels = c(traits = TRUE,
#' fitness = TRUE,
#' country = FALSE,
#' authors = FALSE))
#' mtext('Across', side = 2, line = 6, at = -1,
#' las = 2, cex = 0.8, col = 'black')
#' mtext('studies', side = 2, line = 6, at = -2.2,
#' las = 2, cex = 0.8, col = 'black')
#'
#'
#'
#' #~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~#
#' #### Assess whether responses are adaptive: products of slopes #
#' #### from conditions 1 & 2 vs. the slopes from condition 3 #
#' #### binomial test #
#' #~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~#
#'
#' # assess condition 3: whether response is adaptive?
#' ## extracting slopes
#' all_slopes_phen <- merge_3slopes(meta_obj = mod_T_prcs,
#' meta_obj1 = mod_phen_T_prcs,
#' meta_obj2 = mod_Sel_T_phen)
#'
#' all_slopes_morph <- merge_3slopes(meta_obj = mod_T_prcs,
#' meta_obj1 = mod_morph_T_prcs,
#' meta_obj2 = mod_Sel_T_morph)
#' ## testing adaptation
#' res_phen <- table(all_slopes_phen$slope.x * all_slopes_phen$slope.y > 0)
#'
#' (bc_phenT <- binom.confint(res_phen["TRUE"], sum(res_phen)))
#' (bc_phenF <- binom.confint(res_phen["FALSE"], sum(res_phen)))
#' binom.test(res_phen["TRUE"], sum(res_phen))
#'
#' res_morph <- table(all_slopes_morph$slope.x * all_slopes_morph$slope.y > 0)
#' (bc_morphT <- binom.confint(res_morph["TRUE"], sum(res_morph)))
#' (bc_morphF <- binom.confint(res_morph["FALSE"], sum(res_morph)))
#' binom.test(res_morph["TRUE"], sum(res_morph))
#'
#'
#' par(par_ini)
#' mat <- matrix(c(1, 2, 3, 3), 2,2, byrow = TRUE)
#' layout(mat, widths = rep.int(1, ncol(mat)),
#' heights = rep.int(1, nrow(mat)), respect = TRUE)
#' plot_slopes_products(data = merge_3slopes(meta_obj = mod_T_prcs,
#' meta_obj1 = mod_phen_T_prcs,
#' meta_obj2 = mod_Sel_T_phen),
#' phenological = TRUE)
#'
#' plot_slopes_products(data = merge_3slopes(meta_obj = mod_T_prcs,
#' meta_obj1 = mod_morph_T_prcs,
#' meta_obj2 = mod_Sel_T_morph),
#' phenological = FALSE, xlim = c(-0.1, 0.1),
#' ylim = c(-0.1, 0.1))
#'
#' par(par_ini); par(mar = c(4, 8, 3, 6))
#' theplot_Adapt <- barplot(cbind(c(res_phen["TRUE"]/sum(res_phen),
#' res_phen["FALSE"]/sum(res_phen)),
#' c(res_morph["TRUE"]/sum(res_morph),
#' res_morph["FALSE"]/sum(res_morph))),
#' names.arg = c("Phenology", "Morphology"),
#' las = 1, ylim = c(0, 1.1),
#' ylab = "Proportion of studies",
#' col = c("grey", "black"),
#' cex.lab = 1.3, cex.axis = 1.3, cex.names = 1.3)
#' legend(x = 0.8, y = 1, horiz = TRUE, fill = c("grey", "black"),
#' legend = c("adaptative", "maladaptive"), bg = 'white')
#' arrows(x0 = theplot_Adapt[1], x1 = theplot_Adapt[1],
#' y0 = bc_phenT[11, "lower"], y1 = bc_phenT[11, "upper"],
#' code = 3, angle = 90, col = "grey30", lwd = 4, length = 0.1)
#' arrows(x0 = theplot_Adapt[2], x1 = theplot_Adapt[2],
#' y0 = bc_morphT[11, "lower"], y1 = bc_morphT[11, "upper"],
#' code = 3, angle = 90, col = "grey30", lwd = 4, length = 0.1)
#' mtext('C)', side = 3, line = 0, adj = 0, cex = 2)
#'
#'
#'
#' #~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~#
#' #### Assess whether responses are adaptive ####
#' #### meta-analysis ####
#' #~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~#
#'
#' dat_AllTPhen <- prepare_data(data = dat_All, temperature = TRUE,
#' phenology = TRUE, morphology = FALSE)
#' dat_Adpt <- compute_adaptation(data = dat_AllTPhen)
#'
#' dat_Adpt$slope <- dat_Adpt$WMS * sign(dat_Adpt$prod12)
#' dat_Adpt$SE_slope <- dat_Adpt$se_WMS
#'
#' mod_Adapt_T_Phen <- fit_meta(meta_data = dat_Adpt, fixed = NULL)
#'
#'
#'
#' ## morphology
#' dat_AllTMorph <- prepare_data(data = dat_All, temperature = TRUE,
#' phenology = FALSE, morphology = TRUE)
#' dat_AllTMorph <- dat_AllTMorph[! is.na(dat_AllTMorph$ScaledSE), ]
#' dat_Adpt_Morph <- compute_adaptation(data = dat_AllTMorph)
#'
#' dat_Adpt_Morph$slope <- dat_Adpt_Morph$WMS * sign(dat_Adpt_Morph$prod12)
#' dat_Adpt_Morph$SE_slope <- dat_Adpt_Morph$se_WMS
#'
#' mod_Adapt_T_Morph <- fit_meta(meta_data = dat_Adpt_Morph, fixed = NULL)
#'
#' # using selection dataset only: to identify whether the sign
#' # is in the same direction as phenotypic change
#' dat_SelT <- prepare_data(data = dat_Sel, temperature = TRUE,
#' phenology = TRUE, morphology = FALSE)
#' dat_AllTPhen_Sel <- dat_All[dat_All$id %in% unique(dat_SelT$id), ]
#' dat_Adpt_Sel <- compute_adaptation(data = dat_AllTPhen_Sel)
#'
#' dat_Adpt_Sel$slope <- dat_Adpt_Sel$WMS * sign(dat_Adpt_Sel$prod12)
#' dat_Adpt_Sel$SE_slope <- dat_Adpt_Sel$se_WMS
#'
#' mod_Adapt_T_Sel <- fit_meta(meta_data = dat_Adpt_Sel, fixed = NULL)
#'
#' mod_Adapt_T_Sel_Sens <- fit_meta(meta_data =
#' dat_Adpt_Sel[dat_Adpt_Sel$slope < 0.6, ],
#' fixed = NULL)
#'
#'
#' ## morphology
#' dat_SelT_Morph <- prepare_data(data = dat_Sel, temperature = TRUE,
#' phenology = FALSE, morphology = TRUE)
#' dat_AllMorph_Sel <- dat_All[dat_All$id %in% unique(dat_SelT_Morph$id), ]
#' dat_AllMorph_Sel <- dat_AllMorph_Sel[! is.na(dat_AllMorph_Sel$ScaledSE), ]
#' dat_Adpt_Morph_Sel <- compute_adaptation(data = dat_AllMorph_Sel)
#'
#' dat_Adpt_Morph_Sel$slope <- dat_Adpt_Morph_Sel$WMS *
#' sign(dat_Adpt_Morph_Sel$prod12)
#' dat_Adpt_Morph_Sel$SE_slope <- dat_Adpt_Morph_Sel$se_WMS
#'
#' mod_Adapt_T_Morph_Sel <- fit_meta(meta_data = dat_Adpt_Morph_Sel, fixed = NULL)
#'
#' par(par_ini)
#' plot_adapt_T_Sel <- plot_forest(meta_obj1 = mod_Adapt_T_Sel,
#' meta_obj2 = mod_Adapt_T_Morph_Sel,
#' list_extra_meta_obj =
#' list(mod_Adapt_T_Sel, mod_Adapt_T_Morph_Sel),
#' sort = c('Species', 'Study_Authors',
#' 'Fitness_Categ'),
#' increasing = TRUE,
#' labels = c(traits = TRUE,
#' fitness = TRUE,
#' country = FALSE,
#' authors = FALSE),
#' mar = c(7, 10, 2, 2))
#'
#' mtext('the direction of climate-driven', side = 1, line = 4.2, cex = 1.3)
#' mtext('trait change', side = 1, line = 5.3, cex = 1.3)
#'
#'
#'
#' #~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~#
#' #### Precipitation: Condition 1 ####
#' #~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~#
#'
#' ## for PRCS
#' mod_Prec_prcs <- fit_all(data = dat_Clim, temperature = FALSE,
#' precipitation = TRUE, phenology = TRUE,
#' morphology = TRUE, condition = '1',
#' nb_cores = nb_cores, rand_trait = FALSE,
#' fixed = NULL, digit = digit)
#' ## for PRC
#' mod_Prec_prc <- fit_all(data = dat_Clim_prc, temperature = FALSE,
#' precipitation = TRUE, phenology = TRUE,
#' morphology = TRUE, condition = '1',
#' nb_cores = nb_cores, rand_trait = FALSE,
#' fixed = NULL, digit = digit)
#'
#'
#'
#' #~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~#
#' #### Precipitation: Condition 2 ####
#' #~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~#
#'
#' ## for PRCS
#' mod_phen_Prec_prcs <- fit_all(data = dat_Trait, temperature = FALSE,
#' precipitation = TRUE, phenology = TRUE,
#' morphology = FALSE, condition = '2',
#' nb_cores = nb_cores, rand_trait = FALSE,
#' fixed = NULL, digit = digit)
#'
#' ## for PRC
#' mod_phen_Prec_prc <- fit_all(data = dat_Trait_prc, temperature = FALSE,
#' precipitation = TRUE, phenology = TRUE,
#' morphology = FALSE, condition = '2',
#' nb_cores = nb_cores, rand_trait = FALSE,
#' fixed = NULL, digit = digit)
#'
#' ## fixed effects
#' mod_phen_Prec_prc_Taxon <- fit_all(data = dat_Trait_prc,
#' temperature = FALSE,
#' precipitation = TRUE,
#' phenology = TRUE, morphology = FALSE,
#' condition = '2', nb_cores = nb_cores,
#' rand_trait = FALSE, fixed = 'Taxon',
#' digit = digit)
#'
#'
#'
#'
#' #~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~#
#' #### Precipitation: Condition 3 ####
#' #~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~#
#' ## for phenology, because there is no data on morphological traits
#' ## in the dataset on precipitation
#' mod_Sel_Prec_phen <- fit_all(data = dat_Sel,
#' temperature = FALSE, precipitation = TRUE,
#' phenology = TRUE, morphology = FALSE,
#' condition = '3', nb_cores = nb_cores,
#' rand_trait = FALSE, fixed = NULL,
#' digit = digit)
#'
#' ## fixed effects
#' mod_Sel_Prec_phen_Fitn <- fit_all(data = dat_Sel, temperature = FALSE,
#' precipitation = TRUE, phenology = TRUE,
#' morphology = FALSE, condition = '3',
#' nb_cores = nb_cores, rand_trait = FALSE,
#' fixed = 'Fitness_Categ', digit = digit)
#'
#' mod_Sel_Prec_phen_Gener <- fit_all(data = dat_Sel, temperature = FALSE,
#' precipitation = TRUE, phenology = TRUE,
#' morphology = FALSE, condition = '3',
#' nb_cores = nb_cores, rand_trait = FALSE,
#' fixed = 'GenerationLength_yr',
#' digit = digit)
#'
#' ## supplementary plot of all effects for each condition (as Fig. S3)
#' par(par_ini); par(oma = c(2,1,1,0), mfrow = c(1,3))
#'
#' ## panel a)
#' plot_Prec_Cond1 <- plot_forest(meta_obj1 = mod_Prec_prcs$meta_res,
#' meta_obj2 = NULL,
#' list_extra_meta_obj =
#' list(mod_Prec_prcs$meta_res,
#' mod_Prec_prc$meta_res),
#' sort = c("Coord"),
#' increasing = FALSE,
#' labels = c(traits = FALSE,
#' fitness = FALSE,
#' country = TRUE,
#' authors = FALSE),
#' mar = c(4, 9, 2, 2))
#' mtext('Across', side = 2, line = 5,
#' at = -1.2, las = 2, cex = 0.9, col = 'black')
#' mtext('studies', side = 2, line = 5,
#' at = -1.8, las = 2, cex = 0.9, col = 'black')
#' mtext('a)', side = 2, at = 11.1, cex = 1.5, las = 2, line = 1)
#'
#' ## panel b)
#' plot_Prec_Trait_Cond2 <- plot_forest(meta_obj1 = mod_phen_Prec_prcs$meta_res,
#' meta_obj2 = NULL,
#' list_extra_meta_obj =
#' list(mod_phen_Prec_prcs$meta_res,
#' mod_phen_Prec_prc$meta_res),
#' sort = c('Species', 'Trait_Categ_det'),
#' increasing = FALSE,
#' labels = c(traits = TRUE,
#' fitness = FALSE,
#' country = FALSE,
#' authors = FALSE),
#' mar = c(4, 9, 2, 2))
#' mtext('Across', side = 2, line = 6,
#' at = -1.2, las = 2, cex = 0.8, col = 'black')
#' mtext('studies', side = 2, line = 6,
#' at = -1.8, las = 2, cex = 0.8, col = 'black')
#' mtext('b)', side = 2, at = 12.3, cex = 1.5, las = 2, line = 1.2)
#'
#' ## panel c)
#' plot_Sel_Prec_cond3 <- plot_forest(meta_obj1 = mod_Sel_Prec_phen$meta_res,
#' meta_obj2 = NULL,
#' list_extra_meta_obj =
#' list(mod_Sel_Prec_phen$meta_res,
#' mod_Sel_Prec_phen_Fitn$meta_res),
#' sort = c('Species', 'Study_Authors',
#' 'Fitness_Categ'),
#' increasing = FALSE,
#' labels = c(traits = TRUE,
#' fitness = TRUE,
#' country = FALSE,
#' authors = FALSE),
#' mar = c(4, 13, 2, 2))
#' mtext('Across', side = 2, line = 8,
#' at = -2, las = 2, cex = 0.8, col = 'black')
#' mtext('studies', side = 2, line = 8,
#' at = -3, las = 2, cex = 0.8, col = 'black')
#' mtext('c)', side = 2, at = 13.3, cex = 1.5, las = 2, line = 1.2)
#'
#'
#'
#' #~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~#
#' #### Sensitivity analysis ####
#' #~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~#
#'
#' ## excluding the study by Goodenough et al. (2011)
#' mod_phen_T_prcs_noOutlier <- fit_all(data = dat_Trait[
#' dat_Trait$Study_Authors != 'Goodenough', ],
#' temperature = TRUE, precipitation = FALSE,
#' phenology = TRUE, morphology = FALSE,
#' condition = '2', nb_cores = nb_cores,
#' rand_trait = FALSE, fixed = NULL,
#' digit = digit)
#'
#' ## selection
#' mod_Sel_T_phen_noOutlier <- fit_all(data = dat_Sel[
#' dat_Sel$Study_Authors != 'Goodenough', ],
#' temperature = TRUE, precipitation = FALSE,
#' phenology = TRUE, morphology = FALSE,
#' condition = '3', nb_cores = nb_cores,
#' rand_trait = FALSE, fixed = NULL,
#' digit = digit)
#'
#' ## fixed effects
#' mod_Sel_T_phen_Fitn_noOutlier <- fit_all(data = dat_Sel[
#' dat_Sel$Study_Authors != 'Goodenough', ],
#' temperature = TRUE, precipitation = FALSE,
#' phenology = TRUE, morphology = FALSE,
#' condition = '3', nb_cores = nb_cores,
#' rand_trait = FALSE, fixed = 'Fitness_Categ',
#' digit = digit)
#'
#' mod_Sel_T_phen_Gener_noOutlier <- fit_all(data = dat_Sel[
#' dat_Sel$Study_Authors != 'Goodenough', ],
#' temperature = TRUE, precipitation = FALSE,
#' phenology = TRUE, morphology = FALSE,
#' condition = '3', nb_cores = nb_cores,
#' rand_trait = FALSE,
#' fixed = 'GenerationLength_yr',
#' digit = digit)
#'
#' ## excluding the study on mammal
#' mod_phen_T_prcs_noMammal <- fit_all(data = dat_Trait[
#' dat_Trait$Study_Authors != 'Plard_et_al', ],
#' temperature = TRUE, precipitation = FALSE,
#' phenology = TRUE, morphology = FALSE,
#' condition = '2', nb_cores = nb_cores,
#' rand_trait = FALSE, fixed = NULL,
#' digit = digit)
#'
#' ## phenology
#' mod_Sel_T_phen_noMammal <- fit_all(data = dat_Sel[
#' dat_Sel$Study_Authors != 'Plard_et_al', ],
#' temperature = TRUE, precipitation = FALSE,
#' phenology = TRUE, morphology = FALSE,
#' condition = '3', nb_cores = nb_cores,
#' rand_trait = FALSE, fixed = NULL,
#' digit = digit)
#'
#' ## fixed effects
#' mod_Sel_T_phen_Fitn_noMammal <- fit_all(data = dat_Sel[
#' dat_Sel$Study_Authors != 'Plard_et_al', ],
#' temperature = TRUE, precipitation = FALSE,
#' phenology = TRUE, morphology = FALSE,
#' condition = '3', nb_cores = nb_cores,
#' rand_trait = FALSE, fixed = 'Fitness_Categ',
#' digit = digit)
#'
#' mod_Sel_T_phen_Gener_noMammal <- fit_all(data = dat_Sel[
#' dat_Sel$Study_Authors != 'Plard_et_al', ],
#' temperature = TRUE, precipitation = FALSE,
#' phenology = TRUE, morphology = FALSE,
#' condition = '3', nb_cores = nb_cores,
#' rand_trait = FALSE,
#' fixed = 'GenerationLength_yr',
#' digit = digit)
#'
#' ## excluding both the study on mammal and the study by Goodenough et al. (2011)
#' mod_phen_T_prcs_noBoth <- fit_all(data = dat_Trait[!
#' dat_Trait$Study_Authors %in%
#' c('Goodenough', 'Plard_et_al'), ],
#' temperature = TRUE, precipitation = FALSE,
#' phenology = TRUE, morphology = FALSE,
#' condition = '2', nb_cores = nb_cores,
#' rand_trait = FALSE, fixed = NULL, digit = digit)
#'
#' mod_Sel_T_phen_noBoth <- fit_all(data = dat_Sel[!
#' dat_Sel$Study_Authors %in%
#' c('Goodenough', 'Plard_et_al'), ],
#' temperature = TRUE, precipitation = FALSE,
#' phenology = TRUE, morphology = FALSE,
#' condition = '3', nb_cores = nb_cores,
#' rand_trait = FALSE, fixed = NULL,
#' digit = digit)
#'
#' ## fixed effects
#' mod_Sel_T_phen_Fitn_noBoth <- fit_all(data = dat_Sel[!
#' dat_Sel$Study_Authors %in%
#' c('Goodenough', 'Plard_et_al'), ],
#' temperature = TRUE, precipitation = FALSE,
#' phenology = TRUE, morphology = FALSE,
#' condition = '3', nb_cores = nb_cores,
#' rand_trait = FALSE, fixed = 'Fitness_Categ',
#' digit = digit)
#'
#' mod_Sel_T_phen_Gener_noBoth <- fit_all(data = dat_Sel[!
#' dat_Sel$Study_Authors %in%
#' c('Goodenough', 'Plard_et_al'), ],
#' temperature = TRUE, precipitation = FALSE,
#' phenology = TRUE, morphology = FALSE,
#' condition = '3', nb_cores = nb_cores,
#' rand_trait = FALSE,
#' fixed = 'GenerationLength_yr',
#' digit = digit)
#'
#'
#'
#' #~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~#
#' #### Sensitivity to the inclusion of abundance ####
#' #~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~#
#'
#' dat_PhenT_ab <- prepare_data(data = dat_Trait, temperature = TRUE,
#' precipitation = FALSE, phenology = TRUE,
#' morphology = FALSE)
#' dat_PhenT_ab <- dat_PhenT_ab[! is.na(dat_PhenT_ab$ScaledPop), ]
#'
#' # to exclude studies with < 11 years
#' num_years <- dat_PhenT_ab %>%
#' dplyr::group_by(., id) %>%
#' dplyr::summarise(num = dplyr::n())
#'
#' num_years <- num_years[num_years$num > 10, ]
#' dat_PhenT_ab <- dat_PhenT_ab[dat_PhenT_ab$id %in% num_years$id, ]
#'
#' EfSizes_phen_ab <- extract_effects_all_ids(data = dat_PhenT_ab, condition = '2b',
#' nb_cores = nb_cores)
#' EfSizes_phen_ab$slope <- EfSizes_phen_ab$slope_clim
#' EfSizes_phen_ab$SE_slope <- EfSizes_phen_ab$se_slope_clim
#' meta_phen_abund <- fit_meta(EfSizes_phen_ab, fixed = NULL,
#' rand_trait = FALSE, digit = 3)
#'
#'
#' ## and a model without abundance on this subset of data
#' mod_phen_T_prcs_subs <- fit_all(data = dat_PhenT_ab, temperature = TRUE,
#' precipitation = FALSE, phenology = TRUE,
#' morphology = FALSE, condition = '2',
#' nb_cores = nb_cores, rand_trait = FALSE,
#' fixed = NULL, digit = digit)
#'
#' ## plot as Supplementary Fig. S14 (for all slopes)
#' par(par_ini)
#' plot_abund_effects(meta_data = EfSizes_phen_ab)
#'
#'
#'
#' #~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~#
#' #### Supplementary: Funnel plots ####
#' #~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~#
#'
#' ## funnel plots for PRCS, as Supplementary Fig. S16
#' par(par_ini); par(mfrow = c(2, 3))
#' plot_funnel(meta_obj = mod_T_prcs)
#' mtext('a)', line = 1, adj = 0, cex = 1.7)
#'
#' plot_funnel(meta_obj = mod_phen_T_prcs)
#' mtext('b)', line = 1, adj = 0, cex = 1.7)
#'
#' plot_funnel(meta_obj = mod_morph_T_prcs)
#' mtext('c)', line = 1, adj = 0, cex = 1.7)
#'
#' plot_funnel(meta_obj = mod_Sel_T_phen, model = 'lm')
#' mtext('d)', line = 1, adj = 0, cex = 1.7)
#'
#' plot_funnel(meta_obj = mod_Sel_T_morph)
#' mtext('e)', line = 1, adj = 0, cex = 1.7)
#'
#'
#' ## funnel plots for PRC, as Supplementary Fig. S17
#' ## run only if the models were fit
#' \dontrun{
#' par(par_ini); par(mfrow = c(1,3))
#'
#' plot_funnel(meta_obj = mod_T_prc)
#' mtext('a)', line = 1, adj = 0, cex = 1.7)
#'
#' plot_funnel(meta_obj = mod_phen_T_prc)
#' mtext('b)', line = 1, adj = 0, cex = 1.7)
#'
#' plot_funnel(meta_obj = mod_morph_T_prc)
#' mtext('c)', line = 1, adj = 0, cex = 1.7)
#' }
#'
#' #~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~#
#' #### Supplementary: temperature over years VS ####
#' #### study duration VS first year ####
#' #~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~#
#'
#' ## make sure the model is fitted before running this model
#' \dontrun{
#' plot_eff_dur_firstY(data = dat_Clim_prc, meta_obj = mod_T_prc)
#' }
#'
#' #~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~#
#' #### Prepare Supplementary Tables ####
#' #~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~#
#'
#' # a list of models
#' # can only be run if all the models were fitted,
#' # some of which are time-consuming to fit
#' \dontrun{
#' mod_list <- list(mod_T_prcs, mod_phen_T_prcs, mod_phen_T_prcs_Type,
#' mod_morph_T_prcs, mod_morph_T_prcs_Type, mod_Sel_T_phen,
#' mod_Sel_T_phen_Fitn, mod_Sel_T_phen_Gener, mod_Sel_T_morph,
#' mod_Sel_T_morph_Fitn, mod_Sel_T_morph_Gener, mod_Prec_prcs,
#' mod_phen_Prec_prcs, mod_Sel_Prec_phen, mod_Sel_Prec_phen_Fitn,
#' mod_Sel_Prec_phen_Gener, mod_T_prc, mod_phen_T_prc,
#' mod_phen_T_prc_Taxon, mod_phen_T_prc_Type, mod_morph_T_prc,
#' mod_morph_T_prc_Taxon, mod_morph_T_prc_Type,
#' mod_morph_T_prc_Therm, mod_Prec_prc, mod_phen_Prec_prc,
#' mod_phen_Prec_prc_Taxon)
#'
#' ## A table with the effect sizes and their SE
#' tab_efSizes_ST1(model_list = mod_list)
#'
#' ## A table with LRT statistics
#' tab_LRT_ST2(model_list = mod_list)
#' }
#' # a list of models
#' mod_list_sens <- list(mod_phen_T_prcs_noOutlier,
#' mod_Sel_T_phen_noOutlier,
#' mod_Sel_T_phen_Fitn_noOutlier,
#' mod_Sel_T_phen_Gener_noOutlier,
#' mod_phen_T_prcs_noMammal,
#' mod_Sel_T_phen_noMammal,
#' mod_Sel_T_phen_Fitn_noMammal,
#' mod_Sel_T_phen_Gener_noMammal,
#' mod_phen_T_prcs_noBoth,
#' mod_Sel_T_phen_noBoth,
#' mod_Sel_T_phen_Fitn_noBoth,
#' mod_Sel_T_phen_Gener_noBoth)
#'
#' tab_efSizes_ST4(model_list = mod_list_sens)
#'
#' ## a table for the LRT statistics and the variation due to random effects
#' ## for the sensitivity analysis
#' tab_LRT_ST5(model_list = mod_list_sens)
#'
#' ## a table for heterogeneity
#' ## make sure all models are fitted before running this
#' \dontrun{
#' mod_heterogen <- list(mod_T_prcs, mod_phen_T_prcs, mod_morph_T_prcs,
#' mod_Sel_T_phen, mod_Sel_T_morph, mod_Prec_prcs,
#' mod_phen_Prec_prcs, mod_Sel_Prec_phen, mod_T_prc,
#' mod_phen_T_prc, mod_morph_T_prc, mod_Prec_prc,
#' mod_phen_Prec_prc, mod_phen_T_prcs_noOutlier,
#' mod_Sel_T_phen_noOutlier, mod_phen_T_prcs_noMammal,
#' mod_Sel_T_phen_noMammal, mod_phen_T_prcs_noBoth,
#' mod_Sel_T_phen_noBoth)
#' tab_heterog(model_list = mod_heterogen)
#' }
#'
NULL
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.