TrenchR: An R package for modular and accessible microclimate and biophysical ecology

Much understanding of organismal responses to climate change and variability relies on the assumption that body temperatures are equal to temporally averaged air temperatures high above the ground. However, most organisms experience microclimates near the ground and acute exposure to solar and thermal radiation and thermal extremes can substantially elevate or depress their body temperatures. We introduce the TrenchR package, which aids in Translating Environmental Change into organismal responses. The package includes microclimate models to vertically scale weather station data to organismal heights. Additional functions model and temporally partition air and soil temperatures and solar radiation. TrenchR biophysical modeling tools include both general models for heat flows and specific models to predict body temperatures for a variety of ectothermic taxa. We also offer utility functions to aid in estimating the organismal and environmental parameters needed for bio-physical ecology. TrenchR focuses on simple and modular functions so users can create transparent and flexible models for biophysical applications. The package aims to introduce and enable microclimate and biophysical modeling to improve ecological and evolutionary forecasting. We further this aim through a series of educational modules that introduce the field of biophysical ecology.


Introduction
Responses of organisms and ecosystems to climate change are heterogeneous and thus inconsistent with current predictive models [1].Some of these predictive shortcomings stem from omitting spatial and temporal environmental variation and how it interacts with organismal phenotypes (e.g., size, coloration) [2,3].Many analyses assume that organismal body temperatures are equal to shaded air temperatures inside a screen at weather station height (usually 2 meters).However, air temperatures near the ground where most organisms reside are often considerably warmer, and absorption of solar and thermal radiation can raise body temperatures well above air temperatures, resulting in potential discrepancies of tens of degrees Celsius [4][5][6].The rapid drop in performance and onset of thermal stress at high body temperatures can amplify the biological significance of underestimations [7].
Applications are increasingly demonstrating the importance of accounting for microclimate variation [8][9][10].Acute thermal stress events are often more relevant to the physiology and energy balance of organisms than daily, monthly, or annually averaged environmental conditions.Yet, daily and seasonal environmental variation is often neglected [11,12].We aim to help remedy these shortcomings by introducing the TrenchR package (https://github.com/trenchproject/TrenchR) as a component of the TrEnCh project (https://www.trenchproject.com), which builds computational and visualization tools to Translate Environmental Change into organismal responses.
The mathematical tools for leveraging environmental data to predict organismal conditions have long lingered in books [13,14] and articles [15,16], but adoption of these tools has not kept pace with research on climate change responses.In recent decades, adoption of biophysical approaches is expanding with recognition of the importance of spatial and temporal environmental variation to determining climate change responses [17,18] and the increased dissemination of computational tools.
Estimating how organisms experience their thermal environment generally entails two classes of models [13,14] ( Fig 1).Microclimate models allow scaling conditions from sensors to organism height.They can characterize heat and air transport to estimate vertical air and soil temperature and wind profiles.They can also characterize direct and diffuse solar radiation and longwave radiation emanating from the sky, ground, and surrounding objects like vegetation.Microclimate model output is used as input in biophysical models, which estimate body temperatures using an energy budget to balance heat exchanges between organisms and their environment.
The TrenchR package is intended to complement the NicheMapR package, which includes sophisticated microclimate [18] and biophysical [20] models.Although the Fortran source code was recently released, the complexity of NicheMapR functions can make it difficult to understand and modify the source code and Fortran can be complicated to compile and run.The integrated functions made it difficult to separate and adapt different forms of heat flow, but we note that some modular R functions were made available during the review process for this manuscript that enable comparison with the TrenchR package (S1 Text).Many of the approaches are similar between the packages, but TrenchR provides simple R functions to aid the understanding and accessibility of biophysical approaches.The TrenchR functions are modular and easily adapt to a variety of organisms and research questions.
Although many of the TrenchR functions are general, some components of TrenchR are best suited for small ectothermic animals since body temperatures are generally assumed to be uniform and at steady-state.We omit plant specific biophysical models since they are the focus of several R packages [tealeaves: [21]; plantecophys: [22]].Some energy budget calculation components of biophysical models are also available in the ThermImage package [23].Microclimate models, with an emphasis on describing spatial variation, are also available in the microclima package [24].The microclimc package enables accounting for how forest canopies alter microclimates [25].
The development of increasingly complex microclimate and biophysical models before open and reproducible science was emphasized or feasible has limited their uptake.Yet, such models are crucial to improving understand and prediction of species responses to climate change.Simple and modular functions can be combined and extended in a null modeling approach until the model adequately describes organismal interactions with their environment while remaining accessible and transparent.Many applications will warrant the inclusion of more complex microclimate and biophysical algorithms from R packages such as NicheMapR, microclima, and microclimc.

Methods and features
The TrenchR package (http://trenchproject.github.io/TrenchR/)aims to promote accessibility and reproducibility.We welcome contributions and corrections from users.Our package was built using the devtools methodology (https://github.com/r-lib/devtools)with version control managed in Github.Issues and feature requests can be contributed in Github (https://github.com/trenchproject/TrenchR).The package is available via CRAN and Github.
Validation against sensor data including from physical models of organisms is essential to microclimate and biophysical modelling.We adapted well established and validated functions from biophysical ecology textbooks [13,14] and research articles.Source references are cited in function headers and reference sections.Many of the references, particularly those developing taxa-specific biophysical models, contain extensive empirical validation.S2 Text illustrates validation of biophysical models using a physical model of a grasshopper.S1 Text compares TrenchR functions to NicheMapR functions.We have incorporated default parameters (e.g., organism emissivity or heat conduction rates) when general values are available.Many functions additionally include comments that describe potential parameterizations.The biophysical ecology texts [13,14] include many tables with parameter values and we provide several data tables describing solar and thermal absorptivity in TrenchR.

Components
TrenchR functions are organized into the following categories (Fig 2 ): Utility functions.Calculate environmental metrics that form the basis of microclimate and biophysical models (e.g., zenith and declination angles, which describe the angles of A) Microclimate models scale air and soil temperatures and wind speeds from sensor to organismal heights.B) Biophysical models balance heat exchanges between organisms and their environment to estimate body temperatures [19].A thermal image depicts how grasshopper body temperatures vary substantially from vegetation, air, and ground temperatures due to heat exchange.https://doi.org/10.1371/journal.pclm.0000139.g001incident sunlight), facilitate modeling diurnal temperature variation, and estimate biologically relevant aggregate metrics such as degree days available for organismal development.
Allometric functions.Allow estimating the dimensions of organisms needed for energy balances and other analyses.Available functions can convert between organismal length, mass, surface area, and volume as well as estimate the silhoutte area, which describes the organismal area exposed to solar radiation.
Microclimate functions.Facilitate calculating the environmental conditions experienced by organisms.Temperature and wind profile functions scale environmental conditions from sensors to organismal height.Radiation functions allow estimating incoming solar radiation and partitioning the variation diurnally and across components (i.e., direct, diffuse, reflected).
General biophysical functions.Provide generalized models of heat exchanges between organisms and their environment so that users can build custom biophysical models.The functions allow implementing an energy balance including the following components of heat exchange between organisms and their environment: • radiative heat exchanges of solar and thermal radiation; • convective heat exchanges between organisms and their surrounding fluid (air or water) driven by fluid flow; • conductive heat exchanges between organisms and solid surfaces (generally the ground) due to physical contact; and • evaporative and metabolic heat exchanges associated with organisms' evaporative water loss and metabolic heat generation.
Additional functions aggregate these forms of heat exchange into energy balances and to use the energy balances to predict body temperatures.
Specific biophysical functions.provide biophysical models that have been built for particular organisms based on their physical properties, behavior, environment, and life history.We currently provide published biophysical models for lizards, salamanders, butterflies, grasshoppers, limpets, mussels, and snails.Most models predict operative environmental temperatures, which are the steady-state body temperatures of organisms with specified physical properties in a specific microclimate and assume no heat exchange via metabolism or evaporation [26].However, we also present an analytical function for humid operative temperature that incorporates the effects of evaporative cooling on operative temperatures for wet-skinned ectotherms, such as salamanders.

Vignettes and introductory tutorials
We introduce the functions in several vignettes.A good place to start is the Allometry and Conversions vignette, which provides tools for preparing data, such as estimating additional dimensions of organisms from measured dimensions.The Estimating Microclimates vignette provides resources for estimating the environmental conditions experienced by organisms.This includes estimating solar radiation and its components, diurnal variation in temperature and radiation, temperature and wind speed profiles, and soil temperatures and profiles.Finally, the core biophysical modeling functions are described in a tutorial on Using Energy Balances to Estimate Body Temperatures.Components of an energy budget can be estimated using individual functions and then body temperatures can be solved for using either a generic energy balance or taxa specific biophysical models.We additionally offer a List of Symbols used in equations.
For additional background on microclimate and biophysical modeling, we have updated a series of tutorials entitled Physical Processes in Ecosystems.We intend the tutorials to provide a more contained and accessible introduction to microclimate and biophysical modeling than that included in classic textbooks [13,14].The tutorials provide less detailed coverage of many of the topics contained in the books.The series of 14 tutorials starts with an overview of the calculus and physics principles underlying the modeling.
Thermodynamics and energy budget modeling are then introduced with detailed examples.Tutorials address the climate space concept, operative temperatures, biophysical models for leaves and sheep, and heat flow in soils.We have revised the tutorials (which originated from an NSF training grant in 1979 and lacked broad distribution) to include R code and utilization of TrenchR functions.We expanded the original series to include tutorials contributed by M. Kearney introducing the Microclim environmental data [27] and integrating the Microclim data with biophysical modeling to examine broad scale climatic limits.The tutorials are provided in html form and available for download as a pdf at the bookdown server (https:// bookdown.org/huckley/Physical_Processes_In_Ecosystems/).R markdown files for the tutorials are available in Github (https://github.com/trenchproject/TrenchRmodules).

Results
We illustrate the use of the TrenchR package by estimating an energy budget for a Sceloporus lizard on June 1, 2021 in Santa Fe, New Mexico, USA (35.69˚N, -105.944˚W,elevation: 2121 m).The simplified example, which is designed to be self-contained, is also incorporated in the Using energy balances to estimate body temperatures vignette.We start by generating environmental inputs (Fig 2).Using these inputs, we estimate the energy budget with component functions.Finally, we use an integrated biophysical model to estimate operative environmental temperatures.See S2 Text for a more realistic example examining a time series of microclimate data and applying and testing the biophysical models.
We will use the energy budget to estimate body temperature, T b , which can depart dramatically from the air temperatures due to heat exchange with the environment.Heat energy is exchanged with the environment by way of solar and thermal radiation, metabolic reactions, and evaporation.The organism also exchanges heat with the surrounding air or water via convection and with substrate it is in contact with via conduction.The balance of these heat exchanges (omitting metabolism and evaporation, which are often negligible for ectotherms) can be estimated and often referred to as operative environmental temperature, T e [26].T e is an estimate of T b and the package functions refer to T b for simplicity.Our exemplar estimation of T b assumes steady-state thermal conditions.Additionally, we assume that the lizard's body temperature is homogenous, which is generally reasonable for small ectotherms.Approaches to account for thermal gradients between the animal's core and its skin are available elsewhere [28].
Let us assume the lizard is in an unshaded location where a weather station at standard height (2 meters) reports that the daily air temperature varies from a minimum of 10˚C to a maximum of 25˚C and the wind speed averages 1 m/s.The soil surface temperature varies from a minimum of 15˚C to a maximum of 30˚C.We assume that atmospheric transmissivity τ = 0.7 and albedo ρ = 0.6.

Environmental data
At the first stage, we prepare the environmental data for analysis.We will estimate hourly air and soil temperatures and radiation using a function describing diurnal temperature variation.We start by estimating the day of year and the timing of sunrise and sunset: Although measured solar radiation is preferable if available, we can estimate hourly solar radiation by discounting incoming solar radiation as it moves through the atmosphere as follows.We use the approach from Cambpell and Norman [14], which uses an empirical relation to partition radiation into direct, diffuse, and reflected components.The partition_ solar_radiation() function includes 8 empirical relationships for, and the propor-tion_diffuse_solar_radiation() includes a more complex numerical approximation for, partitioning radiation components as described in the Estimating microclimates vignette.We then calculate hourly air and soil surface temperatures based on daily minimum and maximum temperatures.We select the sine-exponential model for air temperature and the sine model for surface temperature [29]:

# Soil surface temperature (C)
Ts <-sapply(1:24, diurnal_temp_variation_sine, T_max = Tmax_s, T_min = Tmin_s) At the second stage, we use microclimate models to scale air temperature (T r ) and wind speed (u r ) from weather station height (reference height z r = 2 m) to lizard height (organism height z = 0.02 m).We assume a surface roughness of z 0 = 0.2 m, which corresponds to bare sand and determines the turbulence of airflow.We implement free air temperature and wind speed profiles driven by density differences but profiles forced by wind speed are also available.

Energy balance
Finally, we will use our microclimates estimates to solve the following energy balance to estimate T e : where Q net is the net energy exchange with the environment (W), Q abs is the solar radiation absorbed (W), Q emit is the net thermal radiation emitted (W), Q conv is energy exchange due to convection (W), Q cond is energy exchange due to conduction (W), Q met is the energy generated by metabolism (W), and Q evap is the energy generated by evaporative water loss (W).We will estimate each term on the right side of the equation in turn.Estimating Q abs requires the surface area exposed to radiation and the solar absorptivity of the animal surface (a proportion).
We use zenith angle psi to estimate the projected (silhouette) area as a portion of the surface area of the organism, which allows estimating absorbed solar radiation.We model a 10 gram Sceloporus lizard with solar absorptivity a = 0.9 [13].We will initially assume T b = T a + 10 to illustrate the calculations before solving for T b given the environmental conditions.

# Estimate surface area (m^2) and the proportion sihouette area
A <-surface_area_from_mass(mass, "lizard") psa <-sapply(psi_deg, proportion_silhouette_area, taxon = "lizard", posture = "prostrate") We estimate thermal radiation Q emit (W) for the lizard outdoors, where psa dir and psa ref are the view factors, also refered to as configuration factors, that indicate the proportions of surface area A (m 2 ) exposed to the sky and ground, respectively.We assume the surface emissivity of lizards, epsilon s = 0.965 [30].We next estimate convection Q conv (W) and conduction Q cond (W).We will estimate the lizard's heat transfer coefficient, H L (Wm -2 K -1 ) using an empirical relationship for lizards (heat_transfer_coefficient()).We average thermal conductivity and kinematic viscosity across the day for simplicity and since there is not substantial diurnal variation.We also illustrate a function estimating H L using a spherical approximation (heat_transfer_ coefficient_approximation()) and a simplified approximation (heat_ transfer_coefficient_simple()) for cases when taxon specific relationships for estimating heat transfer coefficients are not available.We estimate the characteristic dimension, which determines exposure to convective heat exchange as the cube root of volume, assuming the animal density approximates that of water [31].
These coefficients assume convection is forced by the wind.TrenchR includes approaches for free convection and a function (free_or_forced_convection()) that evaluates whether free or forced convection is appropriate.The function uses dimensionless numbers, which have been developed to describe heat transfer coefficients associated with convection over different geometries and can be estimated using TrenchR (e.g., Grashof, Nusselt, and Reynolds numbers).
# Use DRYAIR from NicheMapR to estimate the thermal conductivity of air and kinematic viscosity.We estimate convective heat exchange between the animal and surrounding air using the following relationship: where an enhancement factor, ef, multiplier can be incorporated to account for increases in heat exchange resulting from air turbulence in field conditions.We implement the function in R assuming that 2/3 of the lizard's surface area is exchanging heat through convection.

}
We estimate conductive heat flow (W) from the lizard to the surface assuming conductance through the animal tissue is the rate limiting step as follows: where K skin in the thermal conductivity of lizard skin (Wm -2K -1).We implement the estimate assuming that conductive heat exchange occurs down to a soil depth of 2.5cm.We use this value rather than skin thickness, which results in rapid conduction and does not readily reach steady state conditions.We assume, as is generally done for lizards, that heat exchange associated with metabolism and evaporation is negligible.However, functions for estimating both forms of heat exchange available in TrenchR.
Qmet <-0 Qevap <-0 The full heat budget can be calculated as follows [13]: We also implement a similar but simplified energy balance [14].The energy balance omits conduction with the ground: We additionally estimate T b using a specialized function for lizards [32], where F d , F r , F a , and F g are the view factors between the surface of the lizard and diffuse solar radiation, reflected solar radiation, atmospheric thermal radiation, and ground thermal radiation, respectively:

}
The microclimate models indicate that air temperatures at lizard height are similar to surface temperatures (Fig 3).The biophysical models indicate that solar radiation can elevate lizard body temperatures far above air temperatures and that the lizard will face thermal stress if it is unable to seek shade (Fig 3).The three biophysical models predict different body temperatures during peak period of solar radiation because they model interactions with the ground differently and users are encouraged to review the details of each biophysical model and perform empirical validations before selection.Differences in estimated body temperatures are accentuated by the high level of solar radiation.Tb_Gates is a general and comprehensive model that is appropriate for many applications.Taxa-specific biophysical models often best account for details of organism environment interactions and have generally been well tested.

Discussion
TrenchR is intended to promote understanding of how organisms interact with their environment and consequences for physiology, energetics, behavior, and demography.Our example implementation highlights the importance of considering organismal body temperatures, rather than air temperatures, when examining thermal stress and other responses to environmental variability and change.TrenchR currently focuses on heat balances but may be expanded to include water balances.Simple functions can be combined as needed to produce comprehensive and transparent models for biophysical ecology and evolution.The resultant models are likely to be sufficiently detailed for many applications, but users are referred to NicheMapR for more detailed biophysical models [18,20].We focus on models that predict steady-state conditions for simplicity (that is, steady-state conditions).Such models generally do not present computational challenges so our models are not optimized for computational efficiency.Making classic biophysical ecology techniques more accessible will allow researchers to take advantage of rapidly accumulating data on environmental conditions and organismal traits to understand and predict ecological and evolutionary responses.Considering how organisms experience their environment is central to understanding responses to variable and changing environments [33].

Supporting information
S1 Text.We compare TrenchR and NicheMapR implementations.We thank Michael Kearney for making modular R functions corresponding to NicheMapR available and providing code comparing TrenchR and NicheMapR functions during the review process.Our adaptation of the provided code corresponds to the operative temperature estimation in the manuscript.The NicheMapR implementation relies on the micro_global environmental data available through NicheMapR.The modular NicheMapR functions are available in version 3.2.1 via GitHub (https://github.com/mrke/NicheMapR/releases).(DOCX) S2 Text.We provide an additional example using microclimate data measured along an elevation gradient in CO, USA to illustrate use of TrenchR.We first use a time series of air temperatures and wind speeds collected at multiple heights to examine profiles and surface roughness.We then use a time series of air and surface temperatures, wind speeds, and solar radiation collected at a single height to implement a biophysical model for grasshoppers and compare estimates to observations.We omit the code to read and process the environmental data here for brevity.R markdown files with the full code and the associated environmental data are available at https://github.com/trenchproject/TrenchRmanuscript/.(DOCX)

Fig 1 .
Fig 1.A) Microclimate models scale air and soil temperatures and wind speeds from sensor to organismal heights.B) Biophysical models balance heat exchanges between organisms and their environment to estimate body temperatures[19].A thermal image depicts how grasshopper body temperatures vary substantially from vegetation, air, and ground temperatures due to heat exchange.

Fig 3 .
Fig 3. Body temperatures (Te) are predicted to drastically exceed air temperature when lizards are exposed to high levels of solar radiation.Air temperatures (Ta, ˚C) at lizard height (0.02 m) are predicted to exceed air temperatures at 2 m and to be similar to surface temperatures (Ts).We estimate body temperatures using two general energy budgets [solid: Tb_Gates(); dotted: Tb_CampbellNorman()] and a lizard specific biophysical model [dashed: Tb_lizard()] that differ in how they model heat exchanges.https://doi.org/10.1371/journal.pclm.0000139.g003