DiscussionTop

The density of young trees was higher than that of adults in all except one plot and reached 961.9 trees ha^-1 in plot 3 (Table 1). The frequency distribution of tree dbh classes showed an overall J-reverse shape reflecting a good regeneration capacity with a prompt decline in the abundance of larger-size trees (Fig. 2). Similar results were found by Linares et al., (2011) and recently by Navarro-Cerrillo et al. (2020) in TNP. Tree density and basal area were higher at high elevation and were associated with increasing rates of mortality. Besides, the percentage of young trees increased as elevation did, while mean tree diameter was relatively higher at lower elevation (Table 1). A similar decrease of tree height as elevation rises has been reported by Benabid (1982). The resulting diversity of tree sizes is lower at low altitude. Regarding species richness, overall diversity of woody species was relatively higher at middle and lower elevations (plots 4 and 6).

MF attained relatively higher mean dbh values in less dense stands (plots 2 and 7) which have lower values of Shannon’s dbh index (H’ dbh, Table 1). Similar relationships between size diversity, wood species diversity and stand density were found by Linares et al., (2011). Compared to Spanish fir forests, MF forests are characterized by higher species diversity (Aafi, 2000; Benabid, 2000; Linares et al., 2011) due to the presence of several tree species mainly C. atlantica, A. granatense, P. nigra and J. oxycedrus. The presence of these species with their structural features descripted in the present study was similar to the results found by Navarro-Cerrillo et al. (2020).

Nevertheless, our results also support other environmental factors affecting MF forest structure and diversity. Hence, stumps frequency, especially of small sized trees indicating human perturbation, mainly at low elevation (plots 8, 7, 6, 5 and 2). In the Rif, forests have deeply suffered from an abusive exploitation since the protectorate period (Benabid, 2002) which was reflected by the absence of very large size classes (Fig. 2). However, before this period, the exploitation by local populations was relatively less important although it has considerably increased over last decades (Benabid, 2002). In a comparison between Moroccan and Spanish fir stands, Linares et al., (2011) reported a widespread logging activity by local inhabitants in Talassemtane forest. Furthermore, anthropogenic pressure, especially the extension of Cannabis crops, remains a major threat to the maintenance of MF forests (Benabid, 2000; Taïqui, 2005). Previous studies reported the predominance of these crops at low altitudes; however, at present, Cannabis crops are increasingly becoming the main driver of deforestation of MF and other mountainous forests in Northern Morocco.

We found a significant spatial aggregation of Moroccan fir trees, at a fine scale of up to 3 m, as showed by the univariate pair-correlation function (Fig. 3a). This result is refined by the bivariate analysis of two size classes; MF young individuals are mainly concentrated around adult trees at smallscale (up to 2 m; Fig. 3b). Additionally, marked point pattern analyses showed that individual tree dbh is independently distributed from location, which indicate that aggregation does not affect tree size (Fig. 4a). As suggested for other Mediterranean ecosystems (De Luis et al., 2008), the observed MF pattern may be produced by a combination of various mechanisms including its limited seed dispersal, shade tolerance and the ecological effects of spatial heterogeneity of resources and conditions. The general context of climatic stress and human pressure should be also considered.

Spatial patterns of tree species usually show a scale-dependent aggregation that is generated by seed dispersal limitation (De Luis et al., 2008; Zhang et al., 2013). Wind dispersal can cause short distance dissemination away from source leading to seed shadow around adult trees (Wiegand et al., 2007). Seeds that did not dispersed far from adult trees are accumulated under the crowns so that few juveniles are found far from the adult stems (LeMay et al., 2009; Raventós et al., 2011; García-Cervigón et al., 2017). Like other Mediterranean firs and pines (García-Cervigón et al., 2017; Abellanas & Pérez Moreno, 2018), MF seeds are wind-dispersed. Melhaoui (1990) estimated that the dispersal distance of MF is almost always less than 5 m. Our finding supports this estimation since the aggregation of juveniles around MF adults occurred up to 2 m (Fig. 3b).

However, seedling establishment at short distance from adults depends on species preferences and it is usually thought to be affected by spatial heterogeneity and its effects on habitat quality (Zhang et al., 2013). For example, shade tolerance or intolerance abilities can play an important role in species spatial pattern (Abellanas & Pérez-Moreno, 2018). But unlike wet sites where seedlings development is strongly influenced by light (Chave, 2000), the limiting factor at dry sites seems to be moisture where large trees offer moister microsites to small trees (LeMay et al., 2009). In our case, like the majority of fir species, MF is a shade-tolerant species that regenerate under canopy layer (Baumer, 1977). Previous studies on MF also found that the moisture factor plays a crucial role in natural regeneration (Melhaoui, 1990). Clumping of MF juveniles under adult trees is due to the favorable moisture conditions below canopy and shade-tolerance ability of juveniles.

Taking separately, almost all plots showed aggregated pattern of MF trees. The exceptions were plots 6 and 8 for the univariate analysis and plots 5 and 8 for the bivariate one, which showed random and independent patterns, respectively (Fig. S2 and S3 [suppl.]). These departures from the general findings might be explained by a higher anthropogenic pressure reflected by the high amount of stumps in these plots. Indeed, many studies reported that anthropogenic activities such as logging can modify the spatial structure of the stand (Motta & Edouard, 2005).

In addition to seed dispersal limitation and environmental heterogeneity, biotic interactions play an important role in the dynamics of an ecosystem and can induce specific structures (Goreaud & Pélissier, 2003). Previous studies have used second order statistics to study the relationships between spatial patterns and biotic processes, and have reported positive, independent and negative correlations between tree positions and sizes. In a study of Pistacia atlantica in a semi-arid region of Iran (Erfanifard et al., 2018), bivariate analysis revealed a spatial aggregation of saplings around females while mark analysis using k_mm(r) showed a negative correlation between dbh and position. In this case, aggregation negatively affects size. In a large study area covering 30 ha in a temperate old-growth forest in China, Zhang et al., (2013) identif ied mostly independence of tree dbh at greater distances using k_mm(r), but highly significant negative correlation at short distances for the majority of species and some cases of significant short-distance attraction of dbh associated with the light-demanding shrub species. Jácome-Flores et al., (2016) found at Martinazo site in Spain a clustered distribution of dwarf palms Chamaerops humilis individuals while their size did not show a spatial correlation.They suggested that this result can be seen when using small sample sizes, which produce noticeably wider simulation envelopes. Getzin et al., (2008) showed no correlation of the studied marks (dbh, crown area and upper crown area) using k_mm(r) at all scales in a low-density deciduous stand in Germany which imply low overall competition effect. Similar results were found for Pinus syivestris in Poland where no correlation was found between trees either using diameter or height of pine individuals (Erfanifard & Stereńczak, 2017). Law et al., (2009) reported that the absence of correlation may suggest that the size is much more equally distributed over space than the individuals are themselves. In our study, mark correlation analysis of Abies marocana trees in Talassemtane forest showed a non-significant negative correlation of tree sizes up to 2 m, as well as at 9-11 m where a negligible negative correlation appeared (Fig. 4a). This non-significant negative correlation may be attributed to the occurrence of a weak competitive effect at small distance. Negative correlation indicates aggregation of young individuals and may be related to clumped dispersal and establishment patterns, which are absent at later developmental stages due to competition and mortality.

At the plot level, we found no correlation of tree sizes except for plot 2 which showed a very weak and non-significant negative correlation up to a very short distance (< 0.5 m) (Fig. S4 [suppl.]). The significant pattern detected by the average function at small scale based on replicated pattern could be reinforced by the additional effect of weak negative correlation tendency observed in plots 5, 7 and 8 (where the k_mm(r) expected value was closed to the lower simulation envelope). However, the facilitative effect recorded among young and aged trees may hide the competitive effect occurring within groups of regenerating individuals.

Furthermore, low density of suppressed trees as well as their spatial pattern confirms weak competitive interaction and self-thinning process. In fact, the result of random labeling either by g₁₂(r) or p₁₂(r) showed a weak but significant repulsion between suppressed and dominant individuals at 2.5-3.5 m (p = 0.005, Fig. 4b and 4c). The highest density of suppressed trees in association with a higher rate of dead trees in the densest plot 3 (Table 1) may indicates a density-dependent suppression process but which could be significantly detected by test functions only at a limited scale. This negative association appears to be a small-scale pattern only emerging within a 2.5–3.5 m range which reflects a fine patchy pattern where small and almost distinct groups of dominant and suppressed trees occur separately within the stands (Carrer et al., 2013). As most of the suppressed trees were small trees, and MF regeneration was abundant in most plots, the spatial negative correlation in tree size could contribute sporadically to suppression of tree growing close to each other in denser stands. Accordingly, MF juveniles clumped in the shelter of the adults and beneficiate from a positive recruitment interaction which seems to be more important than an intra-specific competitive effect.

Here, the random labeling considered for each plot showed the distinct behavior of the four last plots (Fig. S5 [suppl.]), and thus support the previous observations in the three other analysis types (Fig. S2, S3 and S4 [suppl.]). Indeed, a negative association between dominant and suppressed trees is found in plots 5, 8 and particularly 7 which extends to 4.5 m, in addition to a significant positive association at distance of > 14 m shown in this last plot. At the small scale, suppressed trees appear to be separated from dominant ones while they are more grouped at larger scale.

Under current climate conditions and human pressure in TNP (Taïqui, 2005; Linares et al., 2011; Ghallab & Taïqui, 2015), seedlings and saplings of Moroccan fir regenerate and maintain themselves in favorable moisture conditions under the canopy of adult trees thanks to their shade tolerance and the absence of significant intraspecific competitive effect. Spatial intraspecific negative interaction may sometimes be undetectable because MF can withstand competition for light since it is a shade tolerant species (Baumer, 1977; Abellanas & Pérez-Moreno, 2018). This suggests that intraspecific competition did not appear to be a major process that negatively affects MF growth in the studied stands, which is just the opposite of what was found for the Spanish fir counterpart (Linares et al., 2010). In contrast, the facilitative interaction where parents provide shelter from abiotic stress and protection against biotic consumption was more relevant.

Although our study has provided insights into the spatial patterns of Moroccan fir, there are some limitations regarding the relative number of points used within each plot, particularly in the case of analyses involving two pattern types. Indeed, Wiegand & Moloney (2014) recommended the use of 50-70 individuals as a minimum in spatial point analyses. However, Velázquez et al. (2016) reported that many studies used relatively few individuals in point pattern analysis (< 100). For instance, Camarero et al.. (2000) used at least 15 individuals, Zenner & Peck (2009) used 20 individuals. More recently, Janík et al. (2016) and Cordero et al.. (2016) used much less trees (10 and 4, respectively). In our case, we had a sufficient number of trees in univariate analyses, and despite that the available number of trees was split into two classes in the bivariate and random labeling analyses, this limitation was minimized by using replicate analysis with data from the eight sampled plots that were combined into one average function (Illian et al.., 2008; Wiegand & Moloney, 2014). This method is a straight-forward extension of the typical analysis of separate plots (De Luis et al.., 2008) and allows increasing the sample size and the statistical power (Illian et al.., 2008). Although being considered as an extension of the common analysis of isolated plots (De Luis et al.., 2008), replicate analysis was rarely used by ecologists in spatial point patterns (but see Raventós et al.., 2011).

The results of our explanatory study allowed us to describe the overall spatial pattern of the Moroccan endemic f ir. Analysis of the main MF structural characteristics in the sampled plots indicated a remarkable regeneration dynamic. However, the anthropic pressure remains a serious problem which constitutes a constraint to MF community preservation. Moroccan fir showed strong positive intraspecific association and juveniles were significantly attracted to adult trees. We highlight the importance of shade tolerance, seed dispersal and microsite heterogeneity in the observed facilitative intraspecific effect. However, more investigation is needed to test subsequent hypothesis, in particular by integrating other species (e.g. Cedrus atlantica). Understanding spatial structure and dynamics of trees, especially endemic species, is critical particularly in the context of global change.

 

Figure 2. Frequency distribution of diameter at breast height (dbh) of Abies marocana in all sampled plots.

 

Figure 3. (a) Univariate pair correlation function g(r) of the overall Abies marocana. (b) Bivariate pair correlation function g12(r) of two fir size classes: juveniles with dbh ≤ 10 cm (pattern 2) and adult trees with dbh > 10 cm (pattern 1). Dotted line indicates the empirical curve, gray line indicates the expected value under the null model (Heterogeneous Poisson for univariate g(r) and independence for bivariate g₁₂(r)) and black lines indicate limits of 199 Monte Carlo simulation envelopes.

 

Figure 4. (a) Univariate mark-correlation function k_mm(r) of Abies marocana tree size (dbh). (b) Random labeling analysis of firs (living trees vs suppressed trees) using partial-pair correlation functions g₁₂(r) and (c) mark connection function p₁₂(r). Dotted line indicates the empirical curve, gray line indicates the expected value under the null model (independent marking for k_mm(r) and random labeling for g₁₂(r) and p₁₂(r)) and black lines indicate limits of 199 Monte Carlo simulation envelopes.