Inter- and Intraspecific Variation Patterns of Leaf Phenotypic Traits in Two Sympatric Pine Species

Abstract

Individual persistence under changing climate conditions can be aided by phenotypic plasticity. Needle morphology reflects pine species’ adaptation to their habitats, and adaptive plasticity is beneficial to pine survival. As two closely related pines, Pinus massoniana Lamb. and Pinus hwangshanensis W. Y. Hsia sympatrically occur in the subtropical region of China, forming hilly forests (i.e., subtropical evergreen needleleaf forests) at lower elevations and montane forests (i.e., temperate evergreen needleleaf forests) at medium and high elevations. However, little is known about the intraspecific phenotypic variation patterns of the two pine species and their relationship with environmental factors. Here, we investigated the interspecific and intraspecific variation patterns of needle traits in the two sympatric pine species, focusing on six traits—i.e., leaf length, leaf thickness, leaf width, leaf area, specific leaf area, and leaf dry matter content—as well as the relationship between needle traits and environmental factors. We found significant phenotypic trait differences among populations of each species. Except for leaf length and leaf width, variation within species was greater than variation between species in the needle traits measured. Even more intriguingly, the leaf size traits (i.e., length, width, thickness) were more conservative than the leaf economic traits (i.e., specific leaf area and leaf dry matter content). In other words, the intraspecific variability of the former was weaker than that of the latter. The nature of P. massoniana needle traits was mainly shaped by latitude, while the needle traits of P. hwangshanensis were significantly affected by annual precipitation. Therefore, phenotypic plasticity may be an essential mechanism for the two pine species to better cope with changing external conditions. The intraspecific variation patterns found in the two pine species and the relationships between traits and environmental factors can provide substantial scientific data for large-scale exploration of intraspecific phenotypic variation in pine species and their breeding practice.

1. Introduction

Understanding phenotypic variation (i.e., the measurable or observable variation in a trait) or phenotypic plasticity (i.e., the ability of a given genotype to express different phenotypes in different environmental circumstances) is necessary to investigate a trait’s function and fitness, as well as its ecology and evolution [1,2]. Among other hypotheses, phenotypic plasticity is assumed to facilitate evolution via ‘buying time’, ‘plasticity-led evolution’, and ‘non-genetic evolution’ [3]. Environmental phenotypic change precedes evolutionary adaptation, according to the ‘Plasticity-first’ hypothesis, which uses a mechanism similar to ‘buying time’, [4]. However, a thorough review of related literature reveals that plasticity can both support and impede evolutionary change, and that it is not always adaptive but can also be neutral or even maladaptive [5].

Conversely, several eco-evolutionary mechanisms—such as dispersal, selection, plasticity, and (epi)genetic mutations—can influence how phenotypic variation is expressed in unexpected contexts or under environmental stress [2]. Empirical studies, reciprocal transplanting strategies, and standard garden experiments have demonstrated that genetic or epigenetic characteristics [6] and environmental factors (e.g., geography, climate, and soil) [7] affect plants’ phenotypic variation. Despite this, empirical studies on phenotypic variation in plants (i.e., incredible natural populations) are greatly limited.

Pines are a dominant species in forests of the Northern Hemisphere [8] and form large-scale coniferous forests in Asia, Europe, and North America. With 111 species divided into 2 subgenera, 4 sections, and 11 subsections, the genus Pinus enjoys a rich species diversity and an intricate evolutionary history [9]. With the deepening of global changes in climate and biodiversity, pines are facing many challenges, including death from drought stress [10] and disease attacks [11]. Undoubtedly, observing phenotypic variation is critical for delving deeply into pine trees’ strategies and mechanisms for coping with global changes.

Among the pine species, P. massoniana and P. hwangshanensis are two genetically closely related taxa that belong to the subgenus Pinus, section Pinus, and subsection Pinus. Here, P. hwangshanensis is treated as an independent species after being merged with P. taiwanensis—a species found in Taiwan [12]. The evolutionary events of these pines are related to climate and geography. On the other hand, these two pine species are considered to have distinct ecological characteristics and form coniferous forests with different properties. P. massoniana is a pioneer species that forms a subtropical evergreen needleleaf forest [13], whereas P. hwangshanensis is a pioneer or climax species that forms a temperate evergreen needleleaf forest [14]. More interestingly, the horizontal ranges of P. massoniana (27°55′~32°30′ N, 102°40′~109°55′ E) and P. hwangshanensis (22°48′~31°48′ N, 112°30′~122°30′ E) are highly overlapped. However, there are noticeable differences in their elevation distribution. Generally, the boundary is 700 m, with P. hwangshanensis at the medium and high elevations and P. massoniana at the low elevations [15]. It is known that plants’ functional traits are expected to change as altitude changes, with plants at higher altitudes having smaller and thicker leaves than plants at lower altitudes [16].

Phenotypic variation in physiology and growth has been observed in pine species populations for needle traits such as leaf length, leaf width, leaf thickness [17], epidermis thickness with cuticle, resin duct diameter [18], stomatal row density, specific leaf area [19], and leaf area [20]. Meanwhile, existing research has reported that needle length, needle length–width ratio, the biomass of needles, and needle longevity follow a precise longitude and latitude pattern [21,22], and that hydraulic efficiency traits exist in a pattern along the altitudinal gradient [23]. However, the results demonstrated that in garden experiments, conifer traits associated with cold and drought resistance may be affected by a combination of genotype and environment [24]. Moreover, it has been emphasized that multiple analytical approaches have led to the conclusion that the environment has an overwhelming impact on survival and growth [25]. Since the two viewpoints are still inconclusive, we focused on the effects of the environment on the needle traits of P. hwangshanensis and P. massoniana.

When studying phenotypic variation, it is important to note that plants frequently exhibit substantial intraspecific variation in leaf traits—especially in species with heterogeneous environments, along environmental gradients, and with wide ranges [26]. Relevant data show that up to 40% of total variation in some plant traits occurs within species [27]. With increasing altitude, needle length decreases, and leaf dry matter content decreases predictably [28]. Indeed, intraspecific variation arises from heritable differences and plasticity [29], as in the case of Pinus pinaster Ait., which has evolved to adapt to its local environment, resulting in high genetic differentiation among populations, along with geographically highly structured intraspecific variation in quantitative traits and neutral genetic variation [30]. Based on previous studies, we explored the sources of variation in P. hwangshanensis and P. massoniana through empirical studies.

Two species of pine trees from different elevation gradients on five mountains were selected for the present study, and a total of six leaf morphological traits—leaf length (LL), leaf thickness (TL), leaf width (DL), leaf area (LA), specific leaf area (SLA), and leaf dry matter content (LDMC)—were determined for various groups of pine needles. LL, TL, DL, and LA are traits that are related to needle size and are intuitive trait indicators; SLA is indicative of leaf biomass, and LDMC reflects the efficiency of resource utilization and the environmental adaptation strategy [31]. By using a broad survey of needle structure, we tested the following hypotheses: (1) the six needleleaf traits of the two pine species have a wide range of variation, with LDMC plasticity being significantly higher in P. hwangshanensis than in P. massoniana; (2) the total variation in most needle leaf traits of the two pine species has higher intraspecies variation than interspecies variation; and (3) the needle leaf morphological traits of the two pine species are influenced by different environmental factors, among which latitude and annual precipitation are the most significant.

2. Materials and Methods

2.1. Study Area and Sampling

The study areas were located in Jiangxi, Anhui, and Zhejiang, China; these provinces are among the core ranges of P. hwangshanensis and P. massoniana. In order to investigate the impact of elevation on the leaf phenotypic traits of two closely related pines, we identified five mountains as the objects of investigation. P. massoniana and P. hwangshanensis are located at low and high elevations, respectively. Five different populations for each pine—namely, Mt. Sanqingshan (SQS, 1820 m above sea level), Mt. Lushan (LS, 1474 m), Mt. Jinggangshan (JGS, 1779 m), Mt. Huangshan (HS, 1864 m), and Mt. Tianmushan (TMS, 1506 m)—were involved in our investigation (Figure 1). All five mountains are located in the central subtropical zone, belonging to the subtropical monsoon climate; the climate varies significantly vertically due to the influence of elevation, and the soils are yellow loam, red loam, yellow–brown loam, and other types, mostly related to elevation changes [32,33,34].

The geographic locations for each population are shown in

Table 1. The geographic locations and main meteorological factors of the sample populations of P. massoniana and P. hwangshanensis.

Population	Abbreviation	Longitude (°E)	Latitude (°N)	Annual Average Temperature (°C)	Annual Average Precipitation (mm)
Mt. Sanqingshan, Jiangxi	SQS	118°00′~118°06′	28°52′~28°57′	11	1860
Mt. Lushan, Jiangxi	LS	115°50′~116°10′	29°28′~29°45′	12	2024
Mt. Jinggangshan, Jiangxi	JGS	113°39′~114°23′	26°27′~26°40′	17	1890
Mt. Huangshan, Anhui	HS	116°49′~118°59′	29°13′~31°05′	14	2583
Mt. Tianmushan, Zhejiang	TMS	118°36′~120°06′	29°52′~30°55′	14	1870

. According to the elevation gradients of the two pine species, two or three sites (i.e., elevations) were selected for each population of the two pine species. A 20 × 50 m quadrat was set up for each site. The number of sites was 20 and 15 for P. massoniana and P. hwangshanensis, respectively (Table 1 and Table S1).

Needle samples were collected from 60 P. massoniana individuals and 44 P. hwangshanensis individuals. Three mature individuals of P. massoniana and P. hwangshanensis with a similar diameter at breast height (DBH) were selected for each site, and fresh needles were randomly sampled from different directions around the tree crowns. The needles were mixed and sealed in a zip-lock bag and then returned to the laboratory for determination of their morphological characteristics.

2.2. Measurements

Ten fascicles of needles were selected randomly from each sampled tree for the measurement of LL, DL, and TL. LL was measured with a metric ruler with an accuracy of 0.1 cm; DL and TL were measured using an outside micrometer with an accuracy of 0.01 mm. An electronic balance was used to measure the saturated fresh weight and dry weight of the leaves, with an accuracy of 0.0001 g. The term ‘saturated fresh weight of leaves’ refers to the mass of fresh needles that have been immersed in distilled water for 36 h and then extracted after being dried on the surface. Following the aforementioned determination, the needles were dried in an oven at 70 °C until the mass remained constant, and the dry weight of the leaves was then measured. Then, the leaf area (LA), specific leaf area (SLA), and leaf dry matter content (LDMC) were calculated separately.

Referring to the method of Zhang et al. [35], the needle trait indices were calculated as follows:

$$\mathrm{Area} \mathrm{of} \sin \mathrm{gle} \mathrm{cambered} \mathrm{surface} \mathrm{of} \mathrm{a} \mathrm{needle} (S_{\mathrm{SCS}}): S_{\mathrm{SCS}}=\frac{\mathsf{\pi } \times \mathrm{DL} \times \mathrm{LL}}{2}$$

(1)

$$\mathrm{Leaf} \mathrm{area} (S_{\mathit{LA}}): S_{\mathit{LA}}=\frac{\mathsf{\pi } \times \mathrm{DL} \times \mathrm{LL}}{2} + \mathrm{DL} \times \mathrm{LL}$$

(2)

$$\mathrm{Specific} \mathrm{leaf} \mathrm{area} (S_{\mathit{SLA}}): S_{\mathit{SLA}}=\frac{\mathrm{leaf} \mathrm{area}}{\mathrm{leaf} \mathrm{dry} \mathrm{mass}}$$

(3)

$$\mathrm{Leaf} \mathrm{dry} \mathrm{matter} \mathrm{content} (C_{\mathit{LDMC}}): C_{\mathit{LDMC}}=\frac{\mathrm{leaf} \mathrm{dry} \mathrm{mass}}{\mathrm{saturated} \mathrm{fresh} \mathrm{mass}} \times 100\%$$

(4)

2.3. Data Analysis

Means, coefficients of variation (CVs), and standard deviations for each trait among the populations were estimated. Data from each single needle measurement were used to conduct analysis of variance.

We log-transformed the data for all needle traits and performed a nested ANOVA using the lme [36] and varcomp [37] functions in R [38], in the following increasing order: individuals, elevations, populations, and species. The code was as follows (e.g., LL):

LL < - lme (log10 (LL) ~1, random = ~1|species/population/elevation, data = E, method = “REML”)

mLL < - varcomp (LL, TRUE, cum = FALSE)

The same model was applied to TL, DL, LA, SLA, and LDMC, only replacing LL with these traits in the R code. REML refers to a restricted maximum likelihood method in the ‘lme’ function of R (version 4.2.1).

Clustering analysis among populations of P. massoniana and P. hwangshanensis was performed using the average method in R by the unweighted pair group method with arithmetic mean (UPGMA). The average values of each needle trait were used for cluster analysis. In order to investigate the effects of external factors (geographic and environmental) on the needle traits of the two pine species, redundancy (RDA) analysis was used. The environmental and geographic data for each sampling site were collected by searching WorldClim (http://www.worldclim.org/, accessed on 20 September 2022), using three geographic factors (latitude, longitude, and elevation) and two environmental factors (annual average temperature and annual precipitation). The statistical analysis of this study was completed in SPSS 25.0 and R 4.2.1 [39].

3. Results

3.1. Variation of Needle Traits in the Two Pine Species

The mean values of seven needle morphological traits of P. massoniana and P. hwangshanensis were significantly different between the populations (

Table 2. Mean values, coefficients of variance (CV%), and analysis of variance of needle traits in P. massoniana. The population abbreviations are the same as in Table 1. Numbers in bold symbolize statistical significance.

Traits	SQS		LS		JGS		HS		TMS		ANOVA
	Mean	CV%	Mean	CV%	Mean	CV%	Mean	CV%	Mean	CV%	F	p (>F)
LL (cm)	11.305	13.48	13.574	17.89	15.638	11.05	13.658	15.98	14.189	11.09	140.600	<0.01
TL (mm)	0.404	4.11	0.423	9.04	0.416	9.62	0.466	15.07	0.403	15.39	62.402	<0.01
DL (cm)	0.079	7.10	0.088	14.15	0.081	11.95	0.080	8.48	0.075	18.98	56.719	<0.01
LA (cm2)	2.288	16.00	3.075	20.98	3.255	20.11	2.798	19.33	2.731	22.72	87.540	<0.01
SLA (cm2/g)	149.348	20.46	140.607	19.39	164.363	17.06	137.186	15.37	167.568	25.85	50.105	<0.01
LDMC (g/g)	0.443	16.85	0.455	18.91	0.439	18.12	0.532	14.09	0.402	17.52	83.717	<0.01

Abbreviations: SQS, Mt. Sanqingshan; LS, Mt. Lushan; JGS, Mt. Jinggangshan; HS, Mt. Huangshan; TMS, Mt. Tianmushan.

and

Table 3. Mean values, coefficients of variance (CV%), and analysis of variance of needle traits in P. hwangshanensis. The population abbreviations are the same as in Table 1. Numbers in bold symbolize statistical significance.

Traits	SQS		LS		JGS		HS		TMS		ANOVA
	Mean	CV%	Mean	CV%	Mean	CV%	Mean	CV%	Mean	CV%	F	p (>F)
LL (cm)	9.477	20.84	8.623	19.52	11.438	10.31	9.763	17.96	9.671	14.14	58.209	<0.01
TL (mm)	0.427	4.94	0.477	10.36	0.498	9.13	0.466	19.51	0.540	8.89	97.383	<0.01
DL (cm)	0.094	8.66	0.099	10.64	0.104	11.29	0.094	13.63	0.106	9.94	50.904	<0.01
LA (cm2)	2.289	22.94	2.190	20.17	3.066	16.20	2.370	26.23	2.662	20.86	61.234	<0.01
SLA (cm2/g)	142.461	17.61	121.316	20.69	129.870	15.83	120.501	19.83	117.189	16.31	27.646	<0.01
LDMC (g/g)	0.474	21.74	0.485	21.18	0.249	94.80	0.515	16.30	0.461	84.71	22.928	<0.01

Abbreviations: SQS, Mt. Sanqingshan; LS, Mt. Lushan; JGS, Mt. Jinggangshan; HS, Mt. Huangshan; TMS, Mt. Tianmushan.

). The variation coefficients of morphological traits of P. massoniana ranged from 4.11% to 25.85%, and the average variation coefficients of LL, LA, and LDMC were relatively small, while TL, DL, and SLA showed a great range of variation. Additionally, the fluctuation in needle traits of P. hwangshanensis was marked with variation coefficients of 4.94~94.80%. LDMC and TL showed significant variation, while the average variation of the other traits was slight. Interestingly, according to the boxplot, seven needle traits within a mountain were markedly different between the P. massoniana population and the P. hwangshanensis population, and such patterns occurred on all five mountains (Table S1, Figure 2).

3.2. Sources of Needle Traits’ Variation in the Two Pine Species

The four levels of species, populations, elevations, and individuals had different effects on the variation of the seven functional leaf traits (Figure 3). Except for LDMC, the sources of variation for the needle leaf traits consisted mainly of interspecific and intraspecific contributions. Most of the variation between populations was smaller than the variation within populations for most of the measured needle traits. The variation in LL and DL was mainly derived from interspecific variation, with contributions of 61.48% and 55.73%, respectively. Compared with the other four traits, intraspecies variation contributed more to LA and SLA. The variation in TL was mainly due to intraspecific variation (35.84%), followed by interspecific variation (29.50%). Elevation contributed the most to LDMC (63.79%).

3.3. Cluster Analysis

The two species were divided into two groups by cluster analysis (Figure 4), except for JGS clustered into another group. For P. massoniana, the populations in LS and HS were clustered into one group, and the populations in TMS and JGS were clustered into another group. Similar results occurred in P. hwangshanensis, where the populations in LS and HS (i.e., the high-latitude populations) were clustered into one group.

3.4. Redundancy Analysis

The results of the RDA analysis showed that the first three RDA ranking axes included 44.65% of the information on needle traits and their corresponding geographic variables for P. massoniana, and latitude (La) had more influence on the distribution of needle traits than other factors (Figure 5a). On the other hand, the first three ordination axes in Figure 5b explain 62.01% of the information of P. hwangshanensis needle traits and their corresponding environmental variables. Annual precipitation (MAP) had the greatest influence on the distribution of needle traits.

4. Discussion

4.1. Phenotypic Variation and Plasticity

The study of needle traits in P. massoniana and P. hwangshanensis revealed a wide range of phenotypic variation (Table 2 and Table 3). Furthermore, P. massoniana had less variation than P. hwangshanensis—for example, the variation in the leaf dry matter content of P. hwangshanensis was pronounced (up to 98.4%). The latest histological studies also proved that P. hwangshanensis differed significantly from P. massoniana in the formulation phase [40]. Stapley et al. [41] proposed that phenotypic variability was one of the foundations in determining the intrinsic factors that define species and their biogeography, and that adaptation in any species requires phenotypic variation in traits [1,42]. In fact, nearly all traits will be influenced by plasticity [5]. Physiological changes, on the other hand, are thought to be central to phenotypic variation, and because of physiology’s systemic nature and sensitivity to environmental changes, the division of labor between traits could be regulated by these changes—particularly the endocrine system [43]. That is, divergence in phenotypic plasticity could result from physiological mechanisms with a genetic basis. This allows us to have a more in-depth discussion of phenotypic variation and phenotypic plasticity.

In the present study, the coefficient of variation of needle traits was used to characterize the ability of plasticity; the average coefficients of variation of LL, TL, and DL were lower than those of LA, SLA, and LDMC, indicating the low plasticity in LL, DL, and TL. Despite the fact that our study did not calculate the phenotypic plasticity index (PPI) [44], the coefficient of variation produced results that were consistent with it [45]. Based on 20 morphological and anatomical needle traits and 10 geographic populations of Pinus tabuliformis, Zhang et al. [46] reported that the inter-population heritability of LL, TL, and DL was greater than 0.8, indicating good genetic stability. Meanwhile, they had low levels of plasticity and small coefficients of variation when compared to other traits.

4.2. Phenotypic Variation Patterns of Leaf Traits among and within Species

Conifers are widely assumed to be typical wind-pollinated species, with chloroplast DNA that is primarily inherited paternally via pollen [47]. However, evidence has suggested that geographical barriers may result in low levels of differentiation between populations by limiting the pollen dispersal [48]. A similar conclusion was reached by Meng et al. [49], who reported that the genetic diversity of most needle traits was found to be most abundant within the population of Pinus tabuliformis in 18 needle traits. Xie et al. [50] used SRAP markers to study five endemic Pinus taxa from China (Pinus tabulaeformis, P. tabulaeformis var. mukdensis, P. tabulaeformis f. shekanensis, Pinus massoniana, and Pinus henryi), which showed more genetic diversity within taxa, rather than between taxa.

When investigating genetic variability and adaptive differentiation in conifers, it is well known that variation in needle variability traits is a useful visual indicator [51,52]. Within populations of Pinus yunnanensis, the proportion of phenotypic variation in nearly all needle traits was greater than 50%. In other words, trees within populations accounted for the vast majority of total variation [53]. Based on the phenotypic variation of needle traits in Pinus tabuliformis, Zhang et al. [46] found that more variation was expressed within populations. Similar results were obtained in our study (Figure 3), which demonstrated that the intraspecific components were generally higher than interspecific components in most of the measured needle traits. However, different needle traits showed different levels of variation within species; intraspecific trait variation (ITV) was relatively high for morphological traits linked to the leaf economics spectrum (e.g., leaf N and P, SLA, LDMC) [54], but low for traits related to leaf size (i.e., LL, TL, DL), which are typically considered to be independent of the leaf economics spectrum but have been linked with adaptation to broad climatic gradients [55].

Additionally, it has been demonstrated that intraspecific variation could have consequences that are comparable to or even greater than those of interspecific variation [27], and that genetics and habitat divergence are the primary forces behind phenotypic variation within species [56].

4.3. Adaptive Phenotypes in the Context of Changing External Factors

Researchers began to understand as early as the middle of the 20th century that environmental changes in a variety of morphological and physiological traits could precede and facilitate genetic changes [57]. Levis and Pfennig [58] underlined the dual function of environmental change as a generator and then selector of phenotypic variation, which reoriented how we think of the ‘plasticity-led evolution hypothesis’. In other words, by enabling organisms to react appropriately and immediately to environmental variation, plasticity may confer a special evolutionary benefit [3]. As a result, environmental variation is essential in phenotypic plasticity as well as evolutionary adaptation.

Phenotypic trait variation and genetic architecture are highly influenced by geographic and environmental variations—particularly those linked to photoperiod and temperature (e.g., latitude, longitude, and elevation) [6]. Our study’s RDA analysis showed that the two pine species’ needle traits were affected by various environmental factors to various degrees (Figure 5). Elevation had a significant impact on needle traits, and the results of earlier investigations [53,59] were consistent in that LL and LDMC showed a negative correlation with altitude. Since P. hwangshanensis populations are distributed at higher elevations, they are subjected to more severe winter stress than P. massoniana populations, which are found at lower elevations. Thus, their needles are thicker and shorter to withstand the stress of strong winds and snow [60]. In addition, trees that are located at high altitudes face a lack of available water, leading to drought stress [61], and the lower leaf area minimizes stomatal transpirational water loss. Moreover, due to the short growing season at high altitudes and the tendency for soil moisture to freeze in winter—resulting in root uptake—MAP is extremely crucial for this species. P. hwangshanensis, on the other hand, had a higher level of plasticity in LDMC than P. massoniana, which can be attributed to nutrient deficiencies and low-temperature stresses caused by increased elevation. Cluster analysis could support the peculiar latitude patterns that SLA and LDMC show in P. massoniana and P. hwangshanensis (Figure 4). SLA decreased as latitude increased, while LDMC increased as latitude increased. Due to the fact that SLA is associated with plants’ growth rates and ecological strategies, it was discovered that Pinus koraiensis has smaller LDMC and larger SLA at low latitudes [61]. The greater the SLA, the more effectively leaves photosynthesize in order to adapt to changing habitats [62]. The investment in defense against low temperatures increases with latitude [63]. Based on ‘leaf economic spectrum theory’ [54], we could conclude that both P. massoniana and P. hwangshanensis were typical ‘slow investment-return’-type species. However, our findings indicated that the resource trade-off strategies of the two pine species differed markedly.

P. massoniana and P. hwangshanensis share a common geological history across the mainland, which accounts for their similar levels of genetic diversity [64]. Due to the presence of refugees during the ice ages, climatic changes might have low impacts on these populations at mid-latitudes. Therefore, P. massoniana was able to achieve population expansion, while population decline has only recently started, indicating that its population was dealing with serious survival issues [48]. In cold biomes, where the climate is changing most rapidly, the ability of a species to adjust its phenotypes under new circumstances is vital [7].

5. Conclusions

In sum, the needle traits of P. massoniana and P. hwangshanensis displayed significant heterogeneity, whereas the leaf size traits displayed more consistent heritability between populations. Likewise, the variation of all needle traits within populations was greater than that between populations, suggesting that both genetics and environment played a role in the intraspecific phenotypic variation of the needles. The variation of needle traits was affected by different environmental factors (i.e., longitude, latitude, elevation, temperature, and precipitation), reflecting the adaptability of the traits to the heterogeneous environment. Limited by knowledge, this study is the first to systematically report the intraspecific and interspecific variation patterns of needle traits in multiple populations of P. hwangshanensis. Our findings might provide a theoretical reference for the adaptive evolution and breeding practice of Pinaceae.