Evolvability under climate change: Bone development and shape plasticity are heritable and correspond with performance in Arctic charr (Salvelinus alpinus)

Environmental conditions can impact the development of phenotypes and in turn the performance of individuals. Climate change, therefore, provides a pressing need to extend our understanding of how temperature will influence phenotypic variation. To address this, we assessed the impact of increased temperatures on ecologically significant phenotypic traits in Arctic charr (Salvelinus alpinus). We raised Arctic charr at 5°C and 9°C to simulate a predicted climate change scenario and examined temperature‐induced variation in ossification, bone metabolism, skeletal morphology, and escape response. Fish reared at 9°C exhibited less cartilage and bone development at the same developmental stage, but also higher bone metabolism in localized regions. The higher temperature treatment also resulted in significant differences in craniofacial morphology, changes in the degree of variation, and fewer vertebrae. Both temperature regime and vertebral number affected escape response performance, with higher temperature leading to decreased latency. These findings demonstrate that climate change has the potential to impact development through multiple routes with the potential for plasticity and the release of cryptic genetic variation to have strong impacts on function through ecological performance and survival.

can alter the degree of variation available to selection through the release of cryptic genetic variation . Therefore, current issues surrounding environmental change are highly relevant to plasticity research as such developmental changes are likely to precede the shifts in demography or genetic variation that are often the focus of conservation efforts (Campbell et al., 2017;Skúlason et al., 2019).
To understand how plasticity may affect an organism's response to a changing environment, conservation biology should expand its focus beyond assessments of extant phenotypes and genetic variation (Campbell et al., 2017). To enable this, evolvability, the capacity of a population to adaptively evolve, could be adopted as a target for conservation biology to deal with changing environments (Campbell et al., 2017;Hendrikse et al., 2007). Evolvability extends conventional trait heritability to include both genetic variation in a population and how this variation interacts with the environment to produce the phenotype-the ultimate target of selection (Campbell et al., 2017;Mayr, 1997;Parsons & Albertson, 2013).
Climate change will rapidly impose major changes in the environment worldwide. However, temperature increases are likely to be greater at high latitudes (Cohen et al., 2014;Overland et al., 2013;Parmesan, 2007;Pithan & Mauritsen, 2014). Estimates of this "Arctic amplification" effect predict warming in the far north to be three to four times the rate of temperate and tropical regions (Bracegirdle & Stephenson, 2013;Overland et al., 2013;Pithan & Mauritsen, 2014). This suggests that thermal plasticity (i.e., within-generation phenotypic responses to temperature) could play a major role in how populations, especially those from high latitudes, cope. Ectotherms, such as fishes, in these areas are likely to be more susceptible to changes, with evidence indicating that there are already shifts in life-history traits in some fish species (Cai et al., 2014;Crozier & Hutchings, 2014).
Early-stage developmental processes are likely to play an important role in the persistence of populations in variable environments. For example, early developmental events, such as the condensation of bone, frequently occur at life stages that are also subject to strong selection (Johnson et al., 2014;McCormick et al., 2018;Pechenik, 2006;Sifa & Mathias, 1987). Bone development is critically important as it underlies variation in morphology and is a major determinant of adaptive traits such as swimming performance, foraging efficiency, and predator avoidance (Blake, 2004;McCormick et al., 2018;Nesteruk et al., 2014;Webb, 1978Webb, , 1984. Bone development is also highly sensitive to environmental variation partly because it exhibits continuous growth and resorption (mediated by osteoblasts and osteoclasts, respectively). The processes of bone remodeling are often linked to mechanical stress but for natural populations, it is largely unknown how these may also be affected by temperature (Totland et al., 2011;Witten & Huysseune, 2009).
Many adaptive traits are dependent on variations in the expression of bone tissue that affect ecological performance. At least, some of these traits are known to result from plasticity (Parsons & Robinson, 2007;Skúlason et al., 2019). A strong link between the effects of temperature on bone development and potential fitness effects has not been demonstrated. In fishes, temperature is well known to impact bone through the expression of meristic traits such as the number of fin rays and vertebrae which show generally decreased numbers linked with higher temperature (Ackerly & Ward, 2016;Fahy, 1972;Sfakianakis et al., 2011). Similarly, temperature can influence the sequence in which bones become mineralized with later developing bones being more susceptible to change (Mabee et al., 2000). Such temperature-induced alterations in the timing and the rate of bone ossification suggest an environmentally induced heterochrony (Rundle & Spicer, 2016). Therefore, the plastic effects of temperature (thermal plasticity) in bone could provide new insights into heterochrony and how it contributes to evolvability.
Skeletal traits are likely to be especially important for larval fish to escape predation (Domenici & Blake, 1997;McCormick et al., 2018;Pechenik, 2006;Walker et al., 2005). The standard "c-start" escape response of fish consists of a rapid bending of the body and tail to create the force necessary to propel the fish away from a potential threat (Domenici & Blake, 1991, 1997Webb, 1978). Skeletal anatomy is of particular importance for a c-start because the vertebrae provide both a point for muscle attachment but also influence the degree of body and tail curvature (Ackerly & Ward, 2016;Brainerd & Patek, 1998). Temperature is known to influence the number of vertebrae that form during ontogeny (Ackerly & Ward, 2016;Beacham, 1990;Fahy, 1972;Lindsey & Arnason, 1981;Sfakianakis et al., 2011) thus suggesting a direct route through which temperature may affect fitness.
Given that climate change is predicted to have greater effects at latitudinal extremes, this study focuses on a representative northern species; the Arctic charr (Salvelinus alpinus (L.)). Among freshwater fishes, Arctic charr have the most extreme northern distribution, display a wide range of phenotypic variation, and are known to possess high levels of plasticity (Adams & Huntingford, 2002;Alexander & Adams, 2000;Jonsson & Jonsson, 2001;Klemetsen, 2010;Parsons et al., 2010Parsons et al., , 2011. Phenotypic variation in this species is strongly linked to adaptive divergence (Adams et al., 1998;Skúlason et al., 1989;Snorrason et al., 1989), with variation in craniofacial and body-shape traits having frequently been shown to correspond to feeding and swimming performance (Bryce et al., 2016;Kapralova et al., 2015;Skúlason et al., 1989). Therefore, Arctic charr are an ideal model to provide insights into how temperature changes may impact bone development and in turn adaptation in fishes.
The early life stages of fish are more susceptible to predation than older age groups and early life can strongly influence later-life phenotypes (Johnson et al., 2014;McCormick et al., 2018;Pechenik, 2006). Thus, we focused on the developmental and ecological consequences of thermal plasticity for bone in the early-life stages when bone first forms. We combine an investigation of developmental change under different temperatures with ecologically relevant performance trials to provide some of the first insights into how climate change may impact evolvability. We predicted that fish exposed to novel higher temperatures would show greater rates of bone metabolism and ossification, as well as earlier ossification (Balon, 1980). We also predicted that different temperatures would induce plasticity in morphology over ontogeny and induce heterochrony while higher temperatures would increase the degree of morphological variation through the release of cryptic genetic variation. Finally, we predicted that higher temperatures would induce the formation of fewer vertebrae and that this would impede escape response performance (e.g. Ackerly & Ward, 2016;Sfakianakis et al., 2011).

| MATERIALS AND METHODS
To create sibling groups for study, Arctic charr adults of a pelagic ecomorph from Loch Rannoch, Scotland (Adams et al., 1998) were captured by gill net (bar mesh sizes of 30-35 mm) at spawning time in November 2013 and November 2015. Eggs (approx. 800 from each female) were stripped on site and fertilized (one female to one male) with milt from males creating 12 and 10 families in 2013 and 2015, respectively. Eggs from each family were evenly divided between two temperature treatments (5°C and 9°C). Charr collected in 2013 were used for investigations of ossification and morphological variation while those from 2015 were used to investigate bone remodeling, vertebral counts, and escape performance.

| Staining for cartilage and bone
To examine temperature effects, equal numbers of embryos were collected, from five families at both temperatures at two developmental stages; "prehatching" (n = 138) and at "first feeding" (n = 184) both respectively representing key stages of initial charr cartilage formation and subsequent bone mineralization as fish being actively foraging (Balon, 1980). Stages were standardized between temperature treatments by using developmental degree days (pre-hatch were sampled at 350 DD, first feeding at 780 DD). Embryos were euthanized with an overdose of benzocaine and placed into 4% paraformaldehyde for 4 h. To minimize variation in staining, fixed embryos were transferred to CaCO 3buffered 70% ethanol at 4°C until all specimens were collected; allowing all staining to be completed using the same batch of solutions.
To assess bone development, embryos were processed using an acid-free double staining protocol (Walker & Kimmel, 2007). Briefly, this protocol used alcian blue and alizarin red to stain for the presence of cartilage and bone, respectively. To remove excess stain, specimens (pre-hatch n = 138, first-feeding n = 184) were stepped through glycerol concentrations of 25%, 50%, and finally 75% until imaging (~2 months after fixing and staining). The protocol was adjusted from that of Walker and Kimmel (2007) by removing KOH from the glycerol series to prevent the disarticulation of bones.

| Quantification of bone and cartilage staining
Following staining, cartilage and bone were quantified across developmental stages by comparing images of the craniofacial region (defined as the element anterior to the opercular) in ventral and lateral aspects. Photographs were taken of the craniofacial region under standardized lighting conditions using a Leica M165 FC stereomicroscope mounted with a Leica DFC 450C digital camera using associated LAS v4.4 software (Leica Microsystems). Bone and cartilage were quantified from images by taking advantage of the double-staining technique. For the pre-hatch specimens, many did not present alizarin stain (representing bone ossification) indicating that no bone had formed. Therefore, a presence/absence measure for bone was collected. Cartilage at both stages, and bone at first feeding, were quantified using custom color vectors using the "color deconvolution" plugin (Ruifrok & Johnston, 2012) for ImageJ version 2006.02.01 (Schneider et al., 2012). Vectors for the quantification of alizarin red and alcian blue across all specimens were created from a reference specimen. Specifically, from the reference specimen, the anterior area of the mandible that was distinctly stained with alizarin red was selected to create a vector for bone, while an area of the hyomandibular region distinctly stained with alcian blue was selected to create a vector for cartilage. Together, these custom vectors allowed for consistent detection of RGB values specific to bone or cartilage stains across all specimens.
Bone and cartilage can vary in their degree of formation and spatial spread across a developing bone. Therefore, using double-stained images, they were quantified by measuring the area and intensity of stain (to create integrated density) in terms of pixels in ImageJ. The area was defined by the number of pixels in an outlined specimen corresponding to either cartilage or bone. The integrated density of the stain was measured by calculating the average intensity of corresponding pixels and multiplying by the area. A greater area of stain indicated that more bone or cartilage was present in a specific region while a greater integrated density of pixels indicated overall denser bone or cartilage. The area of the stain was calculated by first creating a black and white 8bit mask of the stained image using the "threshold_colour" plugin in ImageJ (Schneider et al., 2012). This plugin isolated the pixels corresponding to either bone or cartilage depending on the vector used. The 8-bit mask was then used to quantify pixels. Notably, in some cases, bone can form and result in the "masking" of some cartilaginous elements. The use of photos from the lateral and ventral aspects is, therefore, limited by this factor, meaning that our analysis, is representative of bone development but not completely comprehensive with regard to every possible bony element.
The effect of temperature on bone and cartilage was analyzed using two-way analyses of variance (ANOVA) in R 3.3.2 (R Core Team, 2016). Each stage and measure (integrated density and area in ventral and lateral aspects) were used with temperature, family, and their interaction as explanatory factors.

| Skeletal remodeling
To examine temperature effects on bone remodeling, samples were stained for osteoblast and osteoclast activity. Specifically, tartrate-resistant acid phosphatase (TRACP) and alkaline phosphatase (ALP) staining were conducted using kits from Takara (MK301; Clontech). ALP and TRACP are enzyme markers of osteoblasts (bone-forming cells) and osteoclasts (bone-resorbing cells), respectively. Due to the potential for the deterioration of TRACP over time, fish were euthanized and stained immediately. A total of 20 fish from each of 5°C and 9°C temperature treatment groups at the first feeding stage (SL = 23 ±2 mm) were euthanized (as above) with half being used for each of the ALP and TRACP stains. Specimens were kept in 4% paraformaldehyde, at 4°C overnight, rinsed in distilled water, and incubated at 37°C in one of the two staining solutions (15 min for AP stain or 30 min for TRACP stain). Specimens were subsequently washed with distilled water and bleached for 80 min in 3% H 2 O 2 at room temperature. Specimens were then stepped through 25%, 50%, and 75% glycerol and photographed as described above.
To quantify TRACP and ALP staining, new vectors were created in ImageJ. Area and intensity were quantified in the craniofacial region from both lateral and ventral aspects, as well as the anal and caudal fins from the lateral aspect. The effect of temperature on bone remodeling was tested using ANOVAs in R 3.3.2 (R Core Team, 2016).

| Measuring shape variation
Craniofacial shape variation was quantified using landmark-based geometric morphometrics. Software belonging to the tps suite (tpsDig264, tpsUtil64, tpsRegr64, tpsRelw64) was used for landmark data collection, analysis, and visualization and is freely available at http://life.bio.sunysb.edu/morph/index. html. Because we aimed to assess changes over ontogeny, and measured pre-hatched embryos before the formation of bony elements, we were limited to collecting four and six homologous landmarks for the ventral and lateral aspects, respectively, along with 15 semi-landmarks for the lateral aspect ( Figure 1). Semilandmarks were aligned using the minimum chord distance using tpsUtil64 in conjunction with tpsRelw64. Landmark coordinates were then adjusted using a General Procrustes Analysis (GPA) in CoordGen8 (available at https://www.animal-behaviour.de/imp/). GPA translates, rotates, and scales all specimens coordinate to minimize squared differences between corresponding landmarks (Rohlf & Slice, 1990). Following this, allometric effects were minimized by regressing GPA-adjusted landmarks against centroid size using standard7 (available at https://www.animal-behaviour. de/imp/).

| Temperature-induced shape change
We tested whether temperature, family, and their interactions could influence larval shape using two factor multivariate analysis of variances (MANOVAs) for each sampling stage and aspect with their partial warp scores as response variables using base functions in R (R Core Team, 2016). Partial warp scores are a rotation of the Procrustes residuals around the Procrustes mean configuration and characterize the location of each specimen in the space of partial warps (Bookstein, 1989).
We conducted two additional tests to assess (i) temperature-induced variation in the rate of shape change over development (heterochrony) (Webster & Zelditch, 2005), and (ii) the effect of temperature on the magnitude of potential shape differences at each of the pre-hatch and first feeding stages. Developmental rate variation differences due to temperature were tested by calculating the partial Procrustes distances (PPD) for each temperature between the prehatch and first feeding stage. The difference between these PPDs was then calculated between fish raised at 5°C and 9°C to provide an observed value. Statistical testing was conducted by bootstrapping individuals within sampling stages 900 times to produce a 95% confidence interval (CI) for each sampling stage separately at 5°C, 9°C, and the difference between them. For stage-specific tests, PPDs were calculated between fish raised at 5°C and 9°C for each of the "pre-hatch" and "first-feeding" stages. Both routines were conducted in TwoGroup8 (available at https://www. animal-behaviour.de/imp/).
To assess changes in the degree of morphological variation, indicating the possible release of cryptic genetic variation, we performed a series of tests focused on landmark disparity. Disparity is a multivariate analog of variance and can be used to quantify how much shape varies within and between treatments . Using this measure, we tested for changes in the degree of disparity using pairwise comparisons between temperatures at both the pre-hatch, and first feeding stages. To assess how disparity was altered over ontogeny, we then compared pre-hatch and first feeding stages under a single temperature. Statistical tests were performed by bootstrapping pairwise differences in disparity 1000 times using the DisparityBox functions within PCAGen8 (available at https://www.animalbehaviour.de/imp/).
To visualize temperature-induced shape changes, we performed a series of steps. Using partial warp scores, we performed a discriminant function analyses (DFA) using temperature as a grouping variable in the MASS package (Venables & Ripley, 2002) in R (R Core Team, 2016). Canonical root scores obtained from the DFA represented shape plasticity in response to 5°C and 9°C temperatures and were used as an independent variable in a multivariate regression on x,y coordinates to model their relationship to shape. Deformation grids visualizing the temperature-induced shape variation during ontogeny were created using the software tpsRegr and magnified 3× to highlight shape changes.

| Escape performance trials
To examine the adaptive relevance of any temperatureinduced changes, the escape performance of fish in response to a simulated predator was assessed at the first feeding stage (SL = 24 ± 2 mm) for fish exposed to both temperature conditions. Individual fish (n = 32 for each of the two temperatures) were randomly selected, and acclimated for 60 min at 7°C in a separate holding area (with the same dimensions as the experimental arena) before moving to the experimental arena for a minimum of 15 min before any stimuli. While this duration was unlikely to enable fish to recover completely from stress induced by netting, we noted whether individuals had returned to previous (pre-netting) rates of respiration (indicated by opercular movement) and swimming activity. If this did not occur within 15 min, the settling period duration was extended. The experimental arena consisted of a section of plastic pipe (31 cm in diameter and 15 cm deep) attached vertically to the bottom of a larger (90 × 40 × 22 cm) tank which was supplied with recirculating water at 7°C (±0.5°C).
Predator attacks were simulated by dropping an object into the tank following the protocol of Killen et al. (2015). The stimulus was an epoxy resin cylindrical object (diameter 30 mm, height 15 mm, and weight 15.5 g). To prevent visual stimulation before impact, the stimulus was dropped through a vertical pipe (295 mm high, 55 mm diameter) ending 0.5 cm above the water surface. Trials were filmed at 240 fps using a Hero4 camera (GoPro) from below. Individual fish were exposed to multiple trials (three per individual, excluding three fish which received two trials) with a minimum of 5 min between trials (following Ramasamy et al., 2015). Finally, to avoid cross-contamination of chemical cues, water in the experimental arena was changed between individuals.
To test for relationships between vertebrae number and escape response performance, specimens were immediately euthanized, cleared, stained, and photographed using the methods described above. These images were then used to count vertebrae for each individual across the entire vertebral column (performed twice for each specimen).

| Escape response video analysis
Individual escape response performance was assessed using a multistep approach. A cellophane 1-cm 2 grid on the underside of the tank was used for calibration in Im-ageJ (Schneider et al., 2012). Head and tail movements were tracked manually and measured frame-by-frame through each response using the ImageJ plugin MTrackJ (available at https://imagescience.org/meijering/software/ mtrackj/). Following Azizi and Landberg (2002), two stages of the escape response were defined. Stage 1 was defined as lasting from the first movement of the fish in response to the stimulus to the moment of maximum body curvature. Duration of Stage 1 was calculated as the number of frames between the beginning and end of the stage. Stage 2 was defined as beginning at the moment of maximum body curvature (the end of Stage 1) to when the spine of the fish was again straight. The duration of this stage was measured as the number of frames elapsed. Maximum curvature was measured using the bending coefficient expressed as the ratio of head-to-tail chord length (CL) to total length (TL) subtracted from 1: BC = 1 − (CL/TL) (Azizi & Landberg, 2002). Performance measures included Stages 1 and 2 duration, maximum curvature, head/tail velocity (maximum and average across Stages 1 and 2) and latency (defined as the time between the stimulus first making contact with the water surface and the first detectable movement of the head to begin a "c-start" escape response). Track measurements were used to determine maximum head and tail velocities for the entire escape response and average head and tail velocities for Stage 1 and Stage 2 individually. To minimize fish size effects, velocities were converted to body lengths per second.
To test the impact of developmental temperature on vertebral number, a t test was conducted. To test the impact of developmental temperature and skeletal changes (vertebral number) on performance, each trait was analyzed using a mixed effects repeated measures design with developmental temperature, vertebral number, and family (with individual nested within family as a random effect) included as factors using functions in the lme4 package in R (Bates et al., 2015).

| Temperature-induced effects on skeletal development
Our ANOVA analyses indicated that differences in bone and cartilage development were caused by temperature treatments (Table 1). Effects differed across regions of the skeleton and across stages of ontogeny. At the pre-hatch stage, fish raised at 5°C exhibited both a greater mean area and the integrated density of cartilage in the ventral aspect compared with fish at 9°C (Table 1). Similarly, all instances where calcified bone was present in the prehatch stage (10 out of 138 observations) occurred in fish T A B L E 1 ANOVA outputs for the effects of temperature, family, and their interaction on the area and density of stain for cartilage or bone at two defined stages of of Arctic charr development (pre-hatch or first-feeding)

Stage
Aspect | 339 raised at 5°C (Pearson's chi-squared test; χ 2 = 8.8; df = 1; p = .003). Also, at the pre-hatch stage, family by temperature interactions were present for both lateral and ventral measures of integrated density. At the first feeding stage, the same trends largely held for temperature effects, with fish raised at 5°C exhibiting greater levels of bone and cartilage than those at 9°C (Table 1). Family by temperature interactions were also present for lateral cartilage area.

| Skeletal remodeling effects induced by temperature
Across regions, bone remodeling was generally not affected by temperature. However, higher osteoblast and osteoclast activity was found in the anal fins of fish raised at 9°C relative to those at 5°C (Table 2). For osteoblasts, activity on the anal fin was higher (~15% increase in intensity and~10% increase in area of stain) in fish raised at 9°C. Additionally, osteoclast activity was affected by temperature, with greater osteoclast activity indicated in the anal (7% increase in the intensity of stain) and caudal (9% increase in the intensity of stain) fins of fish raised at 9°C (Table 2).

| Temperature-induced changes to skeletal morphology of the craniofacial region
The MANOVA models indicated that temperature, family, and their interactions influenced shape (Table 3).
For the ventral aspect, both family and temperature had significant effects on the shape at both developmental stages. For the lateral aspect, family and temperature showed a significant effect on the shape at both stages with an interaction between these factors at the first feeding stage indicating heritable variation in temperature-induced plasticity increased with age (Table 3).
Temperature also induced changes in the rate of morphological change, and magnitude of stage-specific shape differences. Heterochrony was present in both the ventral and lateral aspects (95% CIs did not cross zero) with the 9°C treatment inducing a faster rate of shape development in both cases (13.8% and 19.3% increases in rates respectively based on PPD differences). Stagespecific tests showed significant differences in all comparisons (all p < .01; Figure 2). However, consistent with our MANOVA results, temperature-induced differences were far greater at the pre-hatch stage relative to the first feeding stage (92.5% and 148.2% greater in the ventral and lateral views respectively based on PPDs).
Morphological disparity was broadly impacted by temperature across our comparisons. At the pre-hatch stage, both the lateral and ventral views differed in disparity (confidence intervals did not cross zero) (Table 4). Higher levels of disparity were present in the lateral aspect at 5°C, while the ventral aspect showed a decrease in disparity at this temperature. At the first feeding stage, the ventral aspect also showed higher levels of disparity at 5°C. Over ontogeny, all comparisons demonstrated a change in disparity with fish reared at 5°C displaying an increase in disparity from pre-hatch to first feeding stages, while all other comparisons showed decreases over ontogeny. Notably, and regardless of direction, changes in disparity over ontogeny were larger in fish reared under 9°C. Deformation grids showed that higher temperatures induced notably different shapes. Shape analysis from the ventral aspect showed that at 9°C, fish had a wider jaw at the prehatch stage, but this changed at the first feeding stage where at 9°C, fish showed a slightly narrower and longer jaw and a deeper craniofacial phenotype in the ventral aspect (Figure 2). Ontogenetic comparisons showed that jaw width generally narrowed and lengthened as fish reach the first feeding stage, but that the T A B L E 2 ANOVA outputs for models assessing the impact of temperature on levels of skeletal remodeling (indicated by osteoblast and osteoclast activity) narrowing of the jaw was enhanced in the 9°C temperature condition. For the lateral view, clear shape differences were induced by temperature, with 9°C fish showing a more blunted, shorter jaw region compared with fish from the 5°C treatment. Differences were more apparent at the prehatch stage, especially in the mandible length which was shorter in fish exposed to 9°C (Figure 3).

| Vertebral counts and escape performance
The number of vertebrae in fish ranged between 59 and 62 in both temperatures. However, temperature influenced the number of vertebrae expressed, with fish at 5°C having an average of 0.53 more vertebrae than those F I G U R E 2 Deformation grids representing the temperature-induced differences in the ventral aspect of charr shape across two developmental stages. Each arrow represents a pairwise comparison between the corresponding groups with a measure of partial Procrustes distance imposed upon them (all p < .001).
Comparisons are between (a) pre-hatch fish at 5°C and 9°C, (b) first-feeding fish at 5°C and 9°C, (c) pre-hatch and first-feeding stages at 5°C, and (d) pre-hatch and first-feeding stages at 9°C. Shape variation is magnified 3× to enhance interpretation at 9°C (mean no. of vertebrae 5°C = 60.57, mean no. of vertebrae 9°C = 60.03, t 1, 60 = 2.2; p = .03). Some parameters of escape response performance were also affected by developmental temperature, such as latency, with fish reared at 9°C reacting 53% faster (0.05 vs. 0.08 ms) to the stimulus than fish reared at 5°C (Table 5). Additionally, vertebral number, but not developmental temperature, affected the duration of Stage 1 responses, as fish with fewer vertebrae showed more rapid Stage 1 durations ( Figure 4).

| DISCUSSION
Here, we show that thermal conditions can induce extensive developmental changes that correspond to ecological performance traits. Elevated temperatures generally led to decreased levels of cartilage and bone deposition over ontogeny in charr. Elevated temperatures also increased bone remodeling, decreased the number of vertebrae, and extensively changed morphology in both the lateral and ventral aspects at both of the life stages we measured. Variation in morphology (measured as disparity) was also broadly altered by temperature indicating that development under different thermal conditions could alter the opportunity for selection. Thus, while physiological and life-history traits have been a focus for the effects of climate change in fishes (Crozier & Hutchings, 2014), it is apparent that bone development and functional morphology will also be widely impacted (Mabee et al., 2000;Ramler et al., 2014). We have demonstrated, on the basis of family by temperature interactions, that thermal responses in bone development and shape also possessed heritable variation in some cases. Coupling this with our measured effects on some aspects of escape response performance (latency and Stage 1 duration) suggests that such changes could be targeted by selection under changing thermal conditions and drive evolutionary change. While such effects are likely to vary in impact across species and populations, it is also possible that they will have strong effects on the functioning of ecosystems (i.e., through changes in predator-prey dynamics (Skúlason et al., 2019)).

| Temperature-induced changes in skeletal development
The lower deposition of cartilage and bone under a higher temperature persisted across developmental stages. This suggests long-term effects are possible as the processes involved with bone development are conserved over life history. For example, morphological plasticity could be impacted as more heavily ossified bone is more resistant to remodeling and mechanical stress (Parsons et al., 2014;Witten & Huysseune, 2009). Indeed, prior research on benthic charr ecomorphs, which are more heavily ossified than pelagic counterparts, shows they are less morphologically plastic (Parsons et al., 2011). Therefore, reduced bone ossification might allow for Comparisons are made at the pre-hatch and first-feeding stages between fish reared under 5°C and 9°C, and between stages under a common temperature. The calculated Foote disparities correspond to the first (1) and second (2) groups in order of mention under the comparison heading. Significant differences are indicated with an asterisk next to 95% bootstrapped confidence intervals.
CAMPBELL ET AL. rapid adjustments to changing environments, but could also incur disadvantages such as reduced skeletal strength Peterson & Müller, 2018). This would be important, for example, for fishes that use a biting or crushing mode of feeding that results in a higher mechanical load on craniofacial bones (Cooper et al., 2010). Indeed, limitations on bone could be mediated through temperature-induced changes in its metabolism. While limited to the anal fin, we found that bone metabolism (osteoblast/osteoclast turnover) can increase under higher temperature. This suggests that in some cases, bone would be more amenable to remodeling under climate change through increased morphological plasticity. However, considering that fish also develop as a whole more rapidly under warmer temperatures, bone would correspondingly form at a faster rate (Angilletta et al., 2004;McCarthy et al., 1998). Therefore, changes to bone remodeling and metabolism could underlie the observation of differential bone ossification in response to temperature. Higher temperature may lead to an increased osteoblast activity, but under a shorter developmental window, and based on our findings from measures of bone area and integrated density, is likely to reduce bone formation at a given developmental stage. Therefore, such changes in the developmental processes underlying bone development could change the duration of ontogenetic events that have evolved in coordination with current environmental conditions. Indeed, temperature also influenced charr head shape over early stages of ontogeny. While previous research has shown temperature-induced changes in shape for various fish species, these findings offer new insights due to our ontogenetic approach (Ramler et al., 2014;Rowiński et al., 2015;Sfakianakis et al., 2011). Specifically, we show that temperature had a substantially greater influence on shape earlier in development (i.e., pre-hatch stage), and that temperature-induced heterochronic shifts occurred with faster shape change at a higher temperature. Stronger temperature influences F I G U R E 3 Deformation grids representing the temperature-induced differences in the lateral aspect of charr shape across two developmental stages. Each arrow represents a pairwise comparison between the corresponding groups with a measure of partial Procrustes distance imposed upon them (all p < .001). Comparisons are between (a) pre-hatch fish at 5°C and 9°C, (b) firstfeeding fish at 5°C and 9°C, (c) pre-hatch and first-feeding stages at 5 o C, and (d) pre-hatch and first-feeding stages at 9°C. Shape variation is magnified 3× to enhance interpretation during early ontogeny could be due to weaker selection at this stage when yolk resources are plentiful. However, increasingly strong selective pressures should occur as larval fish enter the high mortality life-history stage around first feeding when yolk resources have been exhausted (Aral et al., 2011;Jørgensen & Holt, 2013;Pechenik, 2006;Sifa & Mathias, 1987). Therefore, environmental robustness in response to temperature variation would be favored to maintain growth and mortality rates (Jones et al., 2003). However, our results suggest that development may not be able to resolve temperature-induced changes in morphology with broad effects on the disparity of ventral and lateral aspects especially at the pre-hatch stage, and over ontogeny in the ventral aspect.
The observed changes in jaw length could, for example, alter foraging success. The shorter, wider jaws induced at 9°C would normally be associated with a biting mode of feeding. However, at the first feeding stage, Arctic charr are reliant on small zooplanktonic prey that require a suction mode of feeding. Future work could assess changes in foraging performance to infer how this might affect survivorship. In spite of our findings of shape change it is also possible that Arctic charr are employing compensatory mechanisms to correct for the early "mistakes" induced by temperature, and indicated by general reductions in disparity over ontogeny (see also Parsons et al., 2010). While costs for compensatory shape development have not been investigated, compensatory growth and its associated costs have been well documented for fish (Inness & Metcalfe, 2008;Metcalfe & Monaghan, 2001). Given that compensatory changes in shape would involve widespread anatomical alterations and the integration of traits, we would expect a potential for costs. Indeed, the integration of plastic responses in different morphological traits appears to be under the influence of selection in fish (Parsons & Robinson, 2006). Such potential costs could be reflected in a range of activities but may play a role in the changes we observed in escape response performance (see below).
Given that our temperature treatments were in line with climate change predictions, such induced changes are likely to have significant consequences in future populations. For example, charr are noted for their propensity to adaptive divergence within lakes, and there could be devastating impacts on the mechanisms of selection and development that maintain benthic/limnetic divergence (Skúlason et al., 2019;Skúlason & Smith, 1995;Smith & Skúlason, 1996). Given that increased temperature induced a more benthic morphology, it is plausible that specialization on the pelagic habitat will be less obtainable in the future due in part to temperature-induced developmental changes.
The evolutionary relevance of these temperatureinduced changes in phenotype are ultimately dependent | 345 on their heritability. Although not present in the majority of our comparisons, our data indicate that heritable variation can be substantial, through our observed family effects, and their interaction with temperature. Such interactions occurred in both bone development traits and in morphological shape indicating that evolutionary potential could be altered under climate change. Nonetheless, effects on heritability are likely to be stage-and trait-specific as family by temperature interactions for lateral shape were only present at the later first-feeding stage. Stage-specific selection has previously been shown to underly adaptive morphological divergence in Arctic charr ecomorphs making these findings especially intriguing . It could be that Arctic charr will simply evolve in direct response to future conditions and the developmental effects they incur. Alternatively, populations could evolve to shift phenology in a way that matches the optimum temperature profile for their early development. This could effectively "shield" heritable variation from selection by maintaining developmental conditions. However, to be successful, this strategy would also rely on shifts in resource availability for emerging charr, especially at the critical first-feeding stage when yolk resources are unavailable. Thus, understanding heritability in phenology, in addition to the traits here, could be vital for predicting how populations will respond to climate change.

| Evolvability of temperatureinduced changes: Escape response performance
Temperature induced changes in the number of vertebrae, with higher numbers occurring at lower temperature (5°C). This was in line with previous research on fish (Ackerly & Ward, 2016;Sfakianakis et al., 2011) and Jordan's rule whereby the number of vertebrae in fishes increases with latitude and correspondingly cooler temperatures (McDowall, 2008). This suggests that Jordan's rule may also correlate with other phenotypic effects (bone development, morphology) and in part be the result of plasticity. Notably, changes to vertebrae number have functional consequences for flexibility (Long & Nipper, 1996;McDowall, 2003), maximum body curvature for escape swimming (Brainerd & Patek, 1998), body elongation associated with ambush predation (Maxwell & Wilson, 2013), burst swimming (Swain, 1992) and C-start velocity (Ackerly & Ward, 2016). Indeed, the number of vertebrae and temperature influenced the escape responses of our charr. Specifically, fish with fewer vertebrae displayed a shorter Stage 1 duration in escape response trials, which should benefit their performance (Table 4, Azizi & Landberg, 2002). In support of this, we also found that temperature, but not vertebral number, influenced latency time of responses with fish raised at 9°C reacting more rapidly than those at 5°C (Table 4). Therefore, increased temperature could benefit escape responses in prey species and likely cause downstream ecological effects through changing predatory/prey dynamics (McCormick et al., 2018;Skúlason et al., 2019). However, it should be noted that preliminary data on vertebrate counts in fish from eggs collected in 2013 actually showed an opposing trend, with 9°C inducing an average of 0.8 more vertebrate (n = 22 for 9°C, n = 19 for 5°C, t 39 = −3.9079, p < .001). This suggests that while temperature can consistently impact vertebrate number, the direction of change may F I G U R E 4 Scatterplot depicting the relationship between vertebral number and Stage 1 duration (F 1, 171 = 13.08, p < .001). The points have been "jittered" using the ggplot2 package in R to better visualize the spread of data over common data points. Red points represent charr reared at 9°C, while blue points represent charr reared at 5°C interact with other factors. In our case, this may due to the year fish were collected and corresponding differences in maternal factors. Future research investigating the effects and interactions of multiple environmental factors will be valuable for addressing this idea.
Predicted increases in temperature from climate change will undoubtedly alter developmental conditions at a broad scale. This investigation explains, in part, how this will impact adaptive phenotypic variation in the most northerly distributed species of freshwater fish and provides a basis for predictions about how populations may respond (Campbell et al., 2017). However, the phenotypic changes observed here could serve as the basis for a set of "early warning" trait changes for conservationists concerned with fishes, and any number of taxa. Conceptually, this is because development is much more likely to elicit phenotypic responses prior to the demographic markers currently used to recognize threats to populations. At its best, conservation should aim to preserve the evolutionary process, therefore, developmental approaches, as we have taken here, could be especially effective for the persistence of biodiversity because they can reveal effects on evolvability.