Phenotypic plasticity of leaf shape in selected and semi-domesticated genotypes as new tool of Argania spinosa L. Skeels breeding

© Slovak University of Agriculture in Nitra Faculty of Agrobiology and Food Resources


Introduction
Changes in plant traits are signals of plant adaptation in different environments. It was reported before that plant traits are determined by a combination of genetic heritage (genotype) and the capacity for responding to environmental changes within a single genotype (phenotypic plasticity) (Royer et al., 2009). In fact, plants react to these environmental changes by biomass allocating among several plant organs to capture optimum light, water, nutrients and carbon dioxide, as strategy to maximize growth rate (Bloom et al., 1985).
This adaptation type often refers to phenotypic plasticity, defined as a change in the phenotype expressed by a single genotype in different environments (Gratani, 2014). These plasticity responses are expressed at different levels such as plant morphology, anatomy, physiology and growth (Pigliucci, 2001). Thus, the understanding of phenotypic plasticity will be crucial for predicting changes in species distribution, community composition and crop productivity under global change conditions (Van Kleunen and Fischer 2007;Lande, 2009).

Phenotypic plasticity of leaf shape in selected and semi-domesticated genotypes as new tool of Argania spinosa L. Skeels breeding
In response to climate factors, most of studies on leaf have included species occupying different environments. The results reported revealed a remarkable leaves morphological variation associated to their environment (Joel et al., 1994;Bruschi et al., 2003;Uribe-Salas et al., 2008). Therefore, the variation in leaf traits plays a key role to understand the ecological and biological evolution. Some of them are considered as indicators to select species with excellent growth performance and high biomass productivity (Bunn et al., 2004).
Argan tree [Argania spinosa (L.) Skeels] is one of species characterized by a wide phenotypic variation. It is a tree endemically marked the south-western region of Morocco and semi-desert regions of the south-west of Algeria, belonging to the monospecific genus Argania. It is the only representative of the tropical family Sapotaceae in Morocco (Pennington, 1991). It is a multipurpose forest tree (forest-fruit-forage), highly important for the country in terms of biological and ecological biodiversity as well as economic and social development (M'hirit et al., 1998;El Aich et al., 2007). It presents high medicinal, cosmetic and therapeutic interests, due to its extracted fruits and leaves oil. In another hand, studies on morphology and leaf variability are very limited up to date. The first descriptions showed that leaves are entire, leathery, alternate and often united in fascicles (Sauvage and Vindt, 1952). They are more or less stalked, glabrous and midrib clearly marked. Leaves shape of young seedlings are mostly elliptical (Cornu, 1897). Nonetheless, the adult tree present lanceolate, oblong-lanceolate, spatulate, obovate or oblanceolate shape (Perrot, 1907;Rieuf, 1962;El Aboudi, 1990;Prendergast and Walker, 1992;Zahidi, 2013).On the same tree, two leaf types coexist. Simple leaves appear alternately on the stem of the year or on portion of growing shoot, and grouped leaves on woody shoots (Rieuf, 1962;M' Hirit, 1989).
To test the important of morphological plasticity in predicting population dynamics and plant adaptations to different novel environments, we used Argania spinosa leaves to: a) Determine the morphological variability of leaf shape and size between selected genotypes under different climatic conditions in three locations of Morocco, b) Identify the response of the studied genotypes under different contrasted environments.

Study sites, plant material and experimental establishment
To study the impact of the environment on leaves, trials were conducted in three experimental locations belonging to the National Institute for Agricultural Research, located in different regions in Morocco and characterized by distinct environmental conditions (Table 1 and Figure 1): a) South-Western experimental location: latitude 30º 30' N, longitude 9º 33' W, altitude 80 m a. s. l. The soil is alkaline (pH = 8.45) with clay texture. An arid climate predominate this region, with a defined dry season, mean temperature varies between 11 and 28.9 °C and annual mean rainfall around 290 mm. b) North-Western experimental location: latitude 35º 08' N, longitude 6º 08' W, altitude 39 m a. s. l., characterized by sandy and alkaline (pH = 8.1) soil. The climate is subhumid, with a mild winter period. Mean temperatures vary from 9 to 28.3 ºC. Annual mean rainfall is around 664 mm. c) Central experimental location: latitude 33º 56' N, longitude 5º 12' W, altitude 550 m a. s. l. The soil is clay and neutral. The climate is warm Mediterranean with dry summer. Mean monthly temperature varies from 2.8 to 37 ºC and annual mean rainfall ranges around 470 mm.
Parental trees of studied Argania spinosa genotypes were selected and sampled in 2014 from forest stand in South-Western region (30° 15' N, 009° 21' W). Twenty trees were selected on the basis of fruit yield and morphological criteria (Table 2). Fruits were collected from each mother tree upon homogeneity of shape, size and number of contained seeds. Obtained seeds from peeled fruits were germinated. Six months after germination, a total of sixty seedlings (20 half-sibling) were used for field trials.
The transplantation was carried out in October 2014. One progeny seedling from each of 20 selected mother trees was performed in the three locations according to the following experiment protocol: The seedlings were coded from A1 to A20 for Central location, L1 to L20 for North-Western location and M1 to M20 for South-Western location. Afterwards, the genotypes were transplanted in 7 × 7 m plot with four lines containing five plants.
Plant maintenance was conducted including irrigation, fertilization, optimal control of pest and disease and intensive weeding.

Data sampling and landmarks
In this study, leaves were collected in spring (April) as their formation starts in autumn (October) and reaches maturity in February. Five leaves per genotype were sampled in April 2018 from the 20 selected trees in each experimental study site (100 leaves in total per study site). Leaves from each tree were sampled at a height of 0.5-1.5 m around the crown from primary branches of each tree ( Figure 2). A total of 300 leaves were characterized morphologically, dried and scanned using HP scanjet G4010 scanner. The procedure of digitization was repeated twice to estimate digitizing error. Scanned images were used to record landmarks on the left and right sides of each leaf. The shape of the whole leaf was characterized by a set of 61 landmarks that indicate the relative sizes and approximate shapes. Then, for each leaf, 7 landmarks were fixed and were located as follows: lm 1 and 7, beginning of the petiole; lm 2 and 6, junction of the blade with petiole; lm 3-5, the largest width of the leaf blade; lm 4, apex of the leaf blade. The remaining points (54 landmarks) were placed between the fixed points to have more precision data for size and shape ( Figure 3).

Table 2
Morphological and yield criteria of the selected tree

Criteria Description
Morphology general aspect of tree: tree height, first leaf height, trunk diameter, number of branches, display of tree, growth habit leaves: leaf length, leaf width, leaf shape, leaf color fruits: fruit length, fruit width, fruit shape, fruit color nuts: nut length, nut width, nut shape Yield 1,000 fruit weight, mother tree production

Figure 1
Location map of the experimental sites Leaves are bilaterally symmetric above their respective central axes (petiole⁄midrib). The symmetry shape analysis method was used to take in consideration this symmetry and to generate information about shape asymmetry (Klingenberg et al., 2002). This method was based on the Procrustes superimposition of the original and reflected copies of each landmark configuration (Dryden and Mardia, 1998). A symmetric component of shape variation was obtained from the averages of original and reflected copies. An asymmetric component was computed from the differences between original and reflected copies (Klingenberg et al., 2002).
To evaluate variations in leaves size and shape, a Procrustes ANOVA (Klingenberg and McIntyre, 1998) were conducted with centroid size and shape used as independent variables. Shape coordinates were computed by standardizing each configuration to unit centroid size and by minimizing differences in translation and rotation of all individuals using, a leastsquare algorithm. Size was measured for each individual as the centroid size of the landmark configuration. The centroid size measures the dispersion of landmarks using a function of their distances from the centroid, which is the 'baricenter' of a configuration (Rohlf, 2000).
Therefore, the effect of argan tree genotypes represents shape individual variation, and the effect of body side represents leaf asymmetry. Analysis was carried out using MorphoJ 1.06d software.
To quantify the similarity of covariance structures of leaves shape, at the different levels, matrix correlations were computed between the corresponding covariance matrices. Matrix correlation is a procedure for making a quick overall comparison of two covariance matrices (Klingenberg, 2011). To test the association statistically, a matrix permutation test was used to evaluate correlation matrix against null hypothesis of no relationship between

Figure 3
Global schema of landmarks measured for the whole leaves by TPSDig2W3

Software applications and statistical analysis
To characterize leaves shape and size software from the TPS Series was used to create files for data manipulation and digitization, with set of programs for two dimensional geometric morphometric analysis (Rohlf, 2010): Firstly TPS file was created in TPSUtil32 version 1.74 with the extension TPS. After the TPS file was opened in TPSDig2W32 version 2.30 to digitize the 61 landmarks in the same order on each leaf, after setting a scale factor. The TPS was converted into NTS file to rearrange a data as a matrix with rows corresponding to samples and columns corresponding to coordinates. the matrices, by permuting the order of landmarks. The matrix correlations between-and within location levels were computed; similarity was analyzed between paired locations for symmetric and asymmetric covariance structures.
The amount of genotypes variation and the shape variation in the entire dataset was treated using Principal Component Analysis (PCA) based on the covariance matrix of the leaf shape variation.
Analysis for environment discrimination was conducted on the average data of each region. Canonical Variation Analysis (CVA) was also performed in MorphoJ to test group differences and to plot their differences. For CVA, the statistical significance of pairwise differences in mean shapes was assessed with permutation tests (10 000 permutations per test) using the Mahalanobis distance (Iwaizumi et al., 1997).
Thin-plate spline interpolation functions were used to at each end of the axes to compute deformation grids and to represent shape variation across genotypes.

Procrustes ANOVA
Effects of genotype (effect of each tree individual), half-sibling and locations were estimated on geometric size using the isotropic model. The results for size and shape of the sixty genotypes are reported in Table 3. Concerning the size, the genotype effect was statistically significant and explained about 68% of the total variance. Furthermore, differences in leaf size among locations and among half-sibling were highly significant (17% vs. 14% of the variance). No significant variation was found for leaves within the same genotype.
Procrustes ANOVA for shape variation showed that genotype was involved in large part of the total variance (ca 60%). Half-sibling and location were statically significant explaining around 22% and 2% of the variance, respectively. The "side" factor, indicating a directional asymmetry, also showed statistical significance at all three biosystemic levels (p <0.0001).
Interaction factors were statistically significant as at the location level, at the tree level as a genotype, and at the leaf level as a photosynthetic organ (in all cases p <0.0001). Finally, fluctuating asymmetry was significant, explaining a 10% for "half-sibling × side effect", ca 2.99% for "genotypes × side effect" and 0.2% for "location × side effect".
In summary, for the shape parameter, the differences between genotypes were statistically significant. The variation of leaf shape was high due to genotype rather than location effect. Procrustes ANOVA within locations showed similar patterns of variance partition for the three regions (Table  4). For leaf size in the three regions, 99% of the variance was explained by the genotype effect.
Regarding a shape, the three sites showed almost the same pattern of variation, which was mainly explained by the genotype factor (79-82%). Variation among leaves (0.1%) was lower than the fluctuating asymmetry (13-15%). Thus, according to each region, directional asymmetry was significant and digitizing error explained a negligible part of the total variance.

Correlations of covariance structures
The global similarity of the covariance structures inside and among regions was measured for sixty genotypes/ individuals. Then, the matrix correlations were computed at inter and intra-location levels between the covariance matrix of fluctuating asymmetry and symmetric component. For the entire dataset, this matrix correlation was 0.12 (P = 0.01). Within each site, the matrix correlations between fluctuating asymmetry and symmetric variation were significant: 0.1 (P = 0.0089) for South-Western location, 0.08 (P = 0.0605) for North-Western location and 0.14 (P = 0.0031) for Central location.
Moreover, for the paired sites study, great and significant (P <0.0001) matrix correlations were detected between the covariance structures for the symmetric factor of variation (South-Western location vs. North-Western location was 0.839; South-Western location vs. central location was 0.651; North-Western location vs. central location was 0.720). The covariance structures of fluctuating asymmetry showed a strong correlation among locations (South-Western location vs. North-Western location reached 0.983; South-Western location vs. Central location reached 0.918; North-Western location vs. central location reached 0.940).

Principal Component Analysis
It was noticed from ANOVA results that, the variability within same genotype is lower than variability between genotypes. The averaged data of leaves per genotype was used to further evaluation of region effect.
The morphological variation corresponding to the first two PCs for each component was represented by the corresponding scatter plots and by the TPS-grids representing the extremes of shape variation along these axes (Figures 4 and 5). Eigenvalues from each PCA reveal the shape features contributing to shape variance along each principal component (PC). Deformation grids for the most negative and positive extremes of the PCs reveal a large variation in the sample for the variation components.

Figure 4
Ordination plots and leaf shape variation of the first two principal components (PCs) obtained from the principal component analyses computed on symmetric component; confidence ellipse at 95%. PC1 and PC2 presented 89% of the total variance. For each PC, the left and right diagrams represent the shapes for PC scores of -and +

Figure 5
Ordination plots and leaf shape variation of the first two principal components (PCs) obtained from the principal component analyses computed on asymmetric component; confidence ellipse at 95%. PC1 and PC2 presented 61.93%% of the total variance. For each PC, the left and right diagrams represent the shapes for PC scores of -and +

Figure 6
Ordination plots of Canonical Variates extracted by Canonical Variate Analysis (CVA) computed on the symmetric component; confidence ellipse at 95%. CV1 and CV2 explaining 91.65% of total variance For symmetric component, PC1 presented 80.12% of the total variance. From negative to positive value, the leaf shape varied from obovate to ovate and the relative length of the petiole decreases. In addition, apical and basal parts were strictly contracted for negative values of PC1. PC2 accounted 9.7% of the total variance concerning changes in the relative length of the petiole and the shape of the basal and apical leaf parts. Results showed that the three sites mostly overlap with a tendency of ovate leaves for the Central location (Figure 4).
The percentage values of PCA for fluctuating asymmetry ( Figure 5) represent the level of variability in the data. 61.93% of overall variation was explained by PC1 (36.73%) and PC2 (25.20%). Resulting scatter plot showed a dispersion of genotypes. PC1 revealed the asymmetry caused by bending of the leaf blade with respect to the left or right orientation of the apical and basal parts. PC2 principally explained the asymmetry related to the relative distal or proximal positions of both left and right sides of the leaf. These involved differences in the relative sizes of the left and right leaf sides. The distribution of the sixty genotypes was almost equal around the central axis. PCA of the asymmetry covariance matrix did not show differentiation of genotypes and location in terms of asymmetry.

Canonical variation analysis
Concerning symmetric component, Mahalanobis distances among locations were 35.53 between South-Western location vs. North-Western location, 21.88 between North-Western location vs. Central location and 18.29 between South-Western location vs. Central location. The shape space occupied by each location was presented by projecting the scores onto the first two CV axes. The CVA indicated that each location was distinct from the others. In fact, the Canonical Variate CV1 explaining 91.65% of total variance, separated each region apart. CV2 (8.34%) differentiates the Central location from South-Western and North-Western ( Figure 6).
It was also showed that, leaves of the Central location genotypes were distinct from those of South-Western location and North-Western location by the strictly oval leaf contour and a short petiole. The genotypes of South-Western location were distinguished by a large petiole, ovate shape and obtuse apex.
For asymmetric component, the Mahalanobis distances were 8.86 between South-Western location vs. North-Western location, 25.62 between North-Western location vs. Central location and 19.40 between South-Western location vs. Central location. All permutation tests indicated that the mean shapes differed significantly between locations (all P <0.0001 in pairwise permutation tests between locations). Outline reconstructions with asymmetric datasets shows a left-right basal asymmetry on CV2 (4.95% of overall variation). The other CVs revealed no observable variations (Figure 7). In another hand, the multivariate analysis was focused on symmetric dataset, because the asymmetric component made little contribution to total leaf outline variation.

Genotypes × Environment interaction of leaf morphology
Morphometric analysis of twenty half-sibs leaves studied under three contrasted environments showed significant differences between genotypes and locations. According to grids deformation describing leaf shape and size

Figure 7
Extreme leaf shape reconstructions using the Warped Outline Drawing along the first two CVs from the asymmetric data changes, three principal groups were observed; each group corresponds to leaves of genotypes in the same location. Leaves morphological aspect of genotypes reflects their response to environment. Indeed, North-Western genotypes were characterized by leaves with a large blade, short petiole, obtuse apex and elliptical shape. However, narrow ovate leaves and sub-obtuse apex were observed in Central location. Concerning South-Western, genotypes presented large leaves, long petiole, obtuse apex and ovate shape. The results indicate a strong variability and plasticity among genotypes from the same mother under different environment.
Leaves are the most important organ for synthesis and production of assimilates during life cycle of plants. They are responsible for regulating final plant growth and yield (Coleman, 1994). In the present research, mature leaf variability was investigated to determine the response of selected genotypes to their environment. Basically, in argan tree, leaf formation begins after the first rains, the observation of our genotypes in the three localities showed that the appearance of first buds started in autumn precisely in September at North-Western and Central locations and in October at South-Western location. Leaves achieved their maturity after around 60 days at North-Western and South-Western location and 125 days at Central location.
Morphometric results show that genotype generates a large difference in studied leaf parameters than environment effect. In fact, there is a significant genetic basis to the variation in leaf shape and size in Argania spinosa, the morphological response to the environment generally overrides the genetic influence. Natural variation in leaf shape has been described as being polygenic, mostly additive, and highly heritable, which explain the important effect of genotype (Langlade et al., 2005;Tian et al., 2011;Chitwood et al., 2013). The same result was observed at other woody species, the foliar morphological variation of plants growing in a different state and several species such as Quercus petraea (Viscosi and Cardini, 2011), Azadirachta indica (Kunduand and Tigerstedt, 1997), Eucalyptus sideroxylon (Warren et al., 2005) and Quercus rugosa (Uribe-Salas et al., 2008) have been studied. These research works showed that the genotypic effect was really compatible with leaf variations under various conditions. These differences between genotypes could be explained by influence of environment condition on expression of leaf genes and also by the considerable genetic variation observed due to argan tree allogamy.
Great part of genotypes from South-Western location was characterized by an ovate shape. However, the same genotypes in the other regions presented a very wide spectrum of variation and their forms vary from obovate to ovate. Genotypes L04, M12, A18, M10 and L12 are specific, close to each other and far from other individuals. They have very narrow leaves and an ovate shape. The morphological variations in leaves probably related to genotypic difference in the case of Argania spinosa trees. Pyakurel and Wang (2013) reported that the genotypic variation and environmentally induced variation in leaf morphological characteristics in Betula papyrifera Marsh. The similar pattern was also reported before by Possen et al. (2014) who studied the leaf morphology and the specific leaf area for 15 genotypes of Betula pendula. They concluded that the presence of the leaves diversity was related to the genetic diversity presence. In the present research study, leaves variation may be due to phenotypic plasticity explained by the interaction of genes and environment. Furthermore, the results indicate a strong variability and plasticity among certain genotypes from the same mother tree. For example, the genotype G10 ( Figure 8) showed various aspects in response of different environment. At the Central location, the leaves (M10) were small and their shape was strictly obovate. However, they were larger and exhibit ovate shape in South-Western location and elliptical in North-Western location.

Figure 8
Morphological variation of genotype 10 leaves in the three studied regions It is notice here that CVA was used to illustrate the strength of the separation of the different ecological groups in shape space, and to detect which morphological changes define the different groups. With three regions, CVA deduced of leaf shapes shows that most sites are clearly separated from each other. The separation of the groups with the CVs evokes that there is a significant environment effect on the studied genotypes. In fact, each station was characterized by distinct environmental conditions (temperature, soil, altitude etc.) which act differently on the leaves of selected trees. This result is in agreement with those of Li et al. (2015) on Ziziphus jujuba Mill., indeed, some important characteristics of jujube leaf morphology are related to climatic factors such as with mean annual temperature and mean annual precipitation. It is found that, with humid regions, Chinese jujube tends to have higher leaf area and leaf perimeter than varieties with drought stress. A similar influence of the environment factors was reported on Populus nigra L. the comparison of nine Populus nigra metapopulations in two experiment sites revealed that leaf morphology and leaf structure exhibit substantial genetic variation and phenotypic plasticity which both represent potentially significant determinants of populations ability to respond to environmental variations (Guet et al., 2015). The influence of different environmental variables seems to play a role in the variation for CVA, for example, leaves of North-Western location were characterized by a large blade, short petiole, obtuse apex and elliptical shape and those of South-Western were presented a large blade, long petiole, obtuse apex and ovate shape. Variability between those two sites lies mainly in the length of the petioles. Indeed, the high length of petiole observed in South-Western seems mainly due to the high temperature in this area. This result concurred with the findings of Ibañez et al. (2017) on Arabidopsis thaliana and Hajlaoui et al. (2007) on Fragaria × ananassa Duch who reported that length of the petiole in these two species increases with temperature rise. Genotypes of the Central location were characterized by a narrow blade compared to other stations. The reduction of leaves size was probably affected by the highest altitude (550 m a. s. l.) and the lowest temperature (2.8 °C). This result is consistent with earlier studies on Eucalyptus pauciflora, Acer rubrum and Quercus kelloggii revealed that leaf size decreased with elevation (Royer et al., 2008;Giles Yong et al., 2017). The effect considered to be particularly attributable to lower air temperatures and the photosynthetic process. It's also reported that, the low temperature acts on photosynthesis via both stomatal and non-stomatal mechanisms. This low temperature associated with stomatal restriction of intercellular CO 2 partial pressure and Rubisco inactivation above 15 °C (Hendrickson et al., 2004). Moreover, soil pH and texture significantly influenced leaf morphology, the result of this study showed that the Central location was characterized by clay and neutral soil (pH = 7.6) in contrary to South-Western and North-Western locations whose soils was sandy and alkaline (pH = 8.1-8.5). The soil pH and texture affect its fertility and had an effect on the mobility of mineral elements particularly P and N which known play a major role in physiological processes and photosynthesis. Furthermore, higher densities of stomata, guard and subsidiary cells recorded in clay soils could be as a result of the ability of these soils to retain water and nutrients as compared to other soils with higher sand proportions which is characterized by larger pores that drain water and leach nutrients (Bengough et al., 2011). The ability of clay soil to retain water for plant uptake in comparison with other soil types could be responsible for increased physiological processes and biomass accumulation (Marenco et al., 2006).
The phenotypic variation of the twenty genotypes was originated from the physiological and morphological plasticity which play different roles in response to environmental changes. Plasticity is a kind of adaptation developed by plants for acquisition of resources (Crickand and Grime, 1987), variations in size and placement of resource acquiring organs such as leaves and roots are of major importance for plant adjustment to resource availability (Dong et al., 1996). The higher leaf plasticity of argan tree can be advantageous in response to varying environmental conditions. Leaf shape and size is highly variable and can be used for identifying genotypes, developmental patterning within and among individuals, assessing plant health, and measuring environmental impacts on tree phenotype.

Conclusions
In conclusion, studied genotypes response under different contrasted environments varies according to genotype and location. Therefore, leaves of the twenty genotypes in the Central location seem to be small and characterized by a narrow blade compared to other locations. However, South-Western location leaves were characterized by an ovate shape and the same genotypes in the other regions presented a very wide spectrum of variation and their forms vary from obovate to ovate. Genotype and environmental variation are important indicators of morphological variability, playing the plants adaptation role to heterogeneous conditions. Furthermore, the geometric morphometric approach produced a clear and observable morphological variation for the sixty studied argan genotypes, being particularly suitable for the detection of shape differences that may reflect the adaptation of argan tree to different locations. This may provide novel insights on the potential adaptation responses of the forest ecosystem to climate change.