Growth and Its Relationship to Individual Genetic Diversity of Mountain Hemlock (Tsuga mertensiana) at Alpine Treeline in Alaska: Combining Dendrochronology and Genomics

Globally, alpine treelines are characterized as temperature-limited environments with strong controls on tree growth. However, at local scales spatially heterogeneous environments generally have more variable impacts on individual patterns of tree growth. In addition to the landscape spatial heterogeneity there is local variability in individual tree genetic diversity (level of individual heterozygosity). It has been hypothesized that higher individual heterozygosity will result in more consistent patterns of growth. In this article, we combine genomics and dendrochronology to explore the relationship between individual genetic diversity and tree growth at a mountain hemlock (Tsuga mertensiana Bong. Carr) alpine treeline on the Kenai Peninsula, Alaska, USA. We correlated average observed individual heterozygosity with average tree-ring width and variance in tree-ring width within individuals to test the hypothesis that trees with higher individual heterozygosity will also have more consistent growth patterns, suggesting that they may be more resilient to climate and environmental fluctuations at the alpine treeline. Our results showed that there was no significant relationship between tree growth and individual heterozygosity. However, there was a significant positive relationship between average tree-ring width and variance in tree-ring width implying that overall, fast growing trees in stressful environments, such as the alpine treeline, grow unstably regardless of the level of individual heterozygosity.


Introduction
Patterns of forest tree growth at landscape and regional scales are often affected by changes in broad scale site quality and climate [1]. However, at local scales heterogeneous environments generally have more variable impacts on individual growth [2,3]. Many characteristic growth patterns in trees, such as leaf size [4], shape [5], height [6] and radial growth [7] are linked to environmental factors like temperature [8], moisture availability [4,5,9], competitive interactions [10,11], and stress [12]. Various morphological and physiological changes can occur when trees are subjected to environmental stress. For instance, changes in foliar growth, leaf/needle length, root growth and the production of defense compounds can be altered across heterogeneous sites [13]. Alpine treelines are unique site that are spatially heterogeneous boundaries between sub-alpine forest and alpine tundra. At global scales treeline position is temperature controlled [14], but at local scales it is structured by microclimatic, topographic, edaphic, and biotic variation [15]. Alpine treeline thus represents a harsh environment with strong influence on spatial patterns of radial growth [7,16,17]. Despite the fact that patterns of growth are influenced by various biotic and abiotic processes, no consensus has been reached regarding the degree to which the stability of growth, or the overall reduction in year-to-year variation in growth, in plants is related to individual genetic diversity. That is, if plants have higher genetic diversity (individual heterozygosity), do they also have more consistent and stable growth patterns (homeostasis) relative to environmental variation [18,19] like that found at the alpine treeline? Understanding how individual trees will respond to changing environmental conditions is important. Conifer trees have high levels of within population genetic diversity and low between population genetic differentiation [20,21]. In addition, conifers have high degrees of phenotypic plasticity [22] allowing them to respond differentially to changing environmental conditions. Climate change has the potential to alter geographic patterns of plant ranges. Often, the most rapid response for a plant is to shift its range to higher latitudes or elevations through seed dispersal [23][24][25]. However, in the case of slow growing long lived species, like conifers, the ability to disperse over long distance to escape changing conditions may be curtailed by the length of time it takes for the species to germinate and grow to reproductive maturity resulting in an adaptational lag [24].
One alternate response to dispersal and range shift is to rely on their plastic responses to withstand variable environmental changes. To this end, individual trees that have less variable radial growth across heterogeneous environments have a greater fitness advantage and may be the subject of natural selection and local adaptation. This idea, known as developmental or genetic homeostasis, is not new [18,19,26,27]. Positive relationships between heterozygosity and reduced morphological variance have been identified in several organisms [28]. Despite the aforementioned relationship, previous results in forest trees, and conifers in particular, have been less consistent. For instance positive associations between heterozygosity and basal growth were identified in quaking aspen (Populus tremuloides Michx.) [29] and height in Douglas fir (Pseudotsuga menziesii (Mirb.) Franco) [30], but in lodgepole pine (Pinus contorta Dougl. ex Loud.), and ponderosa pine (Pinus ponderosa Dougl. ex Laws.) no significant relationship was found [31], while mixed results were identified in pitch pine (Pinus rigida Mill.), where some stands had a positive growth relationship with heterozygosity while others had a negative association [32]. Most of these earlier studies relied on allozyme loci. Because of the inconsistencies in earlier results the subject is worth exploring in more detail under the lens of new genomic approaches which allows for a high resolution investigation into the patterns of individual genetic variation and growth.
A recent study by Babushkina et al. [33] introduced the idea of correlating individual heterozygosity to growth parameters using dendrochronological techniques in Siberian larch (Larix sibirica Ledeb.). The authors of the study compared average individual tree-ring width (AvTRW), variance (VarTRW) and individual heterozygosity (IndHet), using eight highly polymorphic microsatellite markers, also known as simple sequence repeats (SSRs), with the assumption that individual trees with higher heterozygosity would positively correlate with heterosis and have more stable homeostasis reflecting less developmental dependence on their environment. With this in mind, Babushkina et al. [33] tested the hypotheses that there would be a negative correlation between IndHet and VarTRW and a positive correlation between IndHet and AvTRW. Ultimately, their findings were inconclusive. They found non-significant positive correlations between IndHet and both VarTRW and AvTRW. They also found that AvTRW and VarTRW were significantly and positively correlated with each other. Babushkina et al. [33] suggested that perhaps under unfavorable growing conditions most fast growing trees grow unstably regardless of the level of genetic diversity. However, they posit that with only eight microsatellite markers the resolution of their analysis may be limited and recommended an analysis using genome wide genotyping with thousands of molecular markers.

Objective
In this study, we tested the association between tree growth at the alpine treeline, measured as an individual average tree-ring width (AvTRW) along an elevational gradient, and individual genetic diversity, measured as an observed individual heterozygosity (IndHet) averaged over 4665 genomic single nucleotide polymorphism (SNP) loci. We assessed the patterns of individual growth and genetic diversity in a treeline formed by a conifer species, mountain hemlock (Tsuga mertensiana Bong Carr.) on the Kenai Peninsula, Alaska, USA. We followed the approach outlined in Babushkina et al. [33] to clarify our understanding of the molecular basis underlying stability and growth in forest species with a focus on the question: do individual trees with higher heterozygosity, approaching their altitudinal limit, exhibit more stable patterns of growth? We answered our research question by testing two hypotheses: (1) IndHet will be negatively correlated with VarTRW; and (2) IndHet will be positively correlated with AvTRW.

Study Species-Mountain Hemlock
Mountain hemlock is a slow growing, monoecious, and wind pollinated species. Additionally, the species has long-distance seed dispersal capacity [34,35] and is highly outcrossed due to extensive pollen dispersal [36][37][38][39]. The species is found principally in wet cool environments on the Pacific coast of North America, where its range extends from the Alaskan Kenai Peninsula in the north to the Sierras in northern California [37,40], with disjunct populations occurring as far as Idaho, 160 km away from the main distribution [35]. Mountain hemlock is the dominant conifer within the subapline zone of the Kenai Mountains and is an important component of the forest below this zone where it co-occurs with Sitka spruce (Picea sitchensis (Bong.) Carriè) and white spruce (Picea glauca (Moench) Voss.).
Mountain hemlock reaches sexual maturity and begins producing cones between 20 and 30 years of age [37]. Seeds are primarily wind-dispersed, with mast year cone crops capable of producing 215,000 to 4,144,000 seeds/ha [37]. At the alpine treeline in Alaska, mountain hemlock stands can occur in open parkland form. Above the alpine treeline the species can take on a krumholz form. Vegetative reproduction, through layering, does occur in the krumholz form in Alaska [37], though it is unknown how common this form of reproduction is. An analysis of clonality at the mountain hemlock treeline on the Kenai Peninsula, Alaska found no instances of vegetative spread [34]. Mountain hemlock growth is negatively correlated with spring snowpack depth and positively correlated with summer growing season temperature [17,40], with warm July temperatures resulting in increased seed production [41]. Average growing season temperatures on the Kenai Peninsula have increased by 0.03 • C per year over the past 76 years (Figure 1). Warming should result in a shift of mountain hemlock treeline to higher elevations in the future.

Study Area-Palmer Creek Drainage: Kenai Peninsula, Alaska
We based our study along a west-facing elevational gradient, from the alpine treeline down to the valley floor in the Palmer Creek drainage of the Chugach National Forest on the Kenai Peninsula, Alaska. The study area was located in the north central portion of the Peninsula in the Chugach-St. Elias Mountains ecoregion [42] at about 60°47′37″ N and 149°32′11.35″ W ( Figure 2). The vegetation structure is typical of the Kenai Mountains and is composed of shrub communities at the lower elevations of the ecotone, primarily willow (Salix sp.) and alder (Alnus sp.), embedded in a matrix of bluestem (Calamagrostis sp.) with white spruce and mountain hemlock occurring as the dominant conifer species. Mountain hemlock dominates the alpine treeline ecotone, which then transitions into alpine lichen tundra at approximately 800 m a.s.l. The study site was previously described in an analysis of seed dispersal within the alpine treeline ecotone [34].

Study Area-Palmer Creek Drainage: Kenai Peninsula, Alaska
We based our study along a west-facing elevational gradient, from the alpine treeline down to the valley floor in the Palmer Creek drainage of the Chugach National Forest on the Kenai Peninsula, Alaska. The study area was located in the north central portion of the Peninsula in the Chugach-St. Elias Mountains ecoregion [42] at about 60 • 47 37" N and 149 • 32 11.35" W ( Figure 2). The vegetation structure is typical of the Kenai Mountains and is composed of shrub communities at the lower elevations of the ecotone, primarily willow (Salix sp.) and alder (Alnus sp.), embedded in a matrix of bluestem (Calamagrostis sp.) with white spruce and mountain hemlock occurring as the dominant conifer species. Mountain hemlock dominates the alpine treeline ecotone, which then transitions into alpine lichen tundra at approximately 800 m a.s.l. The study site was previously described in an analysis of seed dispersal within the alpine treeline ecotone [34]. The green points represents the study site both at the state scale and on the Peninsula. The transect map shows the geographic position of each tree on the sampled transect with circle size related to the age of the tree and the purple gradient representative of the level of individual observed heterozygosity (IndHet) measured at 4665 genome wide single nucleotide polymorphisms (SNPs); the darker the color purple, the more genomically diverse the individual (see Table 1 for IndHet values). Histogram inset: Age structure of Mountain Hemlock along the study transect representing the chronology length from 1900 to 2012.

Sampling Plant Material
Increment cores and needle tissue of mountain hemlock were collected in July 2013, along a single 830-m long elevational transect extending from the highest mountain hemlock individual on the slope down to the valley floor (between 880 and 610 m a.s.l.) ( Figure 2). We observed no evidence of human disturbance along the transect. We sampled 32 mountain hemlock individuals along the transect for analysis in this study. We collected foliage from each tree (approximately 15-cm long branch tips with young needles attached) for DNA extraction, and placed the needles in a plastic zip-lock bag with silica gel for preservation [43]. We stored the bags at approximately 10 °C in the field until they were shipped back to the laboratory. A basal tree-core was collected from the selected trees with basal diameters >5 cm using an increment borer. We recorded the geographic location of each sample using a handheld GPS (±3 m). Upon return from the field the foliage samples were placed in a −20 °C freezer to await DNA extraction and sequencing. The green points represents the study site both at the state scale and on the Peninsula. The transect map shows the geographic position of each tree on the sampled transect with circle size related to the age of the tree and the purple gradient representative of the level of individual observed heterozygosity (IndHet) measured at 4665 genome wide single nucleotide polymorphisms (SNPs); the darker the color purple, the more genomically diverse the individual (see Table 1 for IndHet values). Histogram inset: Age structure of Mountain Hemlock along the study transect representing the chronology length from 1900 to 2012.

Sampling Plant Material
Increment cores and needle tissue of mountain hemlock were collected in July 2013, along a single 830-m long elevational transect extending from the highest mountain hemlock individual on the slope down to the valley floor (between 880 and 610 m a.s.l.) ( Figure 2). We observed no evidence of human disturbance along the transect. We sampled 32 mountain hemlock individuals along the transect for analysis in this study. We collected foliage from each tree (approximately 15-cm long branch tips with young needles attached) for DNA extraction, and placed the needles in a plastic zip-lock bag with silica gel for preservation [43]. We stored the bags at approximately 10 • C in the field until they were shipped back to the laboratory. A basal tree-core was collected from the selected trees with basal diameters >5 cm using an increment borer. We recorded the geographic location of each sample using a handheld GPS (±3 m). Upon return from the field the foliage samples were placed in a −20 • C freezer to await DNA extraction and sequencing.

Genome-Wide Marker Development
We processed tissue for sequencing by first grinding approximately 30 mg of dry needle tissue, and extracted genomic DNA following a modified Cetyltrimethylammonium bromide (CTAB) protocol [44]. We used a genotyping-by-sequencing (GBS) approach, ddRAD-Seq [45], to genotype sampled individuals. Genome-wide markers were developed as outlined in Johnson et al. [39]. In brief, after digestion with the restriction enzyme (RE) pair SphI-MluCI the paired-end (PE) sequencing libraries with~350 bp long inserts were generated and sequenced with 150 × 2 cycles in a single HiSeq 2000 lane (Illumina, San Diego, CA, USA). The library preparation and sequencing were carried out at the University of Texas Genomic Analysis and Sequencing Facility.
Bioinformatic analysis was performed by the Texas A&M Institute for Genome Sciences and Society. We assessed sequence read quality; trimmed or removed low quality reads, de-multiplexed samples, de novo aligned reads and assembled them into contigs, and identified genomic SNP variants using the dDocent pipeline [46]. SNPs were quality filtered using VCFtools [47]. We filtered the SNP dataset by selecting a 10× minimum coverage depth cutoff and a Phred quality score >30. We selected one SNP per RAD tag. We discarded SNP loci that were not in Hardy Weinberg equilibrium (HWE) at p < 0.05, failed to genotype in greater than 80% of individuals, with a minor allele frequency less than 5%, and SNP loci that were not bi-allelic. Individuals missing more than 15% of the identified SNPs were also discarded. The obtained SNPs were supposedly selectively neutral based on the F ST outlier analysis in Bayescan [48]. Most of them were likely located in the non-coding regions and, therefore, also increased their chances of being selectively neutral and making them and our data comparable to that of [33].

Genetic Analysis
We calculated the level of individual genetic diversity (IndHet) by averaging individual heterozygosity across all loci using GenAlEx [49]. Hardy-Weinberg equilibrium was assessed using pegas [50], and its significance was tested by permutation via 1000 simulations.
We tested our prediction of a negative association between VarTRW and IndHet and a positive relationship between AvTRW and IndHet using Pearson's correlation coefficient and linear regression. Moreover, we tested the degree to which AvTRW, VarTRW, and IndHet were correlated with elevation. We tested the significance of the correlations using a Student's t-distribution and stats in the R base package [51]. Significance was assessed at 95% confidence.

Tree-Ring Processing and Analysis
To assess the age structure, AvTRW, and VarTRW along the transect, we measured tree age and ring-width from our collected increment cores following Lafon [52]. We dried cores in an oven for a minimum of 24 h at 100 • C and then sanded them with progressively less abrasive sand paper (80 to 400 grit) to draw out the cellular structure of the annual rings [53]. We developed a master chronology [54] using a subset of the longest cores collected along the transect. This allowed us to identify significant marker rings which facilitated visual crossdating. We dated the remaining increment cores under varying magnification with a stereomicroscope. For any core not showing the pith, we estimated establishment date from the curvature and width of the innermost ring [55]. In most cases, we did not add more than five years when estimating establishment date. We created histograms to investigate the age structure of the transect binned by decade.
We measured the ring width of all 32 cores to the nearest 0.001 mm using a Velmex tree-ring measurement system. We used COFECHA to visually crossdate the series. Spaghetti plots of individual ring width series were produced using dplR in R [53,56]. We calculated dendrochronology related descriptive statistics, including mean ring width (mm), mean sensitivity (MS), standard deviation (SD), series inter-correlation, and autocorrelation (AR1), to evaluate the quality of our chronology [53,57]. A simple horizontal line (method = "mean") was used to normalize the raw ring widths. The normalized tree-ring widths were then used to calculate mean and variance which were then assigned as AvTRW and VarTRW, respectively, to sampled individuals.
In order to assess if the values of heterozygosity in our sample were representative of the larger treeline population at the study site we compared the distributions of calculated IndHet to that of the entire population. Though we do not have increment cores for all trees along the study transect, we do have exhaustive ddRADseq genotypes for all trees (n = 163) described in [34]. The increment cores used in this study represent one fifth of the entire population. We calculated IndHet for each of the 163 census trees and compared the two distributions using a one sample Kolmogorow-Smirnov test using the mean and standard deviation from the census data as the reference distribution. We also tested for equivalent variances in the two distributions using an F-test for two sample variances.

Genome-Wide Marker Development
Our GBS resulted in 41,057,267 unique sequencing reads and 50,952 de novo assembled contigs with an average per nucleotide read depth coverage greater than 30×. Using the dDocent [46] pipeline 171,019 putative SNPs were identified throughout the mountain hemlock genome. Quality filtering reduced the number of SNPs using VCFtools [47]. We identified 4665 SNPs distributed across the mountain hemlock genome. Two individuals had greater than 15% missing genotype information and were removed from genomic analysis (however they were retained in development of the master tree-ring chronology).

Genomic Diversity
Population genomic diversity was moderate among the 30 individuals (Table 1)  Our analysis of the distribution of IndHet in our samples and that of the surrounding treeline population showed no difference between the distributions (D = 0.130, p = 0.641) or the variances (F = 0.839 (29,164), p = 0.49) ( Figure S1; Table S1) suggesting that despite our limited sample size our findings should be representative of the larger population of mountain hemlock. Table 2. Pearson's correlation coefficients with p-values in parentheses calculated in 30 trees to assess associations between tree growth (average tree-ring width (AvTRW) and variance (VarTRW)), elevation (Elev, m) and individual genomic diversity (IndHet). Bold values are significant at p < 0.05.

Tree-Ring Processing and Analysis
From 32 cores we constructed a 112 year long chronology of Mountain Hemlock at the alpine treeline on the Kennai Penunsula, Alaska (Table 3, Figure 5). Our chronology was qualitatively similar to other studies of mountain hemlock carried out in the region [58,59]. Our histogram of age class, binned by decade, identified minor spikes in establishment in the 1950s and 1980s (Figure 2). Our chronology failed to capture establishment for the decade 1911-1920.

Tree-Ring Processing and Analysis
From 32 cores we constructed a 112 year long chronology of Mountain Hemlock at the alpine treeline on the Kennai Penunsula, Alaska (Table 3, Figure 5). Our chronology was qualitatively similar to other studies of mountain hemlock carried out in the region [58,59]. Our histogram of age class, binned by decade, identified minor spikes in establishment in the 1950s and 1980s (Figure 2). Our chronology failed to capture establishment for the decade 1911-1920. Chronology length (CL), mean ring width (mm) (MRW), mean sensitivity (MS), standard deviation (SD), series inter-correlation between trees (Rbar), autocorrelation (AR1), signal-to-noise ratio (SNR) and expressed population signal (EPS).

Discussion
Our analysis identified no significant relationships between AvTRW, VarTRW and IndHet. We hypothesized a negative trend between VarTRW and IndHet as a function of growth and development being less dependent on the environment at higher levels of heterozygosity (homeostasis). Our results revealed a non-significant weak positive trend in the data. The direction of the trend was opposite of the expectations we laid out for the study (Table 2, Figure 4). AvTRW had a non-significant weak positive trend with IndHet matching the expectation that if growth were to be considered an adaptive trait, it would be positively correlated with individual genomic diversity due to heterosis. One of the most surprising findings was the significant positive relationship between AvTRW and VarTRW which, as Babushkina et al. [33] suggested, supports the idea that under poor environmental conditions fast growing trees grow unstably regardless of the level of heterozygosity. Though it is difficult to draw sound mechanistic explanations about this finding, we can speculate on the causes of this phenomenon. Alpine treelines are temperature limited and generally have short growing seasons constituting a poor growing environment. Because of this trees will generally grow more slowly regardless of the level of IndHet. This growth pattern is often confirmed at treeline where older trees will have small diameters relative to trees of the same age in more productive systems in the continuous subalpine forest [16,60]. However, when relatively favorable conditions do occur at the treeline, trees with higher heterozygosity should respond more favorably by increasing radial growth leading to a positive association between growth and variance. Moreover, the positive association between AvTRW and VarTRW could be the

Discussion
Our analysis identified no significant relationships between AvTRW, VarTRW and IndHet. We hypothesized a negative trend between VarTRW and IndHet as a function of growth and development being less dependent on the environment at higher levels of heterozygosity (homeostasis). Our results revealed a non-significant weak positive trend in the data. The direction of the trend was opposite of the expectations we laid out for the study (Table 2, Figure 4). AvTRW had a non-significant weak positive trend with IndHet matching the expectation that if growth were to be considered an adaptive trait, it would be positively correlated with individual genomic diversity due to heterosis. One of the most surprising findings was the significant positive relationship between AvTRW and VarTRW which, as Babushkina et al. [33] suggested, supports the idea that under poor environmental conditions fast growing trees grow unstably regardless of the level of heterozygosity. Though it is difficult to draw sound mechanistic explanations about this finding, we can speculate on the causes of this phenomenon.
Alpine treelines are temperature limited and generally have short growing seasons constituting a poor growing environment. Because of this trees will generally grow more slowly regardless of the level of IndHet. This growth pattern is often confirmed at treeline where older trees will have small diameters relative to trees of the same age in more productive systems in the continuous subalpine forest [16,60]. However, when relatively favorable conditions do occur at the treeline, trees with higher heterozygosity should respond more favorably by increasing radial growth leading to a positive association between growth and variance. Moreover, the positive association between AvTRW and VarTRW could be the result of tradeoffs in reproductive output [31]. This hypothesis argues that if seed production (reproduction) and growth are negatively correlated, a plausible tradeoff at alpine treeline, increased growth would not by itself increase the fitness of individuals there. This finding would mean that greater variance for the growth trait would occur due to more fluctuations in cone production in fast growing trees. We cannot add additional insight to this hypothesis here, and we encourage further research into this potential relationship.
In their study, Babushkina et al. [33] identified the same three trends using eight highly polymorphic SSR markers across two populations of L. sibirica. Both AvTRW (r = 0.146, p = 0.147) and VarTRW (r = 0.122, p = 0.225) had non-significant and positive trends with IndHet, as well as a significant positive correlation between AvTRW and VarTRW (r = 0.726, p < 0.001). The authors of the aforementioned study concluded that the relationships between heterozygosity and growth were clearly complex and non-linear. The lack of correlation between individual heterozygosity, derived from the SSRs, and growth parameters could have been a result of the markers failure to accurately reflect underlying genomic diversity. Because of this Babushkina et al. [33] recommended the use of genome wide sequence data as a means to better capture underlying individual genomic diversity. Despite their recommendation, using more than 4000 SNPs, our study failed to identify a significant relationship between growth and genomic diversity. Even more surprising was the consistency in results between the two studies suggesting that even when using high resolution genomic sequence data patterns of homeostasis may be difficult to detect and may be non-linear.
In terms of climate and growth as measured by the ring width index (RWI), our chronology statistics were qualitatively similar to earlier studies carried out in T. mertensiana on the Keani Peninsula and in the region broadly [58,59]. Response function analyses of the ring-width chronology identified significant correlations with both July precipitation in the year preceding growth and June temperature in the year of ring formation, as the dominant climate factors influencing radial growth in T. mertensiana (unpublished data). Peterson and Peterson [17] established that at high elevations in the Pacific Northwest radial growth in T. mertensiana was positively correlated with growth-year summer temperature and negatively correlated with spring snowpack.
It must be noted, the results of our study should be viewed as a case study and care should be taken when extrapolating the findings to other systems. Our study was based on a limited number of individuals, though representative of the local population. Future work should focus on a more robust set of samples from several tree species.
It is tempting to suggest that individual trees with higher genomic diversity should be more resilient under the poor growing conditions found at alpine treeline. Moreover, under future climate change, it is important to begin to untangle the alternative strategies that forest trees employ in response to rapid changes in their environment. Previous studies of heterosis and genetic homeostasis have shown that in some systems there is a relationships between higher individual genetic diversity and stable growth patterns [26]. In fact, both our study and that of Babushkina [33] have been unable to find support for this hypothesis in two different conifer species using two different numbers and types of genetic markers, thereby calling into question the validity of the hypotheses. The mountain hemlock treeline on the Kenai Peninsula, Alaska has been characterized by high pollen and seed dispersal into the ecotone [34]. Perhaps, the high degree of gene flow, particularly from pollen, contributes to the lack of association between IndHet and AvTRW/VarTRW. Genetic diversity is often structured according to the center-periphery (central-peripheral or central-marginal) model. In this model higher gene flow from large central populations into small peripheral ones, such as the altitudinal range limit constituting the alpine treeline, maintains genetic diversity [76]. This phenomenon is often associated with reduced adaptive potential of individuals at the range edge. This occurs because the influx of genes adapted to the center of the range counters the impact of selection for traits suitable to the surrounding environment (e.g., gene swamping) [77][78][79][80][81]. If gene swamping is occurring at the alpine treeline, then the patterns of growth may not accurately reflect the underlying genetic diversity represented in the local gene pool, and this may ultimately limit the adaptability of trees at their range limits. The degree to which the level of heterozygosity contributes to homeostasis needs to be further examined using a more extensive spatial sampling across an expanded geographic extent among a variety of treeline types. We recommend that analysis incorporate both selectively neutral and putatively adaptive genomic markers to more accurately assess levels of heterozygosity. The use of dendrochronological techniques in combination with genomic sequence data shows promise, and the development of longer chronologies may help to untangle some of the climatic variability allowing a longer temporal scale assessment of growth.

Conclusions
This study has shown that even when using a high resolution genomic dataset, the relationship between individual genetic diversity and tree growth at the alpine treeline is likely not linear, and the patterns may be difficult to detect. We observed no relationship between IndHet and both AvTRW and VarTRW. However, our findings did show that AvTRW and VarTRW were significantly and positively correlated suggesting that under poor environmental conditions, such as those found at the alpine treeline, fast growing trees grow unstably regardless of the level of heterozygosity. The combination of genomics and dendrochronology as an analytical approach is promising, and future research untangling the relationship between individual genetic diversity and patterns of growth will benefit from using larger populations of trees across multiple species to further address the hypothesis.