Microclimate complexity in temperate grasslands: implications for conservation and management under climate change

As climate change advances, there is a need to examine climate conditions at scales that are ecologically relevant to species. While microclimates in forested systems have been extensively studied, microclimates in grasslands have received little attention despite the climate vulnerability of this endangered biome. We employed a novel combination of iButton temperature and humidity measurements, fine-scale spatial observations of vegetation and topography collected by unpiloted aircraft system, and gridded mesoclimate products to model microclimate anomalies in temperate grasslands. We found that grasslands harbored diverse microclimates and that primary productivity (as represented by normalized difference vegetation index), canopy height, and topography were strong spatial drivers of these anomalies. Microclimate heterogeneity is likely of ecological importance to grassland organisms seeking out climate change refugia, and thus there is a need to consider microclimate complexity in the management and conservation of grassland biodiversity.


Introduction
Climate and climate change are heterogeneous across multiple spatial scales (Loarie et al 2009, Ackerly et al 2010. At fine spatial scales (typically <100 m), microclimates may occur where variation in the surface environment creates decoupling from broader conditions (Geiger et al 2009, Bramer et al 2018. Microclimates are ecologically significant because they may best reflect the conditions experienced by species in their habitats (Potter et al 2013, Hannah et al 2014. Such microclimate variation can have significant ecological consequences influencing the behavior, habitat use, and demographics of species (Grisham et al 2016, Frey et al 2016a. Further, microclimates may potentially act as microrefugia, buffering species from rapid climate change by allowing them to exploit environmental heterogeneity to remain within their climatic niche (Suggitt et al 2018, Kim et al 2022.
Measuring the influence of climate at ecologically relevant scales is a critical step towards understanding the climate vulnerability (or resiliency) of species, but this task is complicated by a mismatch between available climate products and the scale of environmental heterogeneity experienced by organisms (Potter et al 2013, Nadeau et al 2017, Lembrechts et al 2019. This mismatch can be multidimensional, as microclimate conditions may vary both vertically (due to differing biophysical processes near the Earth's surface) and horizontally (due to fine-scale variation in vegetation and surface topography) (Geiger et al 2009, Bramer et al 2018. An increasing number of microclimate studies have been conducted in forests revealing the ecological importance of temperature buffering provided by horizontal canopy structure (Frey et al 2016b, Zellweger et al 2020. In contrast, few detailed studies of microclimate variation have been conducted in open grasslands (but see Bennie et al 2008). This may be because grasslands are perceived as homogenous and topographically simple systems or because grasslands appear to have limited capacity to buffer climate relative to forests (Loarie et al 2009, Suggitt et al 2011. Yet this limited buffering capacity suggests that grassland-dependent species may be more vulnerable to climate change (Jarzyna et al 2016), emphasizing the importance of understanding microclimate variation within grassland systems.
Grasslands are among the most threatened and least protected ecosystems on earth (Scholtz and Twidwell 2022) and provide habitat for a diversity of climate-sensitive species, including pollinators (Hanberry et al 2021) and declining grassland birds (Rosenberg et al 2019). Grassland ecosystems are projected to experience a high velocity of climate change relative to other habitats because of their open structure and low elevation occurrence (Loarie et al 2009, Dobrowski et al 2013. Yet, despite lacking a canopy structure or large elevational gradients, grasslands may still experience spatial variation in microclimate conditions for several reasons. First, the near-surface environment experiences greater climate variability because of heat storage in the ground, reduced wind velocity, and less efficient heat-transfer and mixing processes (Geiger et al 2009, Bramer et al 2018. Second, variation in microtopography may impart differences in microclimate by influencing topographic shading and the drainage of cold air and surface run-off (Bennie et al 2008, Pastore et al 2022. Finally, variation in vegetation structure may also affect temperature and humidity conditions. For example, although grassland vegetation is narrowleafed and primarily vertically oriented, this foliage may provide shading at lower solar angles throughout the day , and importantly, dense vegetation also influences temperature and humidity conditions via evapotranspiration (Bramer et al 2018). Wooded edges can also shade adjacent grasslands-reflecting incoming shortwave radiation during the day and trapping outgoing thermal emittance at night-potentially creating a temperature gradient with proximity to these edges (Latimer and Zuckerberg 2017).
Such variation in grassland microclimate may affect the ecology of grassland-dependent species. For example, there is evidence that microclimates may influence grassland bird species distributions (Jähnig et al 2020), adult survival (Pérez-Ordoñez et al 2022), and nest success and productivity (Lloyd andMartin 2004, Carroll et al 2018). Similarly, grassland insects, such as butterflies, also appear sensitive to microclimate and may select habitat and oviposition in response to these conditions (Scherer and Fartmann 2022). Both temperature and humidity can impose physiological limits on organisms mediating processes like desiccation, heat stress mortality, egg unviability, and behavioral trade-offs (McKechnie and Wolf 2010, van de Ven et al 2019, Hoffmann et al 2021. Thus, microclimates could be of particular importance to species living in exposed environments where organisms may often operate near tolerance limits (Carroll et al 2016).
However, an obstacle to quantifying grassland microclimates is a lack of fine-scale remote sensing data to characterize such environmental variation. Most widely available satellite sources of land cover and habitat data provide imagery at resolutions of 30 m or much greater (Bramer et al 2018), and 1 m data available from some commercial systems may not be sufficient to capture relevant variation in space or time. In open grasslands, environmental heterogeneity may be more fine-grained relative to forested and montane environments (Zellweger et al 2019). Unpiloted aircraft systems (UASs) have recently emerged as a powerful tool to address such challenges. Specifically, UAS are capable of collecting imagery at sub-meter resolutions describing environmental characteristics relevant to microclimates, such as elevation, microtopography, primary productivity, and canopy height .
Our study had three objectives: (1) assess variation and magnitude of near-surface anomalies in temperature and vapor pressure in temperate grasslands, (2) understand the drivers of this microclimate variation, and (3) create spatially explicit models of grassland microclimates suitable for use in future ecological studies. We modeled microclimate using anomalies in near-surface iButton measurements of temperature and vapor pressure as responses predicted by a combination of fine-scale UAS environmental variables and gridded climate products (to characterize background conditions). Our approach provides what, to our knowledge, are the first fineresolution, spatiotemporal models of microclimate variation in temperate grasslands.

Experimental design and study area
We collected fine-scale, near-surface climate data, and environmental covariates (figure 1) from 15 May to 30 July 2021, at four grasslands located in Dane and Iowa counties in southern Wisconsin, USA (42 • 55 ′ 11.45 ′′ , −89 • 50 ′ 5.91 ′′ ). We focused on microclimate during the summer months because this period encompasses the critical reproductive phase for many grassland-dependent organisms. Our study sites ranged in size from 8.7-11.0 ha and were situated at similar mean elevations of 296-320 m, though the total elevational range between the highest and lowest points among all sites was 43 m. Our study sites were in Wisconsin's historically unglaciated Driftless Area (figure S1) and characterized by rolling topography and a mosaic of agricultural land, forest edge, drainages, riparian zones, planted grasslands, and prairie remnant. The four sites were on a mixture of public, private, and non-governmental organization-owned lands. All sites were considered planted grasslands, but two were non-native, cool-season monocultures, and two diverse warm-season sites that resembled native Tallgrass prairie (figure S2). Cool-season sites were dominated primarily by Smooth Brome (Bromus inermis) while warm-season sites were dominated by Big Bluestem (Andropogon gerardi), Indian Grass (Sorghastrum nutans), and a diversity of native forbs. The sites varied with respect to topographic complexity and level of shrub encroachment. We deliberately chose sites that would capture a range of conditions in midwestern grasslands.

iButton data collection
We systematically deployed iButtons (DS1923, Maxim Integrated, San Jose CA, USA) at a 5 cm height and at a 55 m grid resolution across each of four sites to capture the range of environmental variation within grassland sites; the number of iButtons per site ranged from 23 to 30 for a total of 109 iButtons (figures S3-S6). Five units across all sites failed to collect data, leaving a total sample size of 104. It has been demonstrated that unshielded iButtons in high irradiance, low-wind conditions may be biased toward warmer temperatures by as much as 15 • C . To help mitigate this potential bias, we housed iButtons in ventilated PVC capsules coated with reflective foil tape (figure S7). PVC capsules have been shown to reduce error by approximately 50% , and the addition of foil tape to deflect direct solar irradiance likely further reduced this rate. We programmed iButtons to log hourly temperature and relative humidity measurements for the duration of the study period. Upon retrieving data from iButtons, we converted relative humidity to absolute vapor pressure (i.e. partial pressure of water vapor) using Teton's formula (Norman and Campbell 1998). Because the primary objective of our study was to identify relative variation in grassland microclimates (e.g. 'hot spots' and 'cool spots' , 'humid spots' and 'dry spots') rather than attempt to model exact temperature or humidity mechanistically, we calculated anomaly by subtracting hourly measurements for each unit from the mean of all other iButtons at the same site for each logging interval. We then summarized hourly anomaly by day for each iButton by averaging all hourly anomalies for a given day, thus producing daily mean anomalies of temperature and vapor pressure for each measurement location.

Microclimate predictors
We used a combination of UAS-collected light detection and ranging (LiDAR) and multispectral imagery to derive seven spatially explicit layers of vegetation and topography for each site summarized at 60 cm resolution (figures S3-S6). We collected imagery using a Matrice M210 V2 RTK drone (DJI, Shenzhen, China) equipped with two sensors: (1) Sentera 6× (Sentera, Minneapolis MN, USA) five-band multispectral camera with 5.4 cm resolution, and (2) VLP-16 LiDAR sensor (Velodyne Lidar, San Jose CA, USA) with a point density of 190-300 m −2 . The multispectral camera had an additional downwelling incident light sensor that measured incoming radiation and was used to radiometrically correct images to account for changing illumination conditions throughout the flights. To further reduce the bias of illumination conditions on collected imagery, we conducted flights during the hours of 1000-1500 CST in clear and calm conditions and flew at an altitude of 120 m above ground level (AGL). However, we conducted LiDAR flights at a lower altitude of 80 m AGL to increase point density. We collected all imagery at a speed of 8-10 m s −1 with 80% photo overlap for the multispectral camera and with cameras in a nadir position. We conducted three flights at each site throughout the season-once in May, June, and July-to capture seasonal variation in vegetation predictors.
We processed multispectral imagery using structure from motion techniques with Metashape Pro (Version 1.6.5, Agisoft LLC, St. Petersburg, Russia) to produce orthorectified rasters. We processed LiDAR point clouds in R (Version 4.2.1; R Core Team 2022) using the LidR package (Roussel et al 2020). We used a digital terrain model (DTM; Zellweger et al 2019, Duffy et al 2021) created from LiDAR point clouds to characterize bare earth elevation at each of our sites, and we then used these data to derive layers for topographic positioning index (TPI)-where high index values are associated with higher elevation relative to the neighborhood-slope, aspect, and hill shade with the terra package in R (Hijman 2022). TPI was created using a moving window of ∼20 m in each cardinal direction within which the elevation of each 60 cm focal cell was subtracted from the neighborhood mean. Thus, our TPI product was calculated at a 20 m scale but ultimately yielded a 60 cm resolution layer. There is no commonly accepted standard for the most appropriate window size for TPI in microclimate studies, therefore we selected our neighborhood size based on how well it appeared to capture the microtopography of the site according to our field observations. We found that, at our sites, larger neighborhoods tended to overlook smaller depressions and created a variable that was more similar to simple elevation. Hill shade was calculated using values for solar noon at the mid-point of the study period to best represent the potential influence of terrainadjusted solar irradiance. Photogrammetry terrain models can increase error when bare ground is heavily obstructed because imagery cannot penetrate canopy cover (Klápště et al 2020). To address this issue, we used a hybrid approach to create a canopy height model (CHM) by differencing the multispectralderived digital surface models of the canopy from the LiDAR-collected DTM representing bare earth, thus creating a raster layer representing vegetation height (Zellweger et al 2019). Finally, we calculated normalized difference vegetation index (NDVI) using the Sentera red (670 nm) and near-infrared (NIR; 870 nm) channels [NIR − red]/[NIR + red] as an index of primary productivity (e.g. Duffy et al 2021). NDVI is generally considered a proxy for greenness identifying areas of live, photosynthetically active vegetation. In our study area, greening of cool-season plant species may begin in March or April following snowmelt, while warm-season species often begin in late May or June. During our study season, NDVI at cool-season monoculture sites appeared to reach a peak in June, while diverse warm-season sites showed a more diffuse pattern of greening (figures S8 and S9). For NDVI and CHM (e.g. vegetation height)-which were somewhat dynamic throughout the seasonwe associated iButton measurements with the spectral capture from the appropriate month (e.g. iButton measurements from May were associated with NDVI and CHM values from the May imagery collection) to make our models as realistic as possible (figures S8-S11).
In addition to UAS spatial predictors, we created rasters of distance to nearest wooded edge by manual digitization in ArcGIS (Version 10.5.1; ESRI, Redlands CA, USA). We also included a site-level variable, ecotype, to account for potential differences between warm and cool-season sites. Finally, to account for the influence of broader climate context on microclimate anomalies (Wolf et al 2021), we included mesoclimate predictors (horizontal scale of 1-300 km; Bramer et al 2018) in our models. This also allowed us to examine We sourced mesoclimate predictors from Daymet, a daily, 1 km resolution gridded climate product available continent-wide across North America (Thornton et al 2021). Daymet uses daily weather station data, and data describing terrain, large bodies of water, winds, and storms to interpolate local climate conditions; these predictions are then crossvalidated (Thornton et al 2021). Table 1 summarizes all predictor variables included in our microclimate models including UAS-collected predictors, site-level variables, and Daymet mesoclimate variables.

Random Forests models
We used Random Forests to model microclimate at grassland sites using iButton collected temperature and vapor pressure anomalies as response variables predicted by a set of spatially explicit vegetation and topographic features, as well as gridded mesoclimate variables (table 1). Although we included mesoclimate predictors in our models, our approach was not intended as a downscaling, but rather to model our iButton response variables under varying external conditions. Random Forests are a decision tree-based, machine learning algorithm capable of achieving high predictive accuracy without overfitting through ensemble modeling and a bootstrapping technique known as bagging (Breiman 2001, James et al 2013. Random Forests have several additional advantages including the ability to handle highdimensional data, model complex interactions, and a lack of restrictive assumptions. We calculated Pearson's correlation coefficient (r) among all pairs of spatial predictors and found r < 0.7 in all cases (figure S12). We implemented Random Forests in R with the ranger package (Wright and Ziegler 2017) and performed k-fold cross-validation with ten partitions, in which ground-collected iButton measurements were withheld for testing. These control data were then compared with predictions from the trained model to calculate root-mean-squared errors as a metric of model performance. We performed cross-validation and model tuning using the caret package (Kuhn 2022). We tuned our Random Forest models across a range of mtry values, as well as two different split rules, using 1501 trees for each iteration (table S1). The mtry value determines how many predictors are included in each bagging and thus represents a balance between good prediction and overfitting. We provided a range of low and high mtry values relative to the number of predictors in our models and selected the optimal value from the tuning grid models.
We extracted variable importance, individual conditional expectations (ICEs), and partial dependence for each model using the flashlight package (Mayer 2021). Variable importance is calculated by assessing the drop in model performance when a feature is randomly permuted, thus a low variable importance score indicates a lack of influence on model performance (Breiman 2001). In our results, we present variable importance as a percentage of total drop in performance from all features combined. This allows variable importance to be compared across models with different units in the response variable. Partial dependence is a global metric that examines the average effect of a feature on the response variable when all other predictors are held constant (Friedman 2001) and can be interpreted analogously to a marginal effects plot. ICE plots represent the local observations from which partial dependence is produced, and they are suitable for identifying interactions, non-linearity, and other complex effects that may be obscured when observations are averaged as in partial dependence (Molnar 2022). We assessed interaction between features in our models using Friedman's H, also implemented in the flashlight package (Mayer 2021). Friedman's H is calculated from the decomposition of partial dependence for a given feature and can be interpreted as the proportion of variance explained by that feature that is attributed to covariance with another feature (Molnar 2022). Finally, we tested our models for spatial autocorrelation, and those procedures are described in supplementary materials, appendix A.

Spatial prediction of grassland microclimates
We produced spatially explicit, 60 cm resolution, microclimate predictions for both temperature and vapor pressure anomalies at our four grassland study sites by using Random Forests models to predict across raster stacks containing each spatial predictor. For NDVI and CHM, we used mean layers for prediction to represent typical grassland conditions. To examine how microclimates changed under different climate conditions, we produced three sets of spatial predictions for each study site representing cool and cloudy days, mean conditions for the season, and hot and sunny days. We characterized these conditions, respectively, by setting Daymet temperature and solar irradiance variables to their lower 5th percentile (5%), mean, and upper 95th percentile (95%) values.

Results
Hourly iButton summaries revealed similar microclimate profiles among study sites, suggesting few broad differences in climate dynamics between warm-and cool-season grasslands (figure S13). Across all sites, daily mean temperatures and vapor pressures ranged from 5.4 • C-30.7 • C and 781.1-3168.3 Pa, respectively (figure S14). Hourly profiles of temperature and vapor pressure anomalies exhibited deviations from mean conditions of up to 10 • C and 3000 Pa, respectively, during daytime hours (figures S15 and S16). Daily summaries of raw iButton data revealed that, on average, grassland microclimates at our sites tracked mesoclimate conditions described by Daymet variables (figure 2(a)), but experienced lower minima (figure 2(b)), notably higher maxima (figure 2(c)), and ultimately were subject to larger diurnal temperature ranges-particularly on warmer days ( figure 2(d)).
Random Forests models of daily mean temperature (R 2 = 0.71) and vapor pressure anomalies (R 2 = 0.52) had predictive accuracies of 0.38 • C and 55.2 Pa, respectively (table S3). The top spatial predictors in both models included primary productivity, as represented by NDVI, canopy height (CHM), elevation, and relative topographic positioning (TPI); aspect also informed vapor pressure anomaly (figure 3). In addition to spatial UAS predictors, Daymet daily mean temperature and solar irradiance were also modestly informative in our models (figure 3), suggesting that background climate conditions affected grassland microclimate dynamics.
Analysis of partial dependence showed that areas of high primary productivity (NDVI), and to a lesser extent taller vegetation (CHM), supported cooler anomalies (figures 4(a) and (b)). By contrast, higher elevation, as well as higher index values for relative topographic positioning (TPI), were moderately associated with warmer temperature anomalies (figures 4(c) and (d)). Similarly, areas of more productive and taller vegetation were associated with higher vapor pressure anomalies (more humid; figures 4(a) and (b)), while higher elevation positions were much drier (figures 4(c) and (d)). Partial dependence plots for all variables modeled can be found in the supplementary materials, figures S17 and S18. Daily mean temperature and solar irradiance from Daymet characterizing broader mesoclimate conditions had non-unidirectional effects (figures 5(a) and (b)). Specifically, warm mesoclimate temperatures and high irradiance appeared to increase both negative and positive anomalies in temperature and vapor pressure, effectively increasing the magnitude of microclimate decoupling produced by variation in the spatial predictors discussed above. This relationship was confirmed by plots of absolute iButton anomalies against Daymet mesoclimate conditions throughout the season (figure S19). We also assessed pairwise interactions between all predictors in both of our models and found that only the interaction between elevation and vegetation height (CHM) in the temperature anomaly model appeared to be substantive (figure S20).
Finally, spatially explicit predictions of microclimate revealed heterogeneous conditions at grassland sites reflecting underlying variation in vegetation and topography (figure 6; spatial predictions of microclimate for all sites can be viewed in supplementary materials, figures S21 and S22). However, the magnitude of anomalies became more pronounced under hot and sunny conditions (figure 6; 95th percentile). Under these conditions, the range of daily temperature anomalies between the coolest and warmest areas was as much as 3.5 • C at some sites ( figure S23). Similarly, differences in daily vapor pressure anomalies between the most humid and dry areas at some sites were as much as 250 Pa (figure S24). Model error structures showed little directional bias in prediction and were broadly similar among sites although differences in distribution tails suggested some differences in predictive accuracy (figures S25 and S26). Additionally, at several sites, there was evidence that extreme values were predicted with greater uncertainty (figures S27 and S28). Raw iButton microclimate measurements summarized by day of season and averaged across all sites plotted alongside Daymet observations for the same period. While mean microclimate temperatures closely tracked regional conditions, maxima and minima exceeded Daymet estimates suggesting that grassland microclimates experience a larger range of temperature conditions than regional climate.

Discussion
Grasslands are a globally endangered ecosystem (Scholtz and Twidwell 2022) that is likely to experience rapid climate change (Loarie et al 2009). Consequently, many grassland-dependent species are vulnerable to climate conditions (González-Varo et al 2013, Wilsey et al 2019, and microclimates may play an important role in buffering these species from future changes (e.g. Suggitt et al 2018). However, a lack of studies examining the magnitude and drivers of grassland microclimate has limited understanding of this potential. We found that grasslands support substantial variation in microclimate anomalies, and in some cases, the magnitude of temperature variation (3.5 • C mean; 10 • C hourly) may be of a similar order to projected climate change for the central United States (Pörtner et al 2022;3 and 6 • C warmer under 2 • C and 4 • C average global warming scenarios, respectively). The spatial distribution of grassland microclimates was driven largely by vegetation conditions, including primary productivity and vegetation height, as well as elevation and topographic positioning. Microclimate differences also became more pronounced in hot and sunny conditions, suggesting that microclimates could play an important ecological role as microrefugia during extreme events (e.g. heat waves).
Broadly, microclimate temperatures recorded by iButtons in our study documented greater climate variability near the surface in grasslands, as well as substantially greater maxima, and to a lesser extent, lower minima relative to regional conditions described by Daymet. This was expected given that there is greater climate variability near the surface. In addition, our model residuals indicated that in some cases, the most extreme anomalies were predicted with greater uncertainty and thus these estimates should be interpreted cautiously. Much higher maximum temperatures recorded in grasslands could also reflect bias associated with iButton use in the absence of canopy , or a lack of temperature buffering in grasslands relative to other cover types (Suggit et al 2011). It is likely that both mechanisms were at play suggesting that while some maxima we observed may have been exaggerated, grasslands do experience greater exposure to temperature extremes as reported elsewhere (Loarie et al 2009, Suggitt et al 2011. Temperature anomalies in our study were most influenced by primary productivity, and productive areas were correlated with cooler microclimates suggesting that dense vegetation may buffer heat extremes in grassland habitats. This effect was likely driven by evaporative cooling through plant stomatal conductance associated with photosynthetic capacity (Bramer et al 2018. Vegetation height also had a moderate cooling effect-likely from shading during lower solar angles. The importance of vegetation productivity and structure suggests that microclimate buffering of temperature in grassland systems may be threatened by extreme droughts and changes to precipitation patterns associated with future climate change (Cook et al 2022); prolonged drought conditions may cause browning of vegetation and an associated reduction in evapotranspiration, as well as diminished plant structure. However, an important limitation of our study is that climate observations were made at a 5 cm height, typically below the grass canopy. Therefore, microclimate dynamics above or near the top of vegetation could operate differently.
Both elevation and topographic position also influenced microclimate, and low-lying areas relative to neighborhoods, as well as lower elevations in general, were associated with cooler and more humid anomalies-an effect possibly driven by the pooling of cold air and surface water run-off, as well as nighttime humidity inversion (Geiger et al 2009, Bramer et al 2018, Pastore et al 2022. Although cold air drainage is more often associated with mesoclimatic processes in mountainous systems (Dobrowski 2011, Ashcroft andGollan 2012), even depressions with depths as little as 2 m may be sufficient to influence temperature conditions and cold air movement (Mahrt 2022). Another possibility is that lower and less exposed positions tend to receive more terrain shading and reduced wind speeds, thus preserving moisture and leading to greater evaporative cooling; this might explain the interaction of elevation and vegetation height observed in our temperature model. Regardless, our findings suggest that although grassland ecosystems are characterized as low-lying and lacking in topographic complexity (e.g. compared to montane forests), small microtopographic differences can support strong microclimate variability (Mahrt 2006(Mahrt , 2022. Given that vegetation played an important role in mediating grassland microclimates, lack of differences between cool-and warm-season grassland sites was surprising and suggests that site productivity and vegetation structure, and not species composition or functional traits (Zellweger et al 2019), have a greater influence on microclimate conditions. However, we did not collect primary productivity data at a fine temporal grain (e.g. daily measurements), and it is conceivable that differences in C3 and C4 carbon pathways in cool and warm-season plants (Wang et al 2013) may influence microclimate differently throughout the season depending on the timing of greening. Plant functional diversity may also play an important role in supporting the resiliency of grasslands to extreme climate events (Craine et al 2013), thus further study on the subject could be of value.
The magnitude of temperature and vapor pressure anomalies was also influenced by daily temperature and solar irradiance, which drove increasing anomalies in both positive and negative directions (e.g. warmer and cooler, more and less humid). Hot and sunny days supported greater microclimate anomalies, suggesting heat extremes impart greater contrast across grassland landscapes. This effect is likely driven by diurnal temperature and vapor cycles in which extremes occur at the surface where heat stored in the ground is transferred to the near-surface air layer via convection (Geiger et al 2009). On hot and sunny days where there is high energy input into the system, these extremes are likely to be greater, thus creating high anomalies in exposed areas and low anomalies where conditions are buffered by other factors. The ultimate result of these biophysical processes is greater microclimate deviation from mean conditions with heat and irradiance extremes-a phenomenon that is of high ecological importance to grassland species seeking refugia to avoid exposure to conditions beyond their thermal and moisture tolerance limits (Grisham et al 2016, Ruth et al 2020, Scherer and Fartmann 2022. Many grassland-dependent species are experiencing widespread declines and are vulnerable to climate change (González-Varo et al 2013, Wilsey et al 2019). Fine-scale temperature and humidity conditions are important to many of these grassland organisms and can influence their reproductive success, and use of habitat across landscapes. For example, grassland birds may select nest sites that buffer eggs and nestlings from lethal thermal extremes and desiccation (e.g. Grisham et al 2016, Carroll et al 2018. Landscape-level thermal conditions, and the presence of refugia, can also influence the habitat use and survival of these species (Hovick et al 2014, Carroll et al 2016, Ruth et al 2020. Similarly, oviposition and habitat use in grassland insects and pollinators is also influenced by microclimates (Gardiner andHassall 2009, Scherer andFartmann 2022). Consequently, microclimates may become an increasingly important aspect of managing and conserving grassland biodiversity under climate change.
Grasslands are often subject to intensive management actions such as prescribed fire, grazing, and mowing. Understanding how these tools influence microclimates is of critical importance. For example, given the role of primary productivity in fostering cooler microclimates, managers may wish to consider activities that could enhance productivity such as sustainable grazing regimes (to avoid overgrazing) or irrigation to improve primary productivity during extreme droughts (Greenwood et al 2016). By contrast, intensive activities such as frequent mowing may be less desirable for creating microrefugia during the breeding season as they reduce primary productivity and homogenize structure. In some cases, these activities have even been implicated in disrupting microclimates required by insect species (Gardiner andHassall 2009, Thomas et al 2009). Finally, the importance of topography in creating grassland microclimate complexity suggests that prioritizing the acquisition and conservation of grassland areas of high topographic complexity will likely increase the prevalence of microrefugia for climate-vulnerable grassland species (Bennie et al 2008, Suggitt et al 2018.

Data availability statement
The data that support the findings of this study are available upon reasonable request from the authors.

Acknowledgments
We thank the University of Wisconsin-Madison and the Wisconsin Cooperative Wildlife Research Unit for supporting this work. We thank Wisconsin DNR, The Prairie Enthusiast, The Nature Conservancy of Wisconsin, and D and W Marshall for property access. We thank the technicians who helped collect these data including S Faterioun, F Hill, S Klein, B Sellers, and P Thompson. Finally, we thank K McGinn for design of figure 1 and M Griffith for assistance with imagery processing.

Funding statement
Funding was provided by National Institute of Food and Agriculture, United States Department of Agriculture, Hatch Project 1020176, the Wisconsin Cooperative Wildlife Research Unit, and NSF 'ASCEND' Biology Integration Institute (BII), supported by DBI Award 2121898.