Dataset and scripts from: Predicting organismal response to marine heatwaves using dynamic thermal tolerance landscape models
Data files
May 09, 2024 version files 2.51 GB
-
code_data.zip
-
README.md
Abstract
Marine heatwaves (MHWs) can cause thermal stress in marine organisms, experienced as extreme ‘pulses’ against the gradual trend of anthropogenic warming. When thermal stress exceeds organismal capacity to maintain homeostasis, organism survival becomes time-limited and can result in mass mortality events. Current methods of detecting and categorizing MHWs rely on statistical analysis of historic climatology, and do not consider biological effects as a basis of MHW severity. The reemergence of ectotherm thermal tolerance landscape models provides a physiological framework for assessing the lethal effects of MHWs by accounting for both the magnitude and duration of extreme heat events. Here, we used a simulation approach to understand the effects of a suite of MHW profiles on organism survival probability across 1) three thermal tolerance adaptive strategies, 2) interannual temperature variation, and 3) seasonal timing of MHWs. We identified survival isoclines across MHW magnitude and duration where acute (short duration-high magnitude) and chronic (long duration-low magnitude) events had equivalent lethal effects on marine organisms. While most research attention has focused on chronic MHW events, we show similar lethal effects can be experienced by more common but neglected acute marine heat spikes. Critically, a statistical definition of MHWs does not accurately categorize biological mortality. By letting organism responses define the extremeness of a MHW event, we can build a mechanistic understanding of MHW effects from a physiological basis. Organism responses can then be transferred across scales of ecological organization and better predict marine ecosystem shifts to MHWs.
README: Marine Heatwave Simulation Analysis
Code and data to replicate results in Villeneuve and White 2024, Journal of Animal Ecology.
To set correct relative file paths, download the entire folder onto your local machine, unzip, and open the R project mhw_simulation_manuscript2024.rproj.
There are three major R markdown scripts that replicate analyses and figures from our paper:
- mhw_sim_analysis.Rmd - Bulk of analysis. Can knit entire document, or run select chunks. Some chunks are very time and computationally intensive, so we provide .RDs data files that can be read in at checkpoints to speed rendering of the document. These data files can be replicated by running chunks that are not evaluated in the markdown. A rendered .html file of the markdown is also available.
- conceptual_figure.Rmd - Script to produce conceptual figure panels in the manuscript. A rendered .html file of the markdown is also available.
- tdt_premise_extended.Rmd - Script to produce supplementary figure 1, which is an extension of the species comparison analysis from the main text. This is computationally intensive, and while all code to replicate results is present, we provide a single data file checkpoint to greatly speed the knitting of the figure.
Additionally, there are two R scripts containing functions from other works used in our analysis.
- hour_fxns_heatwaveR - modified functions from Schlegel et al. 2018, the heatwaveR package. Functions altered to allow for inputs of data in hourly format. More information within mhw_sim_analysis.Rmd and manuscript text.
- Thermal_landscape_functions.R - modified functions from Rezende et al. 2020 to create thermal tolerance landscapes and produce estimates of mortality using dynamic survival model. More information within mhw_sim_analysis.Rmd and manuscript text.
There is one additional file provided, references.bib, that helps render citations within the rmarkdown documents.
Finally, there are two folders.
- data - .RDs files that are saved during computationally intensive sections of the analysis, so that they can be read in independent of running their origin code. These files are not necessary to replicate our results, but do make it much faster. Note that two .RDs files, merged_temp_tri.RDs and merged_temp_tri_ext.RDs, are not included in this folder due to their large size (5GB). If you want these .RDs files, and cannot render the code due to computational reasons, email Drew at drew.villeneuve@unh.edu
- result_figs - folder for storing figures from analysis. These folders are empty, but will become populated with figures once scripts are run. All figures from the main document and supplement can be replicated using the R scripts (see above) and will be stored here.
Software and package versions
Below is a readout from sessionInfo() detailing the R version and package versions used in this analysis. R was run in Rstudio version 2023.12.1. Version history recorded May 7, 2024.
R version 4.3.3 (2024-02-29 ucrt)
Platform: x86_64-w64-mingw32/x64 (64-bit)
Running under: Windows 11 x64 (build 22631)
Matrix products: default
locale:
[1] LC_COLLATE=English_United States.utf8 LC_CTYPE=English_United States.utf8 LC_MONETARY=English_United States.utf8
[4] LC_NUMERIC=C LC_TIME=English_United States.utf8
time zone: America/New_York
tzcode source: internal
attached base packages:
[1] grid parallel stats graphics grDevices utils datasets methods base
other attached packages:
[1] plotly_4.10.1 simglm_0.8.9 gridGraphics_0.5-1 rerddap_1.0.2 forcats_1.0.0 stringr_1.5.0 readr_2.1.4
[8] tidyverse_2.0.0 SpatialTools_1.0.5 terra_1.7-46 ggpubr_0.5.0 ggside_0.2.2 heatwaveR_0.5.3.9000 reshape2_1.4.4
[15] rlang_1.1.1 doParallel_1.0.17 iterators_1.0.14 cowplot_1.1.1 zoo_1.8-11 data.table_1.15.4 ggforce_0.4.1
[22] foreach_1.5.2 metR_0.14.0 purrr_1.0.2 kableExtra_1.4.0 ggplot2_3.5.0 tidyr_1.3.0 tibble_3.2.1
[29] lubridate_1.9.3 dplyr_1.1.4
loaded via a namespace (and not attached):
[1] magrittr_2.0.3 compiler_4.3.3 systemfonts_1.0.6 vctrs_0.6.5 httpcode_0.3.0 pkgconfig_2.0.3 fastmap_1.1.1
[8] backports_1.4.1 magic_1.6-1 utf8_1.2.2 spBayes_0.4-6 rmarkdown_2.26 tzdb_0.3.0 xfun_0.41
[15] cachem_1.0.6 jsonlite_1.8.8 tweenr_2.0.2 broom_1.0.5 R6_2.5.1 stringi_1.7.8 parallelly_1.37.1
[22] hoardr_0.5.2 car_3.1-1 pkgload_1.3.4 Rcpp_1.0.11 knitr_1.45 future.apply_1.11.1 Matrix_1.5-3
[29] timechange_0.1.1 tidyselect_1.2.1 rstudioapi_0.16.0 abind_1.4-5 yaml_2.3.6 codetools_0.2-18 curl_4.3.2
[36] listenv_0.9.1 lattice_0.20-45 plyr_1.8.9 withr_3.0.0 coda_0.19-4 evaluate_0.23 future_1.33.2
[43] isoband_0.2.7 polyclip_1.10-4 xml2_1.3.3 pillar_1.9.0 carData_3.0-5 checkmate_2.1.0 ncdf4_1.20
[50] generics_0.1.3 sp_2.0-0 dtplyr_1.2.2 hms_1.1.3 munsell_0.5.0 scales_1.3.0 globals_0.16.3
[57] glue_1.6.2 lazyeval_0.2.2 tools_4.3.3 ggsignif_0.6.4 colorspace_2.0-3 Formula_1.2-4 cli_3.6.2
[64] rappdirs_0.3.3 fansi_1.0.3 viridisLite_0.4.2 svglite_2.1.3 gtable_0.3.4 rstatix_0.7.1 digest_0.6.29
[71] crul_1.4.0 htmlwidgets_1.6.4 farver_2.1.1 memoise_2.0.1 htmltools_0.5.8 lifecycle_1.0.4 httr_1.4.7
[78] MASS_7.3-58.1
Methods
The marine heatwave dataset used in this study was simulated using the included code. While the entirity of each script can be run from start to finish, we have silenced strategic chunks in each .rmd document and added data read-in lines to speed up script execution. Running silenced chunks of each .rmd file can result in long computation times. Detailed information about data and script files can be found in the README.md document.
Creation of MHW profiles: We simulated MHW profiles that covered a range of magnitudes (0-8 °C, 0-20 °C in extended analysis script, 0.5 °C interval) and durations (1 hour - 30 days, 1 day interval). We accomplished this by first simulating the MHW component of a temperature time series in the shape of a triangle - that is, a ramp up and ramp down period of equal lengths and no plateau period at the maximum magnitude. We next simulated a climatological time series over the year 2022 using a sine wave with a mean of 12°C and a range of 14°C, thus giving a minimum temperature of -2°C and a maximum temperature of 26°C. The maximum temperature of this sine wave occurred on 17th, and we added the MHW triangle on top of the sine wave so that the MHW maxima occurred directly over the climatological sine wave maxima. This gave us 527 unique MHW profiles of varying magnitude and duration.
Assessment of MHW strength using statistical methods: Once we simulated the MHWs, we used the statistical MHW categorization techniques developed by Hobday et al. 2018 and implemented into the heatwaveR package by Schlegel and Smit 2018. We slightly altered heatwaveR functions to allow for hourly durations. This method relies on historical (30 year) variability at a site to determine the 90th percentile of temperatures, which it then uses as a threshold to delineate MHWs from non-MHWs and assign categorizations based on the magnitude exceedance of this threshold. Since we simulated our temperature data, we set three different 90th percentile values of 0.25, 0.75, and 1.5°C above climatology, from which we built three additional thresholds delineating events that were 2x, 3x, and 4x the difference between the 90th percentile threshold and climatology. For each interannual variability manipulation, we assess MHW category based on the Hobday 2018 framework for each of the 527 MHW profiles.
Simulation of three organism's thermal death time curves: We simulated three thermal death time curves that intersect at three different time and temperature magnitudes and thus give three different adaptive strategies - an acute tolerator, with a high CTmax and high z, a chronic tolerator, with a low CTmax and low z, and a mixed strategy underperforming species with low CTmax and a medium level of z. These parameters are derived from the shape of thermal death time curves, a method of assessing thermal tolerance changes with time that is highly conserved across taxa (Rezende et al. 2014). These TDT curves intersect at values that should result in changes in ranking of the thermal performance of the organisms throughout our MHW simulations. The parameter values of each TDT curve approximates different marine bivalve species' TDT curves.
Assessment of MHW strength using physiological dynamic tolerance models: We 'exposed' each of the three hypothetical organisms to each MHW profile and calculated survival probability using dynamic tolerance models formalized by Rezende et al. 2020. We extracted the final survival estimate from each organism at the end of each MHW profile.
Assessment of the role of MHW timing on survival: As a final manipulation, we simulated heatwaves of 24 day duration but ranging in magnitude from 0-8°C across ten timepoints, starting from being centered exactly over the climatological maxima (sine wave maxima) to 90 days after this summer maxima. We repeated step 4 on these MHW profiles to extract final survival estimates of the MHW for one species, the acute tolerator.