DOI: 10.5061/dryad.8w9ghx3tx
title: Dataset and scripts from: Predicting organismal response to marine heatwaves using dynamic thermal tolerance landscape models

description: 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. 
description: 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. <strong>Creation of MHW
 profiles: </strong>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.  <strong>Assessment of
 MHW strength using statistical methods:</strong> 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.  <strong>Simulation of three
 organism's thermal death time curves:</strong> 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. 
 <strong>Assessment of MHW strength using physiological
 dynamic tolerance models:</strong> 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. 
 <strong>Assessment of the role of MHW timing on
 survival:</strong> 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.
description: # 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](mailto: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