Intraspecific variation in thermal tolerance differs between tropical and temperate fishes

How ectothermic animals will cope with global warming is a critical determinant of the ecological impacts of climate change. There has been extensive study of upper thermal tolerance limits among fish species but how intraspecific variation in tolerance may be affected by habitat characteristics and evolutionary history has not been considered. Intraspecific variation is a primary determinant of species vulnerability to climate change, with implications for global patterns of impacts of ongoing warming. Using published critical thermal maximum (CTmax) data on 203 fish species, we found that intraspecific variation in upper thermal tolerance varies according to a species’ latitude and evolutionary history. Overall, tropical species show a lower intraspecific variation in thermal tolerance than temperate species. Notably, freshwater tropical species have a lower variation in tolerance than freshwater temperate species, which implies increased vulnerability to impacts of thermal stress. The extent of variation in CTmax among fish species has a strong phylogenetic signal, which may indicate a constraint on evolvability to rising temperatures in tropical fishes. That is, in addition to living closer to their upper thermal limits, tropical species may have higher sensitivity and lower adaptability to global warming compared to temperate counterparts. This is evidence that freshwater tropical fish communities, worldwide, are especially vulnerable to ongoing climate change.

The capacity of ectothermic species to cope with ongoing global warming, especially the increasing frequency, intensity and duration of extreme heatwaves, will be influenced by their upper thermal tolerance limits [1][2][3] . Tolerance of acute warming, measured as the critical thermal maximum (CT max ), varies among fish species according to thermal conditions in their habitat 4 . Tropical species live in warm, relatively thermally stable habitats; they have narrow thermal tolerance ranges but higher CT max than species at temperate latitudes. Their warm habitat temperatures are also, however, closer to their limits of upper thermal tolerance, so they have a limited thermal safety margin (defined as the difference between upper thermal tolerance limit CT max of adult life stage and the maximum habitat temperature during summer 5 ) and consequently are considered to be especially vulnerable to global warming [6][7][8][9] . Temperate species have lower absolute thresholds for tolerance of warming, but they have broader tolerance ranges, presumably because they encounter a wide range of habitat temperatures, both seasonally and spatially. This is linked to wider thermal safety margins than in tropical species 4,10 . These patterns of vulnerability to global warming among species at a geographic scale are major issues in projecting impacts of warming. They have a strong phylogenetic basis, which is believed to reflect local adaptation to common ancestral thermal regimes in related species 11 .
Studies of broadscale geographic patterns in vulnerability have, to date, focused upon average values for CT max among fish species. The significance of intraspecific variation in tolerance remains to be explored. The extent of variation in functional traits within species, particularly of physiological tolerances (e.g. CT max , hypoxia tolerance, pollutant resilience, immune resistance) is expected to have a profound influence on their vulnerability to global change [12][13][14][15] . Possessing a broad range of tolerance phenotypes in populations can reduce sensitivity to impacts of environmental stressors, through various proximate ecological mechanisms [12][13][14] . If phenotypic variation is linked to underlying genetic diversity in the species, this can provide scope for adaptability and evolvability, by yielding genotypes for selection in changing environments [12][13][14] . When fish species are challenged by thermal stressors, www.nature.com/scientificreports/ such as increased seasonal temperatures and extreme heatwaves, the population sensitivity and adaptability will be major determinants of their relative vulnerability [13][14][15][16] (Fig. 1). Fish species show intraspecific variation in CT max , which has a component of both phenotypic plasticity and heritable genetic variation 15,[17][18][19] . The CT max varies among populations of fish species, due to local adaptation [20][21][22] , indicating that the trait evolves in response to prevailing thermal regimes. Given the broader thermal range experienced by temperate fish species, within generations and over evolutionary time, we hypothesized that they would exhibit greater intraspecific variation in their thermal tolerance, measured as CT max , than tropical species. We predicted that lower variation in CT max within species might be observed in species with low thermal safety margin, because a small margin might constrain scope to express variation 10 . We also expected the extent of variation in CT max to have a phylogenetic basis, indicating that it reflected evolutionary processes of adaptation.
We used published data 4 and, after a data selection process (see "Methods"), we estimated the extent of intraspecific variation in CT max of 203 species of ray-finned (actinopterygian) fish (n = 127 freshwater, n = 76 marine), based on the standard deviation of the mean. We were well aware that the selected studies in the dataset did not have the same protocol procedures. They did not use the same heating rate (0.0017-1 °C/min) and fish size, both of which can influence CT max and standard deviation of the mean. We choose to not include these variables in our main analysis because of the high variation of heating rate used and for fish size there was insufficient reporting for this data among studies. We performed a supplementary analysis with heating rate in the model on 186 species. In our main analysis we then compared two latitudinal groups, temperate to tropical species, considering the boundary to be 23° latitude. We also evaluated if variation in CT max depended on whether species were from northern or southern hemisphere or whether species were marine or freshwater and their individual CT max . Finally, we used the magnitude of the difference between acclimation temperature (T a ) and CT max , which we denoted delta temperature (ΔT = CT max − T a ), as an indication of the capacity to increase CT max depending on the acclimation temperature, and evaluated if it was linked to intraspecific variation in CT max . All of the results were based on a phylogenetically informed analysis (phylogenetic least squares regression, PGLS, see "Methods"), to establish how patterns in the extent of variation were linked to evolutionary thermal history of the species.

Discussion
Overall, tropical species show a lower intraspecific variation in thermal tolerance than temperate species. Specifically, freshwater tropical species have reduced within-species variation in thermal tolerance compared to freshwater temperate species. Conversely, marine temperate species display lower intraspecific variability in CT max compared to marine tropical species, although with a lower significant level (p = 0.036) than that found when comparing freshwater temperate vs tropical species (p = 0.005). To better understand the difference in effect direction in marine vs freshwater environments, further investigation is required on the driving factors modulating this difference. Nevertheless, if comparatively low intraspecific variability in CT max reflect a reduced capacity for phenotypic plasticity, this will increase sensitivity to warming in the short term. If intraspecific variability in CT max reflects diminished heritable genetic variation, a low value implies decreased adaptability and evolvability to a warmer and more thermally stressful future, over generational timescales. The lower intraspecific variation in CT max in freshwater tropical as compared to freshwater temperate species (Fig. 2) renders the former especially vulnerable to future warming, in particular to extreme events 23,24 (Fig. 1). This will negatively affect the vulnerability of freshwater tropical species living near their upper thermal limits 4,6,7,25 . Interspecific variability in thermal tolerance tends to be higher in freshwater than in marine www.nature.com/scientificreports/ species, particularly in temperate areas of the northern Hemisphere 4 . This is likely because freshwater species are distributed across limited latitudinal ranges, while marine species have wider latitudinal ranges, giving rise to a relatively invariant thermal tolerances at the faunal level 4 . This difference in variability of thermal tolerance among species in freshwater vs marine habitats may reflect in a higher variability within species in temperate freshwater species. Therefore, local thermal conditions experienced by species are determinant in setting the natural individual variation within populations. The fact that variation in thermal tolerance was greater in the northern compared to southern hemisphere could be the result of two phenomena: (1) greater thermal variability in the northern hemisphere 4,6 ; or (2) a relative paucity of data for the southern hemisphere 26 . Nevertheless, there was no effect of hemisphere on intraspecific variation in CT max .
We found a significant interaction effect between latitudinal group and ΔT. Tropical species with high ΔT showed a larger intraspecific variation in CT max . Lower acclimation temperatures allowed to set an extended range of variation in CT max . We suggest that low acclimation temperatures provided a certain thermal plasticity and allowed scope for thermal variation.
CT max was negatively linked to S.D. CT max in fishes. This reveals a ceiling to thermal plasticity capacity in fishes. This might explain why tropical species show lower S.D. CT max as they have higher CT max than temperate species. At the highest upper thermal tolerance limits, fishes are not able to express a large range of thermal resilience variation within species 15 .
In this study, we decided to choose the low or mid-range acclimation temperatures to collect the standard deviation of CT max . The reason why we did not take the highest acclimation temperatures tested over the studies, is that they can be different according to if the species is temperate or tropical. They can be higher for tropical species and display different results compared to temperate species. We chose not to standardize the SD CT max values to a common acclimation temperature, because we were more interested in using the real values from the studies.
The strong phylogenetic signal for the extent of intraspecific variation in CT max is presumably because many families contain species with a relatively common history of thermal adaptation (Fig. 4). That is, they have occupied similar thermal regimes within temperate or tropical habitats. In particular, there is a latitudinal effect on family distributions, with some families only being present in temperate (e.g. Gadidae) or tropical (e.g. Apogonidae) habitats, although some cosmopolitan families have species in both (e.g. Gobiidea, Blennidae) ( Figure S2). In addition to the geographic collinearity that may be occurring with some families, the phylogenetically based differences in intraspecific variation among species may cause evolutionary constraints on evolvability in the face of ongoing warming and exposure to extreme events in freshwaters. The extent of such constraints is not clear and would depend on the exact genes affecting thermal tolerance and how these are represented within each family 11 . Further highlighting how temperature regime may shape evolutionary trajectories within closely related species or those with a common ancestor, with potential consequences for their vulnerability to thermal stress [27][28][29][30] .
This evidence for higher vulnerability of tropical species to climate variability and extreme warming events 31 may have numerous ecological implications beyond simple tolerance thresholds. Freshwater tropical species may be obliged to seek thermal refugia in colder areas if these are available, potentially changing community structures 9,32 ; such distribution shifts could have major ecological consequences 33,34 . Overall, the extent of intraspecific variation in CT max must be considered in models that project impacts of warming on fishes. Intraspecific variation for tolerance in other environmental conditions such as hypoxia and acidification would be the next step for future research. Further research should focus on the mechanisms that underly latitudinal variation in CT max and whether these reflect universal principles across all species.

Methods
Dataset and data selection process. We used the data on CT max in marine, brackish and freshwater fish species (2722 observations unimputed data set) published by 4 . We performed a three-step selection procedure to identify the species for this study. First, we excluded data where CT max was measured using death as an endpoint (1256 observations) as these do not correspond to the accepted definition of CT max (loss of equilibrium but not death) 35 , so the temperatures recorded will have exceeded the critical threshold. Second, we excluded polar species because of the sample size (n = 5) and discarded brackish water species because no indication was given about the nature of the brackish habitat (e.g. lagoon, estuary or others). Third, several species were tested at different acclimation temperatures resulting in multiple CT max measures for the same species. We therefore took CT max values measured at the lowest or mid-point tested acclimation temperature with the largest sample size of individuals used. This data selection procedure produced a dataset of 203 fish species for which we have S.D. of their CT max (standard deviation).
Calculation of delta temperature. We calculated the ΔT The ΔT defines the distance from thermal acclimation (T a ) to thermal tolerance limit (T CTmax ), providing an index of vulnerability to acute heating 10 and of thermal acclimation capacity. In other studies, ΔT is defined as the difference between the highest experienced summer temperature and the CT max and referred to as thermal safety margins. In our study we use thermal acclimation temperature and decided to define ΔT based on the difference between CT max and thermal acclimation temperature. This accounts for the fact that acclimation temperature is often asymptotically linked to CT max 15 www.nature.com/scientificreports/ pleted by AIC values using the AIC function estimating the best model fit with the lowest AIC value (see Suppl. Table 1). The phylogeny of 203 fish species was found and generated from the comprehension tree of life (Fig. 4) 40 using the "rotl" package 41 . A measure of phylogenetic correlation, λ, the degree to which this trait evolution deviates from Brownian motion 42 , was evaluated by fitting PGLS models with different values of λ to find that which maximized the log-likelihood of the best-fitted model. The level of statistical significance was set at alpha = 0.05.

Phylogenetic analysis.
This was performed by PGLS on the 203 species' specific geographical location, habitat, ΔT, CT max and number of individuals measured. As fishes' physiology is dependent on the environmental thermal conditions, hemisphere was incorporated into the model because of the significant differences in thermal variability between the two hemispheres 6 , with the north having higher thermal variation than the south 26 . Due to the effects of local thermal variation on fish thermal physiology, we included an interaction term between latitudinal groups (tropical versus temperate) and the ΔT (suppl. Table 2). Two further interaction terms were included in the model between latitudinal groups with hemisphere and habitat type (freshwater versus marine). We also conducted general linear model (GLM) analysis to exclude the effect of phylogeny on the outcome of the observed variation in log 10 S.D.CT max , testing the individual effects of our variables in the model (suppl. Table 3) and comparing the outcoming results to PGLS analysis. As heating rate can also influence the intraspecific variation in CT max , we ran a supplementary PGLS analysis on 186 species including heating rate in the selected model (suppl. Table 4).