Introduction

Climate models forecast an increase in temperature along with an intensification of extreme weather events that will likely lead to more prolonged and intense drought periods around the globe (IPCC 2007, 2014). This future climate may have serious impacts on ecosystem dynamics in the Mediterranean basin (Resco de Dios et al. 2007; Sardans and Peñuelas 2013). Seasonal drought is one of the main factors affecting Mediterranean forests, and forest tree species use different strategies to respond to changes in water availability. Among those, access to deep water pools is fundamental to survive dry periods (Dawson et al. 2020; Ripullone et al. 2020; Rossatto et al. 2012), because water uptake is dependent on root architecture (Hernández et al. 2010). A number of functional traits related to hydraulic conductivity such as root vessel diameter are determined by root xylem anatomy, which in turn is influenced by root diameter and rooting depth (Kirfel et al. 2017; Wang et al. 2015). Thus, root structure and morphology play a key role on the ability of plants to explore and access the deeper subsurface when shallower soil layers dry out (Andivia et al. 2019; Padilla and Pugnaire 2007). Species-specific and intraspecific differences in rooting patterns are still poorly characterized, however. This is mainly due to the inherent complexity of underground sampling (Alani and Lantini 2020) and the fact that traditional approaches for root monitoring, such as the Auger and the Monolith methods, are destructive and non-repeatable (Krainyukoy and Lyaksa 2016).

The Ground-Penetrating Radar (GPR) is a non-destructive geophysical prospecting device that utilizes electromagnetic wave pulses for subsurface detection based on the different dielectric properties of materials, as well as on the physical properties of the medium (Lorenzo et al. 2010). GPR transmits electromagnetic signals into a surface and records the reflection of the signal transmitted by a receiving shielded antenna. The amplitude of the reflected pulses and its receiving time can be used in order to estimate the position and the size of the material that generates a discontinuous signal (Butnor et al. 2001). GPR is widely used for many purposes, such as the detection of bedrock structure (Luodes and Sutinen 2011; Valerio et al. 2012) or the quantification of soil water content (Klotzsche et al. 2018). Also, it has been used as non-invasive tool for the detection of root features and belowground biomass estimation (Barton et al. 2004; Butnor et al. 2001, 2003; Guo et al. 2013a, b; Lorenzo et al. 2010). Roots have a different dielectric permittivity compared to soil and, as a result, GPR allows defining the number (or frequency of occurrence; frequency herein), position (depth) and size (diameter) of coarse roots through the timing and features of back-reflected signals induced by hyperbolic reflections on radargrams (Hirano et al. 2012, 2009; Wu et al. 2014). The reflected waves are detectable as hyperbolas, where the hyperbolic signatures’ peak (or the maximum amplitude of the hyperbola) corresponds to the center of every object detected, which is recorded on a portable control unit (Hirano et al. 2012). However, the successful identification of roots with GPR depends on different factors, with soil moisture, root water content and soil properties influencing wave frequency and thus root detection (Barton and Montagu 2004; Hirano et al. 2009).

Pinus halepensis Mill. (Aleppo pine) is the most widely distributed conifer species across the Mediterranean basin (Vennetier et al. 2018), occupying vastly different ecological niches ranging from mesic (annual precipitation exceeding 800 mm) to xeric environments (below 350 mm). The widespread distribution of Aleppo pine suggests the existence of intraspecific adaptive divergence (Grivet et al. 2013; Ruiz Daniels et al. 2018). Indeed, there is strong evidence for genetic variability among populations originating from distinct environments to tolerate water stress, which involves traits such as water-use efficiency (Voltas et al. 2008), wood anatomy (Esteban et al. 2010), biomass allocation (Chambel et al. 2007), or use of water sources (Voltas et al. 2015). Deeper and more developed root systems might also be related to the ability of Aleppo pine to cope with water deficit (Andivia et al. 2019; Voltas et al. 2015). As a result, root morphology may be a key characteristic to understand water use strategies and, consequently, drought resistance mechanisms in this species. Although a large investment in roots is often related to enhanced survival under drought (Grossnickle 2005; Matías et al. 2014), and Aleppo pine has a greater ability for colonizing the soil subsurface and produce a more efficient root system in dry environments than mountain pines (Andivia et al. 2019), the existence of intraspecific divergence in root traits for this Mediterranean pine has barely been investigated thus far (Voltas et al. 2015).

In this work, we applied GPR root detection to a common garden of adult Aleppo pine trees where populations representative of the circum-Mediterranean distribution of the species were tested. Our purpose was to evaluate if GPR is a valid tool to infer intraspecific variability for coarse root frequency, depth and diameter in a widespread tree species. We hypothesized that populations originating from xeric environments would present deeper, thicker and more numerous coarse roots than their mesic counterparts as a result of different adaptive mechanisms underlying the use of soil water pools. In addition, we also hypothesized that variability in root traits among populations of the species is related, to some extent, to intraspecific divergence in aboveground growth, which may point to the existence of allocation trade-offs to aerial growth or belowground dry matter. Therefore, we aimed to: (1) test the use of GPR as a high-throughput root phenotyping technique to assess population differentiation in coarse root architecture for Aleppo pine; (2) relate intraspecific differences in coarse root traits (frequency, depth and diameter) among Aleppo pine populations to their climate at origin, which can be informative of different adaptive strategies in the species; and (3) describe the associations between such root traits and aboveground growth and define potential allocation trade-offs at the intraspecific level.

Materials and methods

Plant material and experimental site

Seeds from 56 populations of P. halepensis originating from mainland Spain, Balearic Islands, France, Greece, Italy and Tunisia were used in this study. These populations are representative of most of the current distribution range of the species (Fig. 1). The seeds were collected in 1995 from 20 to 30 adult individuals per population, spaced at least 100 m apart, and were sown in a forest nursery in Spain using standard container practices (Landis et al. 1990) the following year. In 1997, 896 one-year old seedlings (16 seedlings per population) were transplanted at the study site (provenance trial), which is located in Altura, Castellón province, Spain (39°49′29ʹʹN, 00°34′22ʹʹW, 640 m a.s.l.). The experimental design was a Latinised row-column design with four replicated blocks (John and Williams, 1998). The spacing between each consecutive tree at row or column distances was 2.5 m. Each row was ca. 70 m long and was composed of seven linear plots or experimental units consisting of four trees of the same provenance per plot (Fig. S1). The total area of the provenance trial was about 0.8 ha. For this study, we used three of the four available blocks, as some individuals from the fourth block were affected by fire in 2012. The Latinised row-column design allowed to efficiently control for inherent intra-site variability as systematic changes in both column and row directions were incorporated into the fitted mathematical models.

Fig. 1
figure 1

Geographic origin of 56 Pinus halepensis populations (coloured dots) evaluated in a common garden (red star) located in Altura (Castellón province, Spain). The legend shows ecotype codes as defined in Table S2. The dark green area represents the natural distribution of P. halepensis according to EUFORGEN (http://www.euforgen.org/species/pinus-halepensis/). A climograph of the trial site is included (mean climate of 1995–2019 obtained from Worldclim)

The soil at the study site is a calcic cambisol with loam texture (44.2% sand, 24.3% clay and 31.5% silt) and a maximum depth of ca. 40 cm followed by a petrocalcic horizon with vertical fractures (visual inspection). The general features of the soil are outlined in Table S1. The site has climatic characteristics similar to the average climate of the species across the Mediterranean basin (Santini et al. 2019a). The mean annual temperature is 13.8 °C and the mean annual precipitation is 468 mm, with 18% falling in summer.

Climate variables

For each population, climate data at origin for the period 1970–2000 were obtained from the WorldClim database (Fick and Hijmans, 2017) at 1 km2 spatial resolution. The following climate variables were retrieved, based on previous studies on the climatic drivers of ecotypic variation in Mediterranean pines (Climent et al. 2008; Tapias et al. 2004): mean annual temperature (MAT), mean summer temperature (MST), temperature annual range (TAR; Tmax-Tmin), mean annual precipitation (MAP) and summer to annual precipitation ratio (PsP). Moreover, solar radiation accounting for cloud cover (SR) was retrieved from WorldClim and monthly means of day-time vapour pressure deficit for the warmest (summer) quarter (VPDs) were calculated for each population by subtracting the water vapour pressure (VP) from the monthly average of saturation water vapour (VPsat) obtained at daily level. VP was estimated as in Ferrio and Voltas (2005):

$$\mathrm{ln}\left(VP\right)=6.34+0.47{T}_{m}+0.96({P}_{m}/1000)-0.22(Z/1000)$$
(1)

where Tm is mean monthly temperature, Pm is monthly precipitation and Z is altitude (in m).

In turn, VPsat was estimated from day-time temperature (Tday) as (Jones and Corlett 1992):

$${\text{VP}}_{\text{sat}}\text{ = 613.75 }{\text{exp}}\text{ (17.50}\frac{{\text{T}}_{\text{day}}}{\text{240.97+ }{\text{T}}_{\text{day}}}\text{)}$$
(2)

The populations were further categorised into 16 ecotypes (Fig. 1, Table S2) following previously published approaches based on climate attributes (Climent et al. 2008).

Aboveground growth

Tree height (H) and diameter at breast height (DBH), measured in 2013 (at age 17 years), were used as aboveground growth traits. Additionally, tree crown area was estimated from aerial RGB imagery taken in summer 2016 (at age 20 years) through a canopy height model using the R package ForestTools (Plowright 2018), as reported in Santini et al. (2020). We assumed constancy of phenotypic ranking in tree growth from age 15 onwards, as previously reported for P. halepensis (Sbay and Zas, 2018).

GPR data collection

GPR records were obtained in a single day in mid-June 2015, when trees were 19 years old, and the number (frequency), depth and diameter of coarse roots were estimated using a MALÅ RAMAC X3M GPR (MALÅ Geoscience AB, Sweden) equipped with an 800 MHz shielded antenna coupled to an inspection cartwheel. The 800 MHz antenna was used because it has been shown to provide the best possible resolution in calcic soils (Rodríguez-Robles et al. 2017), which are distinctive of the study site. In particular, roots reaching a depth of up to ca. 0.9 m could be identified within the petrocalcic fractures underneath the lower soil layer of 40 cm, and roots with a diameter of ca. 2 cm or higher were also detectable (Fig. 2) after in situ calibration (see subsection “GPR data processing”), as also shown in previous studies (Butnor et al. 2001; Barton and Montagu 2004; Rodríguez-Robles et al. 2017).

Fig. 2
figure 2

Example of root detection along a linear profile section of the trial (marked in Figure S1) using a Ground-Penetrating Radar (GPR). The GPR radargram was generated with a 800 MHz shielded antenna and band pass filters (both high-pass and low-pass filter) to eliminate low and high frequency noise (Rodríguez-Robles et al. 2017). The upper panel shows hyperbolically shaped reflections representing root reflections (dark purple hyperbolae). The lower panel illustrates the roots detected from the hyperbolae signatures (grey dots), and additional roots measured in situ for calibration purposes are depicted by black dots

GPR measurements were taken continuously along three linear profiles (P1, P2, P3) following the trial’s column direction and passing through the two central trees of each experimental unit (Fig. 3a, b). For every experimental unit, P1 and P3 were settled at ca. 0.30 m distance from the main trunk of the two central trees, while P2 was positioned halfway from both trees, that is, at ca. 1.25 m (Fig. 3c). In order to assign each detected root to a particular experimental unit, the position of trees along each linear profile was manually recorded in the RAMAC XV Monitor implemented in the GPR, along with GPR measurements. With this information, relevant data could be retrieved at sub-plot level in the profile (see “GPR data processing” subsection).

Fig. 3
figure 3

a Aerial image of the study site exemplifying the three geophysical GPR profiles (P1, P2, P3) used at plot level following the column direction of the common garden of P. halepensis. Measurements were done between the two central trees of every plot (experimental unit). b Detailed view of the geophysical GPR profiles: P1 and P3 are at a row distance of 0.30 m from one or another central tree, while the position of P2 is at mid-distance (1.25 m) between central trees. c Scheme (top view) showing the geophysical profiles and the three different sub-plot areas (S50, S100, S250) evaluated, where sub-plots correspond to column distances of ± 0.25 m, ± 0.50 m and ± 1.25 m, respectively, from the two central trees of an experimental unit

GPR data processing

Coarse roots were identified through hyperbolic signatures on radargrams. First, GPR radargrams were examined using RadExplorer v1.42 software (MALÅ Geoscience AB, Sweden) through visual inspection of the hyperbolae (“manual dataset”), which allowed to determine trends in the distribution of roots along the linear profiles.

Based on this preliminary information, a second dataset was created in RadExplorer using a built-in detection procedure along predefined sub-plots (“automatic dataset”), which was subjected to several pre-processing filter functions, in order to reduce noise and improve the identification of hyperbolae. Particularly, we applied background removal filters in order to remove parallel bands, often derived from ground surface reflection. Then we eliminated the initial current signal component (DC) from the GPR control unit (RAMAC XV Monitor and ProEx system – adjust signal position) to guarantee that depth estimates started from 0. Afterwards, a high/low bandpass filter was applied in order to remove unwanted noise. Stolt F-K migration was also used in order to improve GPR image and to correct for root position (Rodríguez-Robles et al. 2017; Barton and Montagu 2004). This procedure did not allow for a case-by-case assessment of hyperbolae; however, it was able to detect a larger number of roots than the visual inspection of radargrams while being free of its subjectivity, allowing for a more robust analysis.

The automatic dataset was generated for three different sub-plot sizes. In all cases, their cross-sections extended from P1 through P3 and, therefore, were centred at each experimental unit’s mid-point, but had different lengths of 50, 100 and 250 cm (hereafter S50, S100 and S250 sub-plots, respectively; Fig. 3c). Note that the largest sub-plot length (S250) was equivalent to the mean inter-plot distance (Fig. 3c). This strategy allowed to identify a particular sub-plot area maximising population differences in root traits, as described in the “Statistical analyses” subsection.

Although pre-processing filters were applied to the original GPR information, a calibration was needed since the depth and size of an object cannot be inferred directly by the GPR. Thus, we measured depth and diameter in situ of 12 coarse roots. Coarse roots were excavated through soil digging underneath GPR scan line profiles (P1, P2, P3), in different experimental units across the trial (Fig. S1). Root depth and diameter were then GPR-estimated through the so-called ‘time interval between zero crossing’ (Rodríguez-Robles et al. 2017; Guo et al. 2013a; Tanikawa et al. 2013). Time interval between zero crossing is the time elapsed between signal emission and registration by the shielded antenna, which is manually adjusted to time zero before processing the data to set the instant in which the radar signal leaves the antenna. A linear regression between root diameter measured in situ and time interval between zero crossing was applied, and the resulting equation used for calibration of root diameter (Rodríguez-Robles et al. 2017).

Statistical analyses

The definition of a proper sub-plot size is important for optimal GPR-based root detection in experimental trials. While using a too large sub-plot size could potentially capture coarse roots of trees from neighbouring experimental units, a too small sub-plot could underestimate the actual number of coarse roots for each unit. To handle this issue, the variability in root number (frequency) was evaluated as a function of the linear distance from the trunk along the aforementioned profiles using evenly spaced bins of 10 cm. Overall, a high number of roots was found at small and medium distances, followed by a progressive decrease in root number and a sudden increase after some value, hence denoting overlapping with roots from neighbouring trees (see Results, “Root detection using GPR” subsection). In addition, the ‘optimal’ sub-plot length (S50, S100 or S250) was explored through linear mixed-effects analysis of variance (ANOVA) fitted to the automatic dataset independently for each sub-plot size. The ANOVA included block, column and population as fixed effects, while row within block and column by block interaction were defined as random effects in order to account for systematic variability at the trial site. We retained as optimal sub-plot the one that maximised the F value for population differences in the ANOVAs for coarse root traits (diameter, depth, and frequency), which approximately matched the optimal distance previously identified in the manual dataset (Table 1, Table S3). The mixed-effects ANOVAs were further extended to partition the variability among populations into fixed ecotype and between-population within-ecotype effects, hence testing for ecotypic structure in root traits.

Table 1 Mixed-effects model analysis of variance for root diameter (cm), root depth (cm) and root frequency (number of roots) of 56 populations of Aleppo pine categorized into 16 ecotypes grown in a common garden in Altura (Spain)

As indirect validation of our phenotyping approach, we compared the information obtained through GPR measurements with estimates of the relative contribution of two consecutive soil layers (0–15 cm, upper soil; 15–40 cm, lower soil) to the total amount of water taken up by the different populations, as inferred through the analysis of the isotopic composition (δ18O and δ2H) of soil and xylem water in 2010 (data retrieved from Voltas et al. (2015)). Additionally, simple correlations were calculated between GPR-based root traits and xylem water δ18O at plot level.

Population (or ecotype) least squares means of root depth, diameter and frequency were subjected to linear and quadratic regressions as a function of long-term climate conditions of each population (or ecotype) at origin. The three populations of ecotype 16MC (reforestations of Northern Spain) were not used in the regression analysis because they have uncertain geographic origins. Linear and quadratic regressions of root traits were also performed on Euclidean geographic distances and climatic distances (Gower’s distances) of each population (or ecotype) from the trial site. Gower’s distances indicate dissimilarities among entities, in this case climate dissimilarities between populations (ecotypes) at origin and the trial site, and were calculated following Rutter and Fenster (2007):

$$\mathrm{GD}=\frac{1}{\mathrm{p}}\ {\sum}_{\mathrm{i}=1}^{\mathrm{p}}\frac{\left|{\mathrm{A}}_{\mathrm{i}}-{\mathrm{B}}_{\mathrm{i}}\right|}{{\mathrm{r}}_{\mathrm{i}}}$$
(3)

where p is the number of climate variables, Ai and Bi are the values of each climate variable (MAT, TAR, MAP, PsP) at site A (trial) and B (population or ecotype origin), and ri is the range of each climate variable in the dataset.

Finally, linear and quadratic regressions were calculated for root traits as a function of aboveground growth variables (H, DBH and crown area), both at population and ecotype level. The associations of aboveground growth variables with those climate variables at origin found to be related to population or ecotypic variability in root traits were also evaluated. For model selection, the best fitting model for each explanatory variable (linear, quadratic) was chosen considering its overall significance (F-value) and goodness of fit (adjusted R2).

Results

Root detection using GPR

Manual measurements of root diameter showed a positive correlation (r = 0.95; p < 0.001) with GPR records (Fig. 4a). Measured and GPR-estimated root depth also showed a positive association (r = 0.98; p < 0.001) (Fig. 4b). Following the visual inspection of radargrams, we detected the largest number of roots at distances between 10 and 50 cm from each tree trunk along the trial’s column direction, with a progressive decrease at larger distances (Fig. 4c). However, the number of detected roots increased again starting at 1 m apart from each tree. This observation suggested that GPR was detecting roots corresponding to trees from adjacent experimental units for distances higher than 90 cm from each targeted tree trunk (Fig. 4c). This was confirmed by analyses of variance carried out on the automatic dataset for three different sub-plot sizes (S50, S100, S250). Particularly, the results showed that the sub-plot size that maximised population differentiation in root traits was S100 (Tables 1 and S3). Further analyses were therefore performed for S100. At sub-plot level (S100), GPR estimates of coarse root traits varied between 3.0 cm and 5.0 cm (diameter), between 14.7 cm and 44.3 cm (depth), and between 3 and 19 roots (frequency).

Fig. 4
figure 4

Calibration of GPR mesurements. a Relationship between root diameter obtained in situ through soil digging and time interval between zero crossing (measured in nanoseconds; Δns) obtained using the GPR radargram (n = 12). b Relationship between root depth measured in situ and root depth inferred by GPR. c Total number of roots manually detected from GPR radagrams as a function of distance from the tree trunk

Phenotypic variation in root traits

We detected significant variation in coarse root diameter among populations (p = 0.030), but the effect was only marginally significant among ecotypes (p = 0.091; Table 1). Mean diameter values of ecotypes varied between 3.5 cm (15TU, Tunisia; Table S2) and 4.0 cm (3ALC, Spanish Southern Plateau). Differences among populations in coarse root depth (Table 1) were marginally significant (p = 0.086), while significant differences were detected among ecotypes (p = 0.011). Mean depth values of ecotypes varied between 21 cm (15TU) and 30 cm (16MC, reforestations of Northern Spain). These results suggested a stronger ecotypic structure for rooting depth than for root diameter. We did not observe significant differences in root frequency among populations (p = 0.460), whereas marginally significant differences were detected among ecotypes (p = 0.096). Mean frequency values of ecotypes varied between 8.9 (2MO, Spanish Ebro depression) and 13.1 (7BM, Spanish Betic system). Estimates of coarse root depth and diameter were significantly and positively correlated both at population (r = 0.42, p = 0.001) (Fig. S2a) and ecotype level (r = 0.63, p = 0.009), but root frequency did not correlate significantly with root depth or diameter (Fig. S2b, c).

The variability of root traits was in good agreement with xylem water isotopes records obtained in 2010 for the same trial, with the exception of root frequency. At population level, we observed significant negative correlations between coarse root diameter or depth and δ18O of xylem water in early autumn (r = 0.62, p < 0.001, and r = 0.50, p < 0.001, respectively; Fig. S3), but not in peak summer (mid-July) (p > 0.10, Fig. S3). At ecotype level, a negative association between root diameter and the relative contribution of the upper soil layer to xylem water in early autumn was also observed (Fig. 5a); this relationship was non-significant for root depth (Fig. 5b) and root frequency (Fig. 5c). There was also a negative association between the relative contribution of bottom soil layer in early autumn and root diameter at ecotype level (Fig. 5d). However, this association was non-significant for root depth (Fig. 5e) and root frequency (Fig. 5f).

Fig. 5
figure 5

Simple correlations between GPR-derived population means of root diameter (a, d), root depth (b, e) or root frequency (c, f) and xylem oxygen isotope composition (δ18O) obtained in mid-July (a, b, c) or late-September (d, e, f) of 56 populations of Aleppo pine tested in a common garden in Altura (Spain). Isotopic records were retrieved from Voltas et al. (2015)

Phenotypic associations with aboveground growth

At ecotype level, root diameter and root depth showed marginally significant quadratic relationships (p = 0.064 and p = 0.102, respectively) with tree height (Fig. 6, Table S4), while root frequency was unrelated to vertical growth. These results suggested that shallower and smaller roots were indistinctly found in ecotypes having lower or higher than average height growth. Other aboveground growth variables (DBH, crown area) were unrelated to root depth, diameter or frequency (Table S4).

Fig. 6
figure 6

Regressions of root diameter, depth and frequency as a function of tree height. Data correspond to 56 populations of Aleppo pine (grey dots) grouped into 16 ecotypes (numbers) (ecological region codes are defined in Table S2) and tested in a common garden in Altura (Spain). Each panel shows the regression (linear or quadratic) that better fits the data at ecotype level. Significant (p ≤ 0.05) and marginally significant (p ≤ 0.10) regressions are indicated with continuous and dashed lines, respectively. Non-significant regressions are also accompanied by the coefficient of determination (R2) and associated probability (linear case only)

Associations with geographic variables and climate at origin

Coarse root traits were related to some geographic and climatic variables at origin of populations, either linearly or quadratically (Table 2a). Root traits decreased linearly (albeit weakly) with geographic distance to the trial and also with longitude (Table 2a). These results suggested geographically structured differences between populations located near the trial site (i.e. the westernmost distribution of Aleppo pine) and those from the eastern Mediterranean basin. A significant quadratic dependence of root depth on Gower’s distance (Table 2a) was also observed, indicating that populations with climate at origin either most similar or most dissimilar to the trial site presented shallower roots.

Table 2 Coefficients of determination (R2) of linear and quadratic regressions of root traits (diameter, depth, frequency) at population (a) and ecotype level (b) as a function of several ecogeographical variables at origin of 53 populations grouped into 15 ecotypes of Aleppo pine tested in a common garden in Altura (Spain)

Root diameter showed a quadratic dependence on solar radiation (p = 0.014), with thicker roots occurring in populations having intermediate solar radiation at origin. In turn, root depth showed a quadratic dependence on PsP (p = 0.035) and, marginally, on MST (p = 0.067), with shallower roots observed under low and high values of both climate variables (Table 2a). Root frequency was linearly dependent on solar radiation (p ˂ 0.001) and marginally on PsP (p = 0.060), with more roots observed under high solar radiation and low PsP at origin. In addition, we also detected significant or marginally significant positive linear dependencies of root frequency on VPDs (p = 0.003) and TAR (p = 0.073). Although these relations hinted to the dependence of roots traits on climate conditions of populations at origin, the proportion of explained variance was relatively low (R2 < 0.25).

Similar dependencies of coarse root traits on geographic and climatic factors occurred at ecotype level, with relationships being either linear or quadratic (Table 2b). In particular, we observed a negative linear dependency of root diameter on longitude and a positive linear association of root frequency with latitude. There were also marginally significant quadratic relationships between root depth and latitude and between both root depth or diameter and altitude (Fig. S4, Table 2b). The quadratic relationships were such that northern and southern ecotypes tended to present shallower roots, while deeper and thicker roots were distinctive of ecotypes from either low or high altitudes. Based on Gower’s distances, ecotypes having high or low climate similarities with the trial climate showed shallower and finer roots, but they were not distinctive with regard root frequency (Fig. S4).

Root frequency showed a negative linear dependence on PsP (Table 2b), while quadratic dependencies of coarse root traits on MST were detected for root depth and diameter, but not for root frequency (Fig. 7a,b,c). Both root depth and diameter showed also a quadratic relationship with solar radiation, while root frequency showed a positive relationship with this variable (Fig. 7d,e,f). These results suggested that shallower and finer roots were present in ecotypes having either relatively low and high MST or solar radiation values at origin. Also, root depth showed a marginally significant (negative) linear dependence on VPDs (Fig. 7e,f), denoting that ecotypes subjected to higher summer transpirative demand tended to show shallower roots. As for root depth, tree height showed a negative dependence on VPDs (p = 0.029). On the contrary, tree height was negatively related to SR (p = 0.035) and unrelated to MST (results not shown).

Fig. 7
figure 7

Regressions of root diameter, depth and frequency as a function of climate variables at the origin of populations: (a, b, c) mean summer temperature (MST), (d, e, f) solar radiation (SR), (g, h, i) summer vapour pressure deficit (VDPs). Top panels depict regressions for root diameter, central panels for root depth and bottom panels for root frequency. Data correspond to 53 populations of Aleppo pine (grey dots) grouped into 15 ecotypes (numbers) (ecological region codes are defined in Table S2) and tested in a common garden in Altura (Spain). Ecotype 16 (and three associated populations) corresponds to reforestations of uncertain geographic origin of the Northern Spanish Plateau and has not been included in the analysis. Each panel shows the regression (linear or quadratic) that better fits the data at ecotype level. Significant (p ≤ 0.05) and marginally significant (p ≤ 0.10) regressions are indicated with continuous and dashed lines, respectively. Non-significant regressions are also accompanied by the coefficient of determination (R2) and associated probability (linear case only)

Discussion

GPR applicability to root phenotyping of forest trees

This study showcases the use of GPR devices as high-throughput phenotyping tools to assess intraspecific differentiation in coarse roots potentially related to water use and drought tolerance strategies of forest trees. It broadens the customary use of GPR for tree root detection in ecological studies carried out since the turn of this century (Guo et al. 2013a, b; Hruska et al. 1999). By using an 800 MHz antenna we were able to detect 1,835 roots with diameters equal to or above 2.2 cm at a maximum depth of 88 cm for about 330 trees monitored in a common garden of ca. 0.8 ha.

In situ calibration provided direct evidence of the adequacy of GPR measurements for root detection. An additional (indirect) indication of the potential of GPR for root detection was the observation that the variability detected in coarse root traits for Aleppo pine agreed (at least partly) with the information derived from xylem water isotopes on the use of water sources by the same trees (Voltas et al. 2015). However, we should note that only associations involving δ18O of xylem water in early autumn were significant, whereas no relevant relationships with GPR records were detected in peak summer. This suggests that the existing electromagnetic gradient between roots and soil at the time of GPR measurements (i.e. mid-June) could have identified more roots than those effectively taking up water during the acute summer drought typical of end-July. For an isohydric species such as Aleppo pine, a fraction of coarse roots might have remained in a quiescent condition during soil dry-down in summer as avoidance strategy against hydraulic failure. This condition might be the result of a reversible development of water transport barriers such as suberin or cutin, stimulated by nutrient scarcity and drought stress (Vandeleur et al. 2009; Schreiber 2010; Barberon et al. 2016). Roots might be also subjected to seasonal embolism, since partial cavitation in roots is less deleterious than xylem failure because roots can recover after drought by refilling the embolized conduit (Domec et al. 2004). Changes in root turnover and turgor are some of the prevailing consequences of high soil temperature and water deficit stress in plants (Brunner et al. 2015; Gill and Jackson 2000), and both have been widely described for the case of fine roots (Gill and Jackson 2000; Kitajima et al. 2010; Montagnoli et al. 2019). However, there is no evidence of inactivation of coarse roots in forest tree species thus far. This possibility would deserve detailed investigation through e.g. temporal assessment of GPR signals following changes in soil water status during the growing season.

Evidence for intraspecific variability in root traits of Aleppo pine

It is known that root system architecture varies among Iberian pines, with mountain species (e.g. Pinus sylvestris) usually having shallower roots than typical Mediterranean species (e.g. Pinus pinaster, P. halepensis) (Andivia et al. 2019). This variability correlates with the species’ ecological niches and the intensity of drought stress present in their environments (Andivia et al. 2019). However, there is still a lack of information on the variability of functional root traits at the intraspecific level, which determines the efficiency of water acquisition by trees (Kirfel et al. 2017). A main goal of this work was to fill this knowledge gap for a Mediterranean pine species and, particularly, to characterize intraspecific patterns of coarse root traits potentially related to different strategies of water uptake and use in Aleppo pine.

We found population differentiation for root traits in P. halepensis, which was geographically structured following an ecotypic pattern. Such patterns have been previously reported in P. halepensis for traits related to the tree’s water budget such as water-use efficiency (Voltas et al. 2008), transpiration (Santini et al. 2019b) and the use of water sources (Voltas et al. 2015). Also, we observed greater variability in coarse root traits among populations originating from the Iberian Peninsula than among their eastern Mediterranean counterparts in spite of the loss of genetic diversity that Aleppo pine suffered after the post-glacial recolonization from the eastern Mediterranean (Grivet et al. 2009). This result, however, requires confirmation because comparatively less populations were available from the eastern Mediterranean compared to the western Mediterranean basin, being also distributed within a narrower latitudinal gradient.

Relationships between root traits and aboveground growth

We found quadratic associations of coarse root depth with tree height and diameter at the intraspecific level. Although these relationships were only marginally significant, they suggest predictable allometric patterns linking below and aboveground growth, which are partially consistent with the general allometric scaling theory (West et al. 1999). In particular, proportionality of carbon investment to above and belowground growth is suggested for ecotypes reaching up to ca. 6 m height at age 19 years in the study site (root to shoot biomass is known to stabilize once pines reach the reproductive stage; Peichl and Arain (2007), Cao et al. (2012)). Above this height threshold, a negative association between aerial growth and investment in roots could be observed at the ecotype level. This suggests a preferential mass allocation to aerial carbon stocks for mesic ecotypes having high growth potential (Climent et al. 2008; Voltas et al. 2018), as exemplified by ecotypes 12GR and 14IT originating from the center-eastern Mediterranean basin. These ecotypes experience mild climate conditions at origin (either having MAP > 500 mm or MSP > 60 mm) and are known to exhibit a larger tree height plasticity compared with xeric ecotypes across water availability gradients (Patsiou et al. 2020). This superior plasticity might be associated with a progressively large carbon allocation to roots under harsher (i.e. drier or nutrient-poor) conditions. In any case, our results suggest that these mesic ecotypes allocate relatively more resources to aerial growth than to root development under the relatively favorable conditions encountered in the trial, as compared with xeric ecotypes. This result partially agrees with previous studies stating that forest trees allocate carbon preferentially to the root system as a response to poor nutrient availability (Hermans et al. 2006; Vicca et al. 2012), which precedes synthesis of secondary metabolites (Prescott et al. 2020).

Variability in root traits of Aleppo pine follows geographical and climatic gradients

The assessment of geographic dependencies of roots traits suggests that root diameter decreases across an eastward longitudinal cline in Aleppo pine populations. This cline agrees with previous results being indirectly informative of allocation patterns to roots in P. halepensis, which reported decreases in reproductive investment (Climent et al. 2008) and water-use efficiency (Voltas et al. 2008) in eastern populations. On the other hand, the larger proportion of shallow roots in populations originating from either extreme of the species’ latitudinal range disagrees with our initial hypothesis. However, root frequency decreased as latitude increased, suggesting that southern populations, which generally are exposed to a drier climate, have more abundance of coarse roots than northern populations, regardless of their root depth. This discrepancy could be related to the higher costs of construction and maintenance of deep roots (Schenk 2008a) in ecotypes with low growth potential and high reproductive allotment, leading to carbon allocation shifts towards shallower and more abundant roots in dry and resource-poor habitats. This could be an advantageous strategy in habitats where the uppermost soil layers may present higher nutrients and moisture as a result of discrete precipitation pulses (Schenk 2008a, b).

Previous studies have shown that rooting depth increases with water shortage in pine species (Andivia et al. 2019), an interspecific pattern that only partially agrees with our results at the intraspecific level. In fact, only root frequency increased with aridity at origin of populations of Aleppo pine. On the other hand, intraspecific relationships between root depth or diameter and temperature were hump-shaped, and peaked towards the thermal midpoint of the species distribution range (Albert et al. 2010). This finding may be related to the existence of different adaptive strategies imprinted in a number of life-history traits and their potential trade-offs – in addition to rooting traits – to cope with the conditions encountered by the species along its distribution range, as modulated by stresses (drought, cold temperatures, pests) and disturbances (fire). For example, the associations between aerial growth or rooting traits and incoming radiation might indicate allocation trade-offs related to higher selective pressures for increased competition for light in the case of Greek and Italian mesic ecotypes. This result also suggests that the existence of trade-offs between reproduction, defense and vegetative growth in Aleppo pine (Climent et al. 2008; Santini et al. 2019a; Voltas et al. 2008) may imply a decrease of C allocation to coarse roots in fire-prone, dry and warm habitats. The negative dependency of rooting depth on VPDs has to be carefully interpreted, however, since this relation was mostly driven by the Tunisian ecotype (15TU), which has shallow roots despite being exposed to high summer evapotranspiration. This could be explained as the result of conflicting functional strategies typical of stressed and disturbed environments (Santini et al. 2019a), implying a relatively low allocation to growth (both aerial and belowground) compared with reproduction (Climent et al. 2008; Santos-del-Blanco et al. 2013). Such possibility is supported by the negative intraspecific association found between tree height and VPDs.

Methodological limitations

Despite its high detection capability and non-destructive assessment of root features, GPR is not devoid of technical limitations for tree phenotyping studies. These are related to the existing variability in root orientation, since overlapping roots and roots underneath the tree trunk are underestimated by the device (Butnor et al. 2016; Li et al. 2016). However, GPR have been recognized as a proxy tool sufficiently precise to estimate coarse root frequency (Butnor et al. 2003; Hirano et al. 2009, 2012; Guo et al. 2013a). Also, the existing electromagnetic gradient between roots and soil makes a previous calibration strictly necessary to reduce the background noise produced by the particular physicochemical characteristics of the soil. In particular, it should be noted that GPR detection capacity is seriously impaired in both nearly saturated (Hirano et al. 2009; Rodríguez-Robles et al 2017) and very dry soils (Hirano et al. 2009). In this last case, GPR detection is limited if roots have low water content (e.g. under 20% volumetric water content in Cryptomeria japonica; Hirano et al. 2009). In this regard, our field campaign targeted a period of full vegetative activity of trees while avoiding such extremes in soil water status (e.g. abundant precipitation in April-early May and acute drought in July–August).

Conclusions

This study provides new avenues for the examination of intraspecific variability in root traits of forest tree species using GPR as high-throughput phenotyping method. We observed population differentiation in P. halepensis following east–west (root diameter) and north–south (root frequency) geographical patterns. Conversely, the interspecific pattern described for the genus Pinus, by which species originating from xeric conditions exhibit large allocation to roots, was not strictly found at the intraspecific level for Aleppo pine. This observation may be attributable to the existence of different adaptive strategies under varying environmental conditions for the species, which likely lead to trade-offs involving the use of resources. In this regard, our work contributes to unravel the evolutionary complexity of a widespread conifer of high ecological significance in the Mediterranean basin.