Many plant species present a non-random spatial distribution, forming aggregates in response to abiotic (e.g., soil moisture, sunlight incidence) and biotic (e.g., competition, seed dispersal) factors (^{Nathan & Muller- Landau 2000}). This is the case of Puya coerulea Lindl. (Bromeliaceae), a plant species endemic to Chile, distributed from 30º to 38ºS, corresponding to one of the seven Chilean species within the Puya genus, and comprised within the subgenus Puyopsis (^{Smith & Downs 1974}). This bromeliad is subdivided into four varieties (Smith 1970) within which we found P. coerulea var. coerulea, present at our study area. This variety is distributed from the Valparaiso Region (cerro La Campana, 32º 57’S) to the O’Higgins Region (Cachapoal, 33° 55’S), ranging between 500 and 2000 m of altitude, and mainly growing in rocky outcrops at semi-arid areas. Examining plant species spatial patterns can provide valuable insights on the underlying ecological processes shaping those patterns (^{Stoll & Prati 2001}). This is particularly interesting for those species that are usually found forming densely aggregated and discrete distributions. Therefore, the aim of this study was to assess the abiotic factors determining P. coerulea spatial distribution.

This study was carried out during October 2017, at the Río Clarillo National Reserve (33º 45’S 70º25’ W), located 45 km from the city of Santiago. We collected data along a 2-km transect (at both sides of the road) with suitable conditions for P. coerulea. At this transect, we defined 71 random points separated by ~20 m from each other. At each sampling point we determined presence or absence of P. coerulea and measured the following abiotic variables: UV radiation (using a Vetus UV340B portable UV-meter, placed horizontally at 1.30 m above the ground), percentage of rock cover (visually estimated, as the fraction of a 5x5 m plot covered by rock), and soil inclination and slope orientation (both measured using a digital inclinometer). At each sampling point we took two measurements, one at each side of the road. Additionally, we registered the accompanying plant species present at each sampling point. Then, we obtained geographical coordinates from each P. coerulea cluster (composed by one to five individuals, using a Garmin Map 62s GPS device with error ≤ 5 m) found at the study area.

We performed a Mixed-Effects Generalized Linear Model (GLMM) with a binomial error distribution to determine which abiotic variables (included as fixed predictors) explain P. coerulea presence or absence; road side was included as a random factor (^{Zuur et al. 2009}). GLMM goodness-of-fit was assessed using a Hosmer- Lemeshow test. Then, we performed an Analysis of Similarities (ANOSIM; using a Bray Curtis similarity metric and 1,000 permutations) to compare plant species composition between sampling points with and without

P. coerulea. ANOSIM is a nonparametric test aimed to test community composition differences, it is based on similarity matrices and its significance level is calculated using permutations (^{Clarke 1993}). GLMM and ANOSIM analyses were conducted in R 3.5.0 (R Core Development Team 2018) using the packages ‘vegan’, ‘pscl’, ‘mgcv’, ‘lme4’, ‘lmerTest’ and ‘pbkrtest’. Finally, we performed point-pattern analysis using the geographic coordinates of P. coerulea clusters obtained in the field. We used the Programita software (freely available at www.programita.org) to perform a spatial aggregation analysis using the g(r) function in an irregular area defined by the entire area of rocky outcrops present at the study area (i.e., P. coerulea suitable habitat). To evaluate the significance of the observed pattern, we constructed a null model upon 999 simulations under a model of complete spatial randomness (^{Wiegand & Moloney 2004}, 2014). To improve the visualization of the spatial aggregation pattern obtained, the values obtained were divided by the expected value and transformed to log (x + 1). Original data and R scripts associated to this submission is available from the figshare repository: https://doi.org/10.6084/m9.figshare.7045094

The presence of P. coerulea at the study area was largely determined by the percentage of rock cover (Table 1), indicating that this factor is the main determinant of its presence over the other environmental factors measured. The GLMM model fitted well to our data, according to the Hosmer-Lemeshow test performed (χ² = 6.653, df = 8, P = 0.575). Also, the composition of the accompanying plant species in sampling points with and without P. coerulea presence was significantly different (ANOSIM R= 0.345; P = 0.001), as there are only a few species that occur altogether with P. coerulea (Cistanthe grandiflora (Lindl.) Schltdl. (Montiaceae), Gamochaeta stachydifolia (Lam.) Cabrera (Asteraceae), Hypochaeris thrincioides (J. Rémy) Reiche (Asteraceae), Leucheria glandulosa D. Don (Asteraceae)), which are able to coexist with P. coerulea as they have similar ecological strategies. We found two large clusters of P. coerulea along the suitable habitat area (Fig. 1), showing a non-random spatial distribution. These results are corroborated by the point-pattern analysis conducted, which show a greater aggregation than expected by chance between 5 and 20 m (Fig. 2).

FIGURE 1 Distribution map of the georeferenced Puya coerulea clusters (yellow dots) at the Río Clarillo National Reserve. Suitable habitat (i.e., rocky outcrops) area is denoted with a red line. / Mapa de distribución de grupos de Puya coerulea georreferenciados (puntos amarillos) en la Reserva Nacional Río Clarillo. El área de hábitat apropiado (i.e., afloramientos rocosos) está delimitada por la línea roja. 

Steep slopes and rocky outcrops are rarely colonized by woody species (i.e., trees and shrubs), leaving open colonization opportunities for herbaceous plant species (^{Callaway & Walker 1997}). This could be the case of P. coerulea, being favored by adaptations to colonize and thrive in habitats that many other plant species cannot use. This is confirmed by the scarce accompanying flora co-occurring with P. coerulea at the rocky outcrops. Also, we might expect P. coerulea to present some photosynthetic adaptations for overcome water stress conditions (^{Jabaily & Sytsma 2010}). In this sense, the entire subgenus Puyopsis (to which P. coerulea belongs) presents only C3-type photosynthesis, whereas many bromeliads present a C3-CAM flexible strategy to deal with arid conditions (^{Rundel & Dillon 1998}). Comparing P. coerulea to one of the most studied species of the Puya genus, P. raimondii Harms (Jabaily & Sytsma 2010), it is possible that P. coerulea inhabits rocky outcrops because they provide microhabitat conditions with soils rich in organic matter, adequate humidity, aeration and radiation that provide protection to the strong temperature oscillations during the germination and sapling emergence stages (^{Vadillo & Suni 2006}). Another non-mutually exclusive explanation is the adaptation to these particular conditions to reduce inter-specific competition with other species that are unable to tolerate water shortage stress as Puya species do.

FIGURE 2 Spatial aggregation pattern according to the g(r) function for the mapped Puya coerulea clusters. / Función g(r) de agregación espacial calculada para las agrupaciones de Puya coerulea mapeadas. 

Puya coerulea is an important species of the Andean ecosystem, confined to the Coastal Range (^{Zizka et al. 2013}). Therefore, we may expect that the distribution of this species responds not only to rock cover but also to geological features of these rocks. Further, intraspecific spatial aggregation in plant individuals can be considered as a process dependent on the dispersion of propagules, vegetative reproduction (clonal growth) and environmental patchiness on those abiotic factors that influencing plant reproduction and establishment (^{Stoll & Prati 2001}, ^{Tirado & Pugnaire 2003}, ^{Syphard & Franklin 2004}). Puya coerulea have very small elliptical seeds, which could be useful to be mobilized via anemocoria (^{Varadarajan 1990}), being a fairly conserved trait within the Puya genus. This type of dispersal allows seeds to travel long distances, increasing colonization capabilities. However, limiting to the framework of observations made for P. coerulea it is possible that the effectiveness of this dispersal mechanism is affected by the substrate in which this species grows, certainly the seeds could be dispersed to distant sites, but few of them will reach suitable sites to germinate (^{Chesson 1991}). Based on this, already established individuals may be skewing seed rain to suitable areas with high survival probabilities, potentially creating a positive feedback on the patchy spatial distribution pattern, increasing spatial aggregation over time.