Individual- based plant– pollinator networks are structured by phenotypic and microsite plant traits

1. The biotic and abiotic context of individual plants within animal- pollinated plant populations can influence pollinator foraging behaviour. Pollinator movements regulate pollen flow among plant individuals, and ultimately determine individual plant reproductive success. Yet the underlying drivers of this context dependency of interactions at the population level and their functional consequences for indi viduals remain poorly known. 2. Here we used a well- characterised population of Halimium halimifolium (Cistaceae), a Mediterranean shrub species, in combination with exponential random graph models (ERGMs) to evaluate how the intrapopulation variation in plant attributes configures individual- based plant– pollinator networks and determines their re productive outcomes.


| INTRODUC TI ON
The foraging movements of pollinators influence how pollen transfer occurs among animal-pollinated plant individuals within a population (Morris et al., 1995), ultimately influencing their reproductive outputs. Within a generalist pollination system, a plant population is usually composed of individuals differing in their level of pollinator attraction, and therefore, the pollinators visiting a given plant are a subset of the pool of potential pollinators (Gómez et al., 2007;Herrera, 2005). These inter-individual differences in the pollinator assemblage of individual plants partially contribute to generate nonrandom patterns of pollen transfer events (Valverde, 2017), which may affect population dynamics and the local structuring of genetic diversity through a network of assortative mating events.
Individual plants within most populations differ in their phenotype and this variation among plants (e.g. in plant and flower size or flowering phenology) may cause differences in their pollination niche due to distinctive preferences by different pollinator guilds (Gómez et al., 2014). Therefore, individual plants displaying similar phenotypes are expected to interact with a similar pollinator assemblage, promoting mating events among those plants (Gómez et al., 2011). Furthermore, dissimilarities in phenological patterns among pollinator species may create temporal variation in pollinator availability, leading to differences in the pollinator assemblage interacting with plants diverging in flowering phenology (Fox, 2003;Primack, 1985). On one hand, individual plants with higher flowering synchrony within their population are expected to share pollinators with more conspecifics; but, on the other hand, they may also interact with fewer pollinators due to intraspecific competition. Thus, a higher flowering synchrony among conspecifics can influence the level of pollinator sharing and therefore, increase intrapopulation assortativity in mating (Weis, 2005).
Variation in the pollinator niche among individual plants may also be influenced by their spatial position in the population (Rodríguez-Rodríguez et al., 2017). The spatial position defines the ecological context of plants, such as the microtopographical conditions and the local plant community to which they are exposed (e.g. intraspecific and interspecific competition), generating different aspects of context dependency of interaction outcomes (Chamberlain et al., 2014).
Because pollinator species may tend to forage in specific microsites, the fine-scale location of individual plants within a population is expected to affect their interactions with pollinators and consequently, the patterns of pollinator sharing among conspecific plants (Rodríguez-Rodríguez et al., 2015). Most plant-pollinator interaction studies do not address this spatial component, which might entail important constraints given the fine-scale division of resources among insect pollinators in a plant population (Dupont et al., 2014) and the spatial heterogeneity in neighbourhood composition (Janovský et al., 2013).
Plant-pollinator communities have been widely studied during the last decades using bipartite interaction networks (Bascompte & Jordano, 2007), where pollinator species are connected to the plant species they visit. However, species-based networks overlook the existing intraspecific variation by pooling all the individuals in a single group (Dupont et al., 2014;Tur et al., 2015). In contrast to the community-level approach, individual-based networks can be used to unveil processes at the population level (e.g. Rodríguez-Rodríguez et al., 2017;Soares et al., 2021;Valverde et al., 2016) and convey a direct link to analyse demographic and evolutionary questions.
In a bipartite, individual-based network within a plant population, pollinator species are connected to the individual plants they visit.
When it is projected into a unipartite network, individual plants are visualised as nodes and the links connecting those nodes depict the level of pollinator sharing, which may serve as a proxy of mating probabilities (Fortuna et al., 2008). Thus, the unipartite network generated with this procedure is undirected, in contrast with real pollen-transfer networks. This undirected unipartite network allows us to elucidate how interactions and hence, potential mating events, are structured within a plant population.
Because intrapopulation variation in phenotype, phenology and microsite may determine how individual plants interact with pollinators, it is expected to also affect the relative position of those plants within the unipartite network. While there are clear hypotheses suggesting that plants occupying central positions may benefit from a higher level of pollinator sharing with conspecifics, increasing the probability of mating events and hence their reproductive success (Gómez & Perfectti, 2012;Soares et al., 2021), it is not trivial how to translate network topologies into functional implications. Recently developed exponential random graph models (ERGMs) allow us to evaluate how plant attributes influence the configuration of complex networks, advancing their study from descriptive metrics into a more cohesive predictive framework (Harris, 2014). Such analyses are useful to test hypotheses about processes that shape networks by modelling the probability of interaction establishment as a function of the characteristics of individuals (nodes) and the value of other links within the network (Snijders et al., 2010). ERGMs differ from traditionally used descriptive approaches in that it encourages and supports identifying underlying mechanisms that explain the whole structure of the network. Although this analytical tool has been widely used in social network studies in recent years, it predictive research on the evolutionary and ecological processes that give rise to complex ecological networks at the population level.

K E Y W O R D S
Doñana, ecological networks, Halimium halimifolium, intraspecific variation, mutualism, plant fitness, plant-pollinator interactions has rarely been considered for ecological network analysis (but see Miguel et al., 2018).
Here we aim to investigate how the intrapopulation variation in intrinsic (i.e. phenotype and phenology) and extrinsic (i.e. microsite) plant attributes shapes the structure of individual plant-pollinator networks and influences their functional outcomes (Figure 1). To address this question, we recorded individual plant-pollinator interactions in a population of the Mediterranean shrub Halimium halimifolium (Cistaceae). First, we assessed how the intrinsic and extrinsic plant attributes influenced the emerging configuration of the bipartite plant-pollinator network and the unipartite plant-plant network derived from pollinator sharing. Second, to evaluate how this association between plant attributes and network structure translates into functional outcomes, we analysed how individual female plant fitness (i.e. the total seed weight per plant) was affected by the combination of plant attributes and the topological position of individual plants within the pollination network.

| Study site
The study was performed in Doñana National Park (37°07′52.2″N 6°31′40.9″W, 12 m a.s.l.) within a 1.2 ha area bounded on its four sides by a stream, pine forests and semi-natural grasslands. Doñana National Park is located on the Atlantic coast of southwestern Spain, in an area with a Mediterranean-type climate where the vegetation is composed mainly of sclerophyllous shrublands. Halimium halimifolium (Cistaceae) is a Mediterranean shrub that occupies the slopes of stabilised sand dunes (Díaz-Barradas & García-Novo, 1987) and usually blooms in late spring (May-June; Herrera, 1988) when all other co-occurring plant species have finished their flowering season. The study site is largely dominated by this plant species, which represents nearly 60% of the total vegetation cover. Although other plant species are also present in our study area ( Figure S2), H. halimifolium was the only species in flower during the study period. Hence, we did not consider competition or facilitation for pollinators among the different plant species. Each H. halimifolium plant can produce up to 1,000 hermaphrodite yellow flowers ( Figure S1) at once and flower opening occurs synchronously each day in the population.
Flowers usually exhibit a conspicuous dark brown spot at the base of each petal, being the presence and size of these spots highly variable among individuals within populations ( Figure S1). Although these dark spots have not yet been proven to act as floral guides in H. halimifolium, this kind of petal marks is widely recognised to enhance pollinator attraction in multiple plant species from a wide range of families (de Jager et al., 2017). As this plant species is selfincompatible and cannot reproduce asexually, cross-pollination promoted by insect visitation is necessary to complete seed production and allow population persistence (Talavera et al., 1997).

| Sampling
We conducted surveys to record pollinator visitations in the study site during the peak flowering period of H. halimifolium (25 days be- with pollinator species were recorded in 60 of those individuals, which were selected using a stratified random sampling where both isolated and aggregated plants were well represented ( Figure S3).  (Talavera et al., 1997). We assessed the completeness of sampling effort by calculating Chao asymptotic estimators (see Appendix S1 and Figure S4). Video cameras were placed approxi-  (Table S3) Two plants that are strongly connected in the unipartite network, and thus share most pollinators, will be more likely to be involved in mutual pollen transfer, and a potential mating event, than two plants with very different pollinators, connected by a weaker link. We thus considered the link weights estimated in this way as a proxy for the probability of mating events mediated by a given pollinator species (see Gómez et al., 2011;Rodríguez-Rodríguez et al., 2015).

| Plant attributes and fitness estimation
We obtained a set of intrinsic and extrinsic attributes for each sur-  Table S1. Extrinsic variables and plant size were estimated using aerial images, an orthomosaic and a 3D surface model obtained with the performance of drone flights in the study population (see Appendix S1 for detailed methods of spatial data acquisition).
To quantify the female reproductive success per individual plant, we first calculated the fruit set as the proportion of flowers setting fruit. To do this, we collected 10 inflorescences per plant and counted the total number of flower buds produced initially and the number of fruits at the end of the season. Second, we sampled five fully developed fruits and counted the number of seeds per fruit to obtain the average seed production per fruit in each plant, and also measured the average seed mass with a precision scale. An overall female reproductive success ('fitness' hereafter) per plant was estimated by weighing the total number of seeds per plant by seed mass (Appendix S1).

| Statistical analyses
2.4.1 | Ecological correlates of the bipartite plantpollinator network structure To analyse the ecological variables producing the overall structure of the individual-based pollination networks we built ERGMs (Lusher et al., 2013). An ERGM is a true generative statistical model of network structure and characteristics (i.e. it allows us to describe how data are generated in terms of probability): Node traits and local structural properties can be used to predict properties of the entire network (e.g. diameter, degree distribution, etc.). The probability of a given network, which can be re-expressed in terms of the conditional log-odds of a single link between two nodes, is modelled as a function of terms that represent network, node or link features.  ing variation in fitness was assessed using the relaimpo r package (Grömping, 2006). This package uses the fitted linear model explained above with fitness as the response variable and intrinsic and extrinsic plant attributes as predictor variables and estimates the R 2 contribution of those variables to explaining variation in fitness averaged over resampled orderings among regressors.

| Functional consequences
All analyses were performed using R software version 3.5.3 (R Core Team, 2019).

| Ecological correlates of the bipartite plantpollinator network structure
We found a 'sum' effect in the ERGM fitted with the bipartite network as the response variable (Table 1) The model coefficients for exogenous variables (Table 1) Table 1  Overall, intrinsic (i.e. phenotype and phenology) and extrinsic (i.e. microsite) plant attributes accounted for 9.54% and 8.01%, respectively, of the variation in plant-pollinator interaction odds, while the pollinator functional group explained 74.42% of this variation.

| Ecological correlates of the unipartite plantplant network structure
The 'sum' effect was significant in the ERGM with the unipartite network as the response variable (Table 2)

| Functional consequences
Plant fitness and intrinsic plant attributes were not spatially autocorrelated (Table S5). The GLM with individual plant fitness as the response variable revealed that the topological position of individual plants within the unipartite network (i.e. closeness centrality) was positively associated with plant fitness (Table 3; Figure 3a). We also found that individual plant fitness was associated with intrinsic and extrinsic plant attributes, after controlling for the network-mediated effects (Table 3)

| D ISCUSS I ON
Our findings revealed how the individual variation in plant attributes existing in natural populations gives rise to the overall structural patterns of individual-based plant-pollinator networks and its functional implications. We found that intrinsic (i.e. phenotypic traits and flowering phenology) and extrinsic (i.e. microsite) plant attributes had comparable influences on structuring both, the plant-pollinator bipartite network and the plant-plant unipartite network.
Furthermore, our results also showed that network structure, and specifically the position of individual plants within this network, together with the direct effects of plant attributes, had considerable consequences for individual plant fitness. Overall, these findings suggest that the ability of individual plants to participate in shared interactions, and therefore to increase their reproductive success, may ultimately be determined by the intraspecific variation in microsite characteristics, phenotype and phenological traits.

| Ecological correlates of the bipartite plantpollinator network structure
Variation in traits of both plant and animal partners is increasingly recognised as forces structuring ecological networks at the level of interactions among species (Bartomeus et al., 2016). Similarly, intraspecific trait variation can play an important role in shaping individual-based networks, yet this aspect remains virtually unexplored in the literature (Miguel et al., 2018). By using ERGMs we were able to further recognise specific plant attributes (individual phenotype, phenology and microsite) that shape the observed individual plant-pollinator network structure. We thus assessed the relative contribution of extrinsic and intrinsic plant attributes to explaining the odds of interactions to increase for individual plants.
Our findings showed that intraspecific trait variation in floral and vegetative phenotypes among plants influenced the establishment of interactions with pollinator species, which is in accordance with previous work (Dupont et al., 2014;Gómez & Perfectti, 2012;Weber & Kolb, 2013). We also demonstrated that this influence upscales to affect the structure of the individual-based pollination network. Specifically, we found that the probability of interacting with more pollinator species increased with the size of the flower guides, which is widely known to positively affect pollinator attraction (Waser, 1983). Although plant size is expected to increase plant attractiveness for pollinators, we found a negative effect of plant height on the probability to interact with pollinator species. This  increase their chances of receiving better quality pollen in comparison to asynchronous plants (Antonovics & Levin, 1980).
Although the effect of the microsite, especially neighbourhood composition and density, has previously been considered for plantfrugivore (e.g. Carlo, 2005;Guerra et al., 2017;Miguel et al., 2018;Morales et al., 2012) and plant-pollinator interactions (e.g. Delesalle & Mazer, 2002;Lázaro et al., 2009), it has rarely been included in network studies. Previous studies support the notion that the local microclimatic variation affects pollinator composition, activity and behaviour at flowers (Herrera, 1995;Potts et al., 2005). In light of all of these results, the context dependency previously found in species-level networks (Poisot et al., 2015)  attributes, similar individual plants were likely to share more pollinator species, which may translate into potential assortative mating in the study plant population (Waser, 1983). Besides, the similarity in traits consistently associated with greater pollinator sharing showed larger effect sizes than the similarity in traits more associated with divergence in interacting pollinator species. These results also suggest that pollinator preferences for plant attributes are distinct among pollinator species even when foraging within a monospecific flowering stand, scaling up to affect the entire network configuration. It must be pointed out that the between-plant spatial distance was the variable that most influenced the structure of the unipartite network of pollinator sharing, highlighting the importance of accounting for the spatial component when studying plant-animal interactions (see also Dupont et al., 2014;Pasquaretta et al., 2017).
Our interpretations of the results outlined above provide evidence that distinct plant attributes operate in different ways when structuring bipartite and unipartite individual-based pollination networks.

| Functional consequences
The actions. By contrast, previous work found that the plant topological role had a higher effect than plant attributes on plant fitness (Gómez & Perfectti, 2012;Rodríguez-Rodríguez et al., 2017;Soares et al., 2021), but this could be strongly dependent on the spatial context and the species-specific characteristics.
The use of centrality measures is not as powerful as the ERGM approach to describe network topology; however, the ERGM design in its current form does not allow us to assess the fitness consequences of the overall network structure. While ERGMs offer real potential to test hypotheses about network configuration, developing an integrative modelling framework to improve our understanding of the linkage between network structure and functionality is still a current challenge in ecology.
Altogether, our findings indicate that the effect of both intrinsic and extrinsic plant attributes on structuring the pollination networks can ultimately determine differences in plant fitness, even at the scale of small populations. The emergence of the pattern of nonrandom, assortative mating events found may create a spatial genetic structure within the population as a result of the spatially restricted gene flow generated (Epperson, 1993). Our approach, bridging the complexity of the whole pollination network and its direct effects on among-individual fitness variation, would contribute to explaining the higher levels of genetic structure found in most animalpollinated plant populations (Vekemans & Hardy, 2004), where the mobility and preferences of different pollinator groups may drive a high frequency of assortative mating events (Valverde, 2017).

| CON CLUS IONS
Our study opens up new analytical approaches to further understand how complex ecological networks are configured at the population level by disentangling the underlying drivers of the context dependency of interactions, that is, plant traits versus microhabitat attributes. Individual-based networks and ERGMs provide novel insights into how individual-plant variation in both intrinsic and extrinsic attributes predictably shifts the functional outcomes of mutualistic interactions. This close relationship between the inter-individual variation in plant attributes and their fitness might be especially relevant in habitats with patchiness in biotic and abiotic variables, which affect animal foraging behaviour within a plant population. Such a predictive framework has potential to move forward ecological network analyses from descriptive metrics to more powerful forecasting approaches. Therefore, the analytical tools we used would help to better assess how the pervasive context dependency of species interactions influences the evolutionary and ecological processes taking place within populations, and to predict the potential responses of these populations to environmental changes.

ACK N OWLED G EM ENTS
We thank Ricardo Díaz-Delgado and ICTS-RBD for technical support with drone flights and spatial data acquisition, as well as logistic support in the study area. We are grateful to Ana Benítez, Irene Mendoza, Lisieux Fuzessy, Elena Quintero, Jorge Isla and Daniel Pareja for helpful comments during fieldwork and data analysis. We also thank Pietro Maruyama and an anonymous reviewer for insightful comments and suggestions which improved this manuscript.
Doñana National Park authorities approved permissions to work in the study area. Financial support was provided by CSIC (JAE Intro fellowship JAEINT18_EX_0080 to BAC) and the Spanish Ministry of Science and Innovation (FPU fellowship to B.A.-C., projects CGL2017-082847P to P.J. and CGL2017-92436EXP to I.B.).

PEER R E V I E W
The peer review history for this article is available at https://publo ns.

DATA AVA I L A B I L I T Y S TAT E M E N T
Data and r code for the analyses are available at the Zenodo Digital