Beyond habitat: effects of conspecific and heterospecific aggregation on the spatial structure of a wetland nesting bird community

Nest location is a key factor influencing reproductive success in birds, and habitat choice is considered the main way in which birds select nest sites. Less attention has been devoted to the demand for proximity to other bird nests, which can provide additional profit, namely defense against predators. Here we analyzed the contributions of habitat, and conspecific and heterospecific aggregation to the spatial arrangement of breeding birds in a model bird community. We surveyed a pristine Siberian wetland bird community with the aim to locate all bird territories or nests, in 1993 and 2013. Habitat explained much of the nest site choice, but the nests were aggregated both intra- and inter-specifically more than the spatial pattern of the habitat could explain. In particular, ducks, grebes and some waders bred nearby the most abundant active nest defenders, such as gulls and terns. Heterospecific associations were particularly pronounced in 2013, when the community was impoverished and one common active defender (white-winged black tern Chlidonias leucopterus ) was replaced by a less numerous but aggressive predator (Mongolian gull Larus mongolicus ). The results suggest that spatial pattern in bird nests may be influenced by the (dis)appearance of one or a few species, which can play a role as umbrella or predator species. Integration of factors supporting the breeding of umbrella species, such as gulls, may became key targets for comprehensive conservation measures in large wetlands.


Introduction
Bird parents actively choose suitable nest sites to avoid nest loss due to predation and harsh weather (Lloyd et al. 2000, Caro 2005). Since a direct habitat choice has been considered the main way in which birds select nest sites, most studies dealing with species occurrence and numbers have emphasized the role of available habitats (Sutherland and Hill 1995, Davis 2005, Baschuk et al. 2012, Żmihorski et al. 2016. However, individuals of many species may select their nest sites also on the basis of the presence or absence of other animals (Caro 2005, Sebastián-González et al. 2010a. For example, nesting aggregation of conspecifics may inform individuals on potential breeding partners (Dubois et al. 1998), suitable food resources (Bayer 1982) and probability of breeding success (Danchin et al. 1998). In addition, the risk of predation can be reduced in such aggregations through defensive or dilution mechanisms Moldsvor 1992, Larsen andGrundetjern 1997). The species passively protecting their nest against predators and relying more on nest crypsis can benefit from shared nesting near active defenders (Cramp and Simmons 1983) who can increase their nest success (Larsen and Grundetjern 1997, Šálek and Šmilauer 2002, Sládeček et al. 2014, Rocha et al. 2016). All these reasons may increase the attractiveness of nesting in both conspecific and heterospecific aggregations.
The particular contribution of the nest site selection based on inter-individual interactions becomes less clear in large multi-species nesting aggregations, as different tactics of various species can confound the multiple effects prevailing in such systems. Accordingly, most studies refer to or have been focused on describing the spatial arrangement of bird nests in simple communities of two or a small number of species (Eriksson and Götmark 1982, Burger 1984, Larsen and Moldsvor 1992, MacDonald and Bolton 2008, Bentzen et al. 2009, Cunningham et al. 2016, Kubelka et al. 2019). Interindividual relationships have been used as an explanation of spatial distributions of birds at coarser landscape scales (Boone and Krohn 2000, Drapeau et al. 2000, Hobson et al. 2000, Titeux et al. 2004, McIntire and Fajardo 2009, and in less species-rich communites inhabiting boreal forests (Heikkinen et al. 2004) and farmland (Freemark and Kirk 2001). On a local scale, comprehensive studies separating the effects of conspecific and heterospecific aggregation from the sole habitat effects are lacking.
Wetland bird communities are highly diverse in terms of species richness and habitat requirements (Paracuellos 2006(Paracuellos , Żmihorski et al. 2016. Many species of these communities breed in groups and use different anti-predator tactics (Gochfeld 1984, Sládeček et al. 2014, Cunningham et al. 2016. Distinguishing the multiple drivers of nest site choice, possibly extending far beyond the exclusive selection of preferred habitat, can be particularly important from a conservation perspective. Decline of some wetland species may be due to the presence or absence of conspecifics or other species rather than due to the decline or the alteration of their preferred habitats. Examples are the loss of umbrella species providing antipredatory defense (Pöysä et al. 2019) or, on the contrary, the emergence of new settlers that act as potential nest predators (Thomas 1972). Distinguishing these drivers of nest site choice can strengthen the conservation implications of primarily habitat-focused observational studies (McIntire andFajardo 2009, Dray et al. 2012). For example, we may ask whether the species which breed in an impoverished community will disperse more with an emphasized selection of specific habitats, or will aggregate under the anti-predator umbrella regardless of specific habitat. In this case, the umbrella species may even play the role of keystone species in the community (Boogert et al. 2006).
To examine the importance of conspecific and heterospecific aggregation in addition to the effect of habitat choice, we analysed the spatial arrangement of nests within a diverse bird community inhabiting a large natural wetland in two breeding seasons separated by a period of twenty years. First, we tested unique and shared effects of habitat and other predictors of spatial aggregation on nest spatial arrangement of active nest defenders. We predicted shared nesting beyond the mere habitat choice in these species, because these breeders should benefit from joint defense against predators attracted to their aggregations. Second, we tested unique and shared effects of habitat, presence of active nest defenders and other spatial effects on nest pattern of cryptic breeders, which usually lack the ability to actively defend their nests. Because the active nest defenders can provide a protective umbrella for cryptic breeders, we predicted that cryptic breeders seek the proximity of active defenders in addition to the effects of habitat and other spatial effects. We discuss the findings in the context of different structure of the bird communities between the two years of study.

Study site
This research was conducted on the isthmus of the Svjatoj Nos Peninsula, Lake Baikal, Russia ( Fig. 1, Supporting information), between the coast of Lake Baikal and the vast marshland east of the lake. The marshland, with an area of approximately 300 km 2 , is covered mostly by a mosaic of open wetland habitats, and is one of the key areas for wetland birds breeding in a wide region around Lake Baikal (Mlíkovský et al. 2002, Mlíkovský 2009 Sládeček et al. 2014). Pine taiga provides refuge for generalist predators of bird nests, such as carrion crows Corvus corone, ravens Corvus corax and red foxes Vulpes vulpes (Sládeček et al. 2014).

Data collection
All fieldwork took place in the same 1.64 km 2 study area between 4 June and 13 July in both 1993 and 2013. The area sub-divided into 100 × 100 m squares was repeatedly searched to mark the location of birds displaying territorial behavior and map their nests. In 1993, the fieldwork was carried out by one person, who recorded all observations onto field maps and marked the nests with a colored ribbon tied to vegetation 5 m away. In 2013, the fieldwork was carried out by 3-6 observers moving at a distance of 10 m from each other and the positions were stored using GPS. The surveys included slow walking in shallower sections, while inflatable boat was used to reach vegetation patches on deep water with floating islets.
Four main habitat types were distinguished (  [water]. The proportions of these four habitats were recorded within the 100 × 100 m squares in a similar way in both years (Fig. 2). In addition, the water depth in the centres of the squares (in the field using a long rod) and the nearest distance to the marshland edge (using a map) were measured.

Bird species categories
We defined two bird species groups with distinct anti-predator defense tactics. The first group, 'active nest defenders' (listed in Supporting information), encompasses species that Figure 1. Location of the study area on the isthmus of the Svjatoj Nos Peninsula, Lake Bajkal, Russia. (A) Lake Bajkal, (B) Isthmus of the Svjatoj Nos Peninsula, (C) Position towards the lake coast, forest (pine taiga) and marsh (modified from Sládeček et al. 2014). are willing to actively attack approaching predators and drive them away from nest territories (Caro 2005). This active behavior can also provide protection for nearby nesting species that do not display such aggressive behavior. The second group, referred to here as 'cryptic breeders' (Supporting information), includes all remaining species that do not use aggressiveness as a major anti-predator tactic. Instead, they apply many behavioral alternatives, from discreetly leaving the nest or remaining motionless and cryptic on the nest to a distraction display before the predator approaches the nest (Gochfeld 1984, Caro 2005, Humphreys and Ruxton 2020. The species were sorted into the groups on the basis of published knowledge (Cramp 1977, Larsen and Moldsvor 1992 that was verified by our own field observations. Some duck and grebe nests found without an incubating adult or depredated in early incubation stages were difficult to determine reliably to species level. We therefore merged all fowl ducks into the genus Anas, both diving ducks (pochards Aythya fuligula and A. ferina) into the genus Aythya and the grebes into the genus Podiceps. These three genera of cryptic breeders were included in the analysis together with the remaining identified species. We singled out two species of gull -common gull Larus canus and Mongolian gull Larus mongolicus -as a specific category of nest defenders as well as potential nest predators due to their ambiguous nature (Kubetzki andGarthe 2003, Skórka et al. 2014). Passerines were excluded them from the analysis as they are not directly associated with wetland habitats.
The incubation stage and derived laying date were measured by floating the eggs (van Paassen et al. 1984) to test our prediction that if active defenders are to influence the nest site choice of cryptic breeders, they generally have to start nesting earlier (Larsen and Moldsvor 1992). We compared the timing of breeding between active defenders and cryptic breeders using the Mann-Whitney test of medians.

Response data and predictor groups
Individual 100 × 100 m squares were used as units in the analysis. Abundances of particular species or species merged into genera Anas, Aythya or Podiceps breeding in particular squares were included as response variables. As only part of the nests could be found within the abundant bird community in 1993, both confirmed and probable breeding events were included. Probable breeding events were assigned to the birds with territorial behavior observed more than once in the same square with an interval of at least five days. In contrast, it was easy to locate the nests in the impoverished community in 2013, so only confirmed breeding attempts were considered. The spatial nesting pattern of active nest defenders and cryptic breeders in 1993 and 2013 were analysed using four separate models. In all models, two groups of predictors in each square were distinguished: 1) habitats, including the percentage cover of the four habitat types, water depth and distance to marshland edge (sandy shore of the isthmus) and 2) spatial descriptors to explore the spatial effects within the marshland area, which included square centroid coordinates (latitude/longitude) and the spatial eigenvectors, calculated using the distance-based Moran eigenvector maps (dbMEM, Legendre and Legendre 2012) method. In addition, a third predictor group was added in the models for cryptic breeders, namely 3) presence and abundance of active nest defenders.
The distance-based Moran eigenvector maps (dbMEM) method processes adjusted matrix of the distances among sampling sites using principal coordinate analysis (PCoA), and in this way it creates a set of predictors enabling to describe the variation at any spatial scale. This full set of predictors must be reduced to a parsimonious subset using stepwise selection (see Legendre and Legendre 2012 for additional details). Although the spatial coordinates, which may represent linear 'global' trends across the study area, can also be reconstructed by combining the spatial eigenvectors of the dbMEM method, we followed the suggestion of Legendre and Legendre (2012, p. 868) to model these global trends separately. We therefore used the 'global' trends (if ascertained as significant predictors) as covariates when calculating the spatial eigenvectors. Selected spatial coordinates and spatial eigenvectors formed together the group of spatial predictors in our analyses. Significant spatial descriptors define the places where certain species are more or less abundant (concentrated) and therefore refer to a conspecific or heterospecific aggregation.

Community spatial analysis
To analyse the spatial pattern of breeding birds (spatial autocorrelation among the squares) in 1993 and 2013, we calculated Moran's I using the abundances of all nests across taxa in the squares (ape package, R ver. 3.6.2, <www.r-project.  Table 1) in the squares under study within a terrain map projection and their proportions (means ± standard errors) in 1993 (first column in graph) and 2013 (second column in graph) org>). For the habitats, we used a correspondence analysis (CA) of data on the proportions of the four habitat types in each square, separately for each year, to summarize numerically the most important gradient in habitat compositionthe first CA axis. We used the scores for individual squares on this axis (representing 67% and 73% of the total habitat variation in years 1993 and 2013, respectively) to estimate the spatial autocorrelation of the habitats, using Moran's I.

Multivariate analysis of predictor effects
We evaluated the compositional variation of the bird assemblages recorded within individual squares using redundancy analysis (RDA, Legendre and Legendre 2012), mostly in the form of partial RDA, including a priori covariates, to separate the effects of individual predictor groups (Supporting information). The original nest counts were log(x + 1)-transformed to achieve homogeneity of the unexplained (residual) variation, and very rare species, occurring in just one or two squares, were excluded from the analyses. All analyses also included the squares without any nest. To separate the effects of various predictor groups, we used the variation partitioning approach (Legendre and Legendre 2012). As we were testing for a hypothesized effect of active defenders on the occurrence of cryptic breeders, we performed the variation partitioning separately for each of these two groups of species, distinguishing just two groups of predictors for the community of active defenders, but adding the third group of predictors (representing the effect of active defenders) to the variation partitioning for the sub-community of cryptic breeders. To disentangle whether the cryptic breeders chose their nest sites on the basis of heterospecific attraction, we focused on the test of the b fraction in the variation partitioning for the cryptic breeders. This b fraction represents the variation explained uniquely by the active defender predictor, excluding the shared effect with habitat (d and g fractions) as well as other unmeasured variables (potentially reflected by the spatial predictors in the e and g fractions; Supporting information).
Care should be taken when interpreting results of variation partitioning when the compared groups have different counts of predictors, as the additive (and therefore partitionable) explained variation is expected to increase with the increasing number of predictors. Therefore, we interpret our results by using an adjusted form of the explained variation (estimated in the same way as R 2 adjusted in a linear regression), and by using mean square statistics, where the explained variation is divided by the number of degrees of freedom of the model. Note that these mean square statistics can be computed unequivocally only for the unique contributions of individual groups, not for the fractions representing overlaps in the explained variation. Percentage of explained variation and mean square statistics represent two complementary characteristics of effect size for analysed predictor groups. To interpret the direction of predictor effects, one must compare the positions of their arrows with the arrows of response variables, in an ordination biplot diagram (Supporting information).

Stepwise selection of predictors
For each set of response data (active defenders and cryptic breeders in 1993 and 2013 in four separate models), we first chose a minimum adequate subset of predictors in each group, using stepwise selection. Because this method presents an inherent danger of overestimating the size of the set of required predictors, we applied measures limiting such a bias based on the suggestions of Blanchet et al. (2008): (a) before selecting from a group of predictors, a permutation test of significance was performed using the whole group and, if the joint effect was not significant (p > 0.05), no predictors were selected; (b) the primary stopping criterion in the predictor selection was based on a partial permutation test (of the candidate predictor effect in addition to the already selected predictors), and the estimated p-values were adjusted in the case of spatial predictors (where a large pool of candidate predictors substantially inflates the type I error) by transforming them into false discovery rate (FDR) estimates (Verhoeven et al. 2005); (c) an additional stopping criterion was based on calculating R 2 adjusted for constrained ordination using the whole group of predictors, which served as a reference value. When this reference value was exceeded by R 2 adjusted of the selected predictor subset, the selection was stopped. As an open water area reduced the nesting opportunities for most species (e.g. on plant deposits on floating islets), we preferably selected predictors other than 'water', if they emerged among the significant candidate predictors.

Testing predictor effects and visualizing results
All tests of hypotheses concerning multivariate data (in the constrained ordination framework) were performed in Canoco 5 software (ter Braak and Šmilauer 2012), using Monte Carlo permutation tests based on pseudo-F statistics, because their traditional, more parametric counterparts are not available here (ter Braak and Šmilauer 2012, p. 72). We presented the pseudo-F statistic values with a simplified 'F' label. As we were comparing just two sampling times separated by 20 years and with a slightly different way of collecting data, we analyzed each sampling time separately rather than focusing on an (inappropriate) interpretation of the time × predictors interaction terms in a dataset pooled across the two years. In addition, we evaluated the community patterns separately for the two guilds, based on a priori assumed hierarchical nature of their relationship, with both groups of species affected by spatial and environmental predictors, but with cryptic breeders further responding to the presence of active nest defenders. Because a significant effect of a predictor on community composition does not necessarily mean that all species are correlated with that predictor, we emphasize in the ordination diagrams (where appropriate, i.e. for constrained ordinations with a single or few predictors) the responsive taxa, based on the suggestions derived from t-value biplots (ter Braak and Šmilauer 2012, p. 226).

Validation of multivariate analyses
To examine the robustness of our conclusions obtained by the variation partitioning with the multivariate RDA (above), we evaluated these multivariate data using an alternative approach, based on the modelling of community variation as described in Jamil et al. (2013). Their models were defined to study relationship between species functional traits and site environmental descriptors, but they could be easily simplified for the present purpose. We have used generalized linear mixed-effect models (GLMMs) with assumed Poisson distribution of random variation and random effects of both sampling square identity and species identity, accounting respectively for the dependency among the observations of individual species at the same location as well as the dependency among observations made on a particular species across different squares. We have used only linear terms to describe the effects of environmental variables to match the modelling approach of the multivariate RDA.
We estimated GLMMs separately for the two sampling years, for the two species groups (active nest defenders versus cryptic breeders) and also for individual predictor groups or their combinations. The predictors were a priori separated into two to three groups, as described in the earlier section Response data and predictor groups. In the case of spatial predictors, however, we used an alternative, more simple approach, where the description of spatial variation is based on the selection of linear or smooth terms (cubic splines with two or three degrees of freedom) of the latitudinal and longitudinal coordinates of sampling squares. Similarly to our primary approach, predictors for GLMMs were selected in each group separately (here using parametric χ 2 tests, with type I error estimates adjusted by their transformation to false discovery rates, Verhoeven et al. 2005) and then we performed variation partitioning (using individual, pre-selected groups of predictors or their combinations) based on adjusted R 2 estimates of the fitted GLMM (Nakagawa et al. 2017). The GLMMs were fitted using lme4 package (Bates et al. 2015) and the adjusted R 2 quantified with MuMIn package (Barton 2020) in the R software (<www.r-project.org>).

Habitat structure
The most dominant habitat in both years was moss, followed by open water and herb cover (Fig. 2). While the proportion of moss remained stable from 1993 to 2013, the water area decreased by 14% on average (median 6%), and the herb cover increased by 11% (median 0%). The observed Moran's I for the habitat scores (I 1993 = 0.190, I 2013 = 0.318) was significantly higher (p < 0.05) than expected for a random pattern [I = −0.006 ± (SD) 0.007 for both years], suggesting the aggregation of particular habitats throughout the area. In general, the water depth gradient from the marshland edge to open water reflected the habitat transition from waterlogged forest to the west through herb cover and moss plains towards the open water with islets of vegetation to the south and east. This was evident in both years, but in 2013 the herb cover was spread out more to the west and north at the expense of moss, whereas the moss plains shifted at the expense of decreasing water covered areas to the southern and central parts of the study area.

Community composition and the spatial pattern of the nests
In 1993, we identified 25 species with 665 confirmed or probable nesting attempts, while in 2013 the number of species dropped to 19 with 188 confirmed nests (Supporting information). Out of eight species of active nest defenders in 1993, only four remained in 2013, with greatly reduced abundances. The most notable differences between the two years included the disappearance of three previously common active defenders: white-winged black tern Chlidonias leucopterus, little gull Hydrocoleus minutus and common gull Larus canus. In contrast, the Mongolian gull was the only active defender and potential nest predator that had newly settled to breed in the area. In addition, the number of 17 species classified as cryptic breeders in 1993 decreased to 14 species in 2013. Most cryptic breeders declined (8 species) or disappeared (4 species), while only one species was new.
Active nest defenders started to breed significantly earlier than cryptic breeders in both years (1993: 29 May as the median laying date for 107 measured nests of active defenders, and 6 June for 45 nests of cryptic breeders, Mann-Whitney test: W = 1531, p = 0.0004; 2013: 29 May for 44 nests of active defenders, and 15 June for 63 nests of cryptic breeders, Mann-Whitney test: W = 410, p < 0.0001).

Active nest defenders in 1993
After controlling for latitudinal trend, three habitat descriptors contributed significantly to the spatial nest pattern of active defenders (test of significance on all canonical axes: F = 3.9, p = 0.002): herb cover, moss and water depth (Table 2). Specifically, active nest defenders avoided breeding on floating islets above deep water and some waders preferred herb cover and moss (Supporting information). However, the unique contribution of habitat predictors to the total explained variation was < 0.1%, and a substantial part of their explanatory power was shared with spatial predictors (7.7% of the total variation; Table 3, Fig. 4). We found a significant effect of latitude (F = 12.4, p < 0.001) explaining 6.5% (R 2 adjusted ) of the total variation in species composition. In addition to the latitudinal gradient, there were 28 significant spatial eigenvectors (test of significance on all canonical axes, F = 9.4, p = 0.001; Supporting information), which uniquely contributed 51.6% of the total variation in the community composition of active nest defenders in 1993 (Table 3, Fig. 4). The effect of space (expressed as a mean square statistics for unique contributions) was > 8 times larger than the effect of habitat. We obtained consistent results using the variation partitioning in the alternative GLMM approach, with larger spatial effect and substantial overlap between both predictor groups (Supporting information). The spatial predictors (indicating conspecific aggregation) played a similar role in the nest clustering of three larids -black-headed gull Chroicocephalus ridibundus, little gull and common tern Sterna hirundo, while the nests of whitewinged black tern were spread across the area (Supporting information).

Active nest defenders in 2013
After controlling for latitudinal trend, the choice of nest site was explained by habitat (test of significance on all canonical axes: F = 4.0, p = 0.006), in particular by moss and herb cover (Table 2) exclusively explaining 4.6% of the adjusted total variation in species data (Table 3, Fig. 4). The appearance of non-colonial Eurasian curlew Numenius arquata correlated with the amount of herb cover, while semi-colonial northern lapwing Vanellus vanellus showed a positive response to moss cover. Larids such as common tern, black-headed gull and Mongolian gull avoided both these habitats (Supporting information), preferring floating islets above deep water. However, the shared effect (3.7%, Table 3) of the selected habitat descriptors and latitude (the only selected spatial predictor -Supporting information) suggests a considerable overlap of these two explanatory components (Fig. 4). Latitude (F = 7.8, p = 0.001) explained (R 2 adj ) 4.0% of the total variation in species data and was included in the further analysis. The unique explanatory power of spatial predictors (indicating conspecific aggregation) was negligible in comparison with the effect of habitat. This was also supported by the mean square statistic, where the average strength of the habitat predictors exceeds the selected spatial predictors more than three times. Our alternative analytical approach with GLMM (Supporting information) also supports substantial overlap between the effects of habitat and spatial predictors and selects the herb cover and moss cover as effective predictors, but leave no substantial unique (partial) effect of the habitat descriptors. Table 2. Partial RDA comparing the effects of habitat descriptors on nest abundance of active defenders in 1993 and 2013 using forward selection procedure with latitude used as a covariate. The 'Explains %' values show the explanatory contribution of selected predictors to total variation. 'F' (Pseudo-F) and 'p adj ' values refer to a partial test performed at the moment of (considered) predictor entry into the set of selected predictors.

Predictor
Explains % F p adj  Figure 3. Occurrences of nesting birds in the squares under study within a terrain map projection. Small, medium and large circle sizes indicate 1-2, 3-9, ≥ 10 nests in total, respectively, with slices representing proportion of active defenders (A), active defenders and potential nest predators (AP) and cryptic breeders (C).

Cryptic breeders in 1993
Variation partitioning with three groups of predictors -habitats, association with active defenders, and spatial predictors -showed that habitats contributed significantly to explaining nest site choice (test on all canonical axes: F = 3.0, p < 0.002).
Moss was the only significant habitat selected in the forward selection procedure (F = 9.7, p < 0.001, Supporting information), moderately preferred by some waders, e.g. common snipe Gallinago gallinago and wood sandpiper Tringa glareola, whereas it was clearly avoided by grebes (Supporting information). Active nest defenders contributed significantly to the nest site choice of cryptic breeders (test on all canonical axes: F = 2.6, p = 0.036), with particular importance of their abundances (Table 4) seen in a tight association with greebes and Aythya ducks (Fig. 5). Anas ducks and waders, such as wood sandpiper and common snipe, were positively associated with the presence of common gulls (Fig. 5). Among spatial predictors indicating conspecific aggregation, three other spatial eigenvectors (Supporting information) were chosen as significant (test on all canonical axes: F = 1.5, p < 0.001) in addition to latitude. Spatial predictors of cryptic breeders had a larger unique contribution (6.6%) to the total (adjusted) variation than for active defenders predictors (2.5%) and habitat predictors (0.5%) (Table  5). Moreover, spatial predictors shared a large part of their explanatory power with habitats (4.5% of the explained variation as a sum of f and g fractions) but shared only a small part of their explanatory power with active defenders (1.4% as the sum of the e and g fractions) (Fig. 4). The unique effect of active defenders (fraction b) was significant (explaining 2.5% of the total variation, F = 3.3, p = 0.013; Table 5). Overall, the average strength of the unique effects of active defenders and spatial predictor groups was about twice as great as the strength of habitat predictors, as judged by the mean square statistics. By analyzing the variation partitioning using the alternative approach with GLMM, we obtained similar results for spatial effects and their overlap with the effect of habitat. However, the effect of habitat (water depth and herb cover) was stronger in this analysis while the contribution of active defenders was non-significant (Supporting information).

Cryptic breeders in 2013
Variation partitioning revealed significant unique effects of habitats and active defenders. The effect of habitats (test on all canonical axes: F = 4.6, p = 0.004) included moss Table 3. Variation partitioning results comparing adjusted variation in the community of active defenders explained by unique effects of habitats, spatial predictors and the fraction shared by habitats and spatial predictors (according to scheme in Supporting information and Fig. 4). Note that the 'MS' (mean square) statistic uses the nonadjusted estimate of explained variation, divided by corresponding degrees of freedom ('DF'). For independent effects of the fractions combined see Supporting information.   and herb cover as two significant predictors selected in the forward selection procedure (both F ≥ 6.9, p ≤ 0.003, Supporting information). These habitats were preferred by some waders and were avoided by ducks and grebes (Supporting information). Three predictors describing the effects of active defenders (test on all axes: F = 34.8, p < 0.001) -nest count for active defenders, presence of nests of an active defender and the presence of Mongolian gull nests -were selected as highly significant (Table 4). In particular, the unique effect of active defenders (fraction b) was significant and explained 32.4% of the total community variation (F = 30.1, p < 0.001 in 2013; Table 5). Ducks and grebes, but not waders, built their nests preferably near to active nest defenders, with the exception of Mongolian gulls (Fig. 6). Habitat predictors had a small (yet significant) unique contribution (1.3% of the total explained variation) compared with the unique effect of active defenders, which is further documented by the mean square statistics for those unique effects (Fig. 4). Shared effects were important among all three groups of predictors, with the highest contribution from the overlap between all three predictor groups (3.7%). Latitude was chosen as the only significant spatial predictor (F = 9.4, p < 0.001, R 2 adj = 4.9% of the total community variation). Using the alternative modelling approach with GLMM, we obtained consistent results including significant contribution of the three predictors describing the effect of nest defenders together with a significant effect of habitat (moss and wood) as well as the stronger explanatory effect of active defenders compared with the year 1993 (Supporting information).

Discussion
The diverse habitat mosaic of large wetlands supports the formation of rich bird communities usually dominated by colonial birds like gulls and terns, accompanied by other less conspicuous species, such as ducks, waders, grebes, divers and rallids (Baschuk et al. 2012, Pagel et al. 2014. We have demonstrated that habitat, association with other species, and additional spatial effects indicating intraspecific aggregation contributed effectively to the spatial pattern of nests within the wetland bird community in both years of research. Because our study could not be performed experimentally with appropriate system manipulation, the observation-based results cannot confirm causality of our findings. Therefore, we used two robust analyzes to separate the effects of defined predictors from other spatial effects. We performed variation partitioning by two alternative approaches, the redundancy analysis and GLMM. Both these approaches provided fundamentally consistent results in terms of relative effects of habitat, spatial predictors and sharing (overlaps) of these effects, although the detailed outcomes inevitably differed to some extent. This is likely affected by the changed nature of spatial descriptors, where in the GLMM approach we selected linear and smooth terms focusing at the description of more global trends. However, the comparison between the years remains consistent, with a clear contribution of heterospecific aggregation and reduced habitat effect in 2013.

Habitat effect
Habitats, thought to be a 'firsthand' predictor explaining species distribution in the area, played a role in both groups of active defenders and cryptic breeders. Detailed inspection revealed that waders preferred herbs away from deep water, while the opposite was true for larids, grebes and ducks, which is in concordance with previous studies (Duebbert et al. 1983, Frederick andCollopy 1989). The general avoidance of moss by most birds was probably associated with excessively low and uniform cover, which provided little opportunity to conceal nests, thus enabling easy access for both avian and mammalian predators. Only the northern lapwing, which has cryptic eggs, prefers a good view from the nest (Cramp and Simmons 1983), and is active in deterring predators from the nest vicinity (Elliot 1985) selected moss as a suitable habitat for nesting.
Although the habitat effect was evident in active defenders in 2013, there was a large overlap with spatial effects, which might be interpreted in two mutually non-exclusive ways. First, bird nests showed spatial aggregation on the same scale as did the habitats, and, second, the birds clumped the nests due to conspecific attraction in addition (or with a lesser regard) to the effect of habitat. There was a disparity in the relative importance of habitat descriptors for active defenders in 2013, when comparing the results from the RDA-based variation partitioning and the alternative GLMM approach. This can be explained by the different way the spatial effects (with larger relative effect in the GLMM approach) were modelled. The spatial configuration of nests of active defenders in 2013 was more effectively modelled by the smoothing splines in the GLMM approach (using 4 DFs), as in the dbMEM approach none of the spatial eigenvector predictors but only a linear trend in latitudinal direction (using 1 DF) was selected.
On the other hand, for cryptic breeders the contribution of habitat was lower than the clumping with active defenders especially in 2013. Habitat explained only some general patterns, such as occurrence of some waders in moss and grebes and ducks in water-associated vegetation. Therefore, despite the obvious habitat effect on nest site selection in many species, the results indicate that the role of habitat should always be extended with the effect of presence of other birds, either conspecifics or other species.

Conspecific aggregation
The unique effect of spatial descriptors indicated that nesting preferences went beyond what was explained through the habitat composition, particularly in 1993. The latitude alone could not have a biological significance on the local scale relevant to this study. It was probably a proxy of the most pronounced gradient in the study area, which included a transition from taiga edge with moss and herbs to open water with islets of floating vegetation and the increasing depth of the wetland further from the shore and closer to the open lake. We suggest that much of the effect of other spatial predictors in active defenders probably arose from conspecific attraction resulting from the need for colonial breeding (Rolland et al. 1998). Active defenders behave conspicuously during the incubation period Simmons 1983, del Hoyo et al. 1996), thus the easier detectability of their nests by predators, directly or through incubating parents (except of well camouflaged Eurasian curlew, own observations) requires to be compensated, e.g. by shared protection (Gochfeld 1984, Šálek andŠmilauer 2002). The three common larids -black-headed gull, little gull and common tern -shared habitat as well as nesting space in 1993, as a result of heterospecific aggregation. As we have shown, they placed their nests preferably on islets surrounded by water, away from marsh edges, and this prevented the access of predators from the taiga forest (Sládeček et al. 2014). Although the aggregating behavior was conspicuous in the abundant gulls and terns, it was less obvious in uncommon active defenders  like Eurasian curlew and northern lapwing. Both of these less common species combine active nests defense with nest crypsis, either through incubating parents (Eurasian curlew; Nethersole-Thompson and Nethersole-Thompson 1986), or through camouflaged eggs (northern lapwing; Simmons 1983, Šálek andCepáková 2006), allowing them to safely inhabit more exposed habitats e.g. near marshland edges, without the conspicuous clustering tactic.

Heterospecific aggregation
In addition to the conspecific aggregating behavior of most active nest defenders, our results showed that the cryptic breeders tended to share colonies with active defenders in addition to their habitat choice. Such behavior was previously suggested for particular species, e.g. for bar-tailed godwits Limosa lapponica seeking association with aggressive whimbrels Numenius phaeopus (Larsen and Moldsvor 1992) or for ducks (Anas spp. and Aythya spp.) breeding nearby actively defending black-headed gulls and little gulls (Väänänen et al. 2016, Pöysä et al. 2019. Species lacking a pronounced active nest defense may benefit from the protective umbrella against predators provided by the actions of active nest defenders. This may result in generally higher nesting success, as has been described in black-necked grebes Podiceps nigricollis sharing colonies with black-headed gulls (Fiala 1991), or in white-tufted grebe Rollandia rolland and silvery grebe Podiceps occipitalis within colonies of brown-hooded gull Chroicocephaus maculipennis (Burger 1984). Our study extends these previous findings and indicates that clumping of breeding birds can be a much more general phenomenon across different habitats occupied by the wetland bird community. The results suggest that particularly the abundance of active nest defenders in addition to their presence can influence nest site choice of cryptic breeders, especially ducks and grebes. In the mixed species colonies with abundant defenders the advantage for cryptic breeders may comprise their lowered risk of nest predation due to the diluting effect of incidental predation (Ringelman 2014), strengthened by a higher predation risk of more conspicuous nests of active defenders (Sládeček et al. 2014).

Possible role of umbrella and predator species
The spatial predictors (indicating the conspecific attraction) played a more important role in 1993 while the effect of active nest defenders (explained by the heterospecific attraction) dominated in 2013. Specifically, the distribution of active defenders predicted the distribution of cryptic breeders in both years, but this relationship was more pronounced in 2013 within a less diverse and less abundant bird community.
A complete disappearance of colonially breeding whitewinged black terns and a new settlement of Mongolian gulls in 2013 might have an effect on spatial nest patterns of cryptic breeders. White-winged black tern is a common active nest defender with no apparent habitat preferences, able to breed on swampy standing water and transitional or fluctuating marginal flooding areas regardless of water depth (Cramp 1985, Goławski et al. 2016. Its wide occurrence across the wetland in 1993 could attract cryptic breeders to various places with less emphasis on habitat and affect the role of all tested predictors in 1993. On the other hand, the recently established population of Mongolian gull in the study area might have had a negative effect on cryptic breeders. The Mongolian gull is colonially breeding nest defender and prefers similar habitats like black-headed gull and common tern. However, despite this similarity it discouraged rather than attracted other species. A probable reason for this opposite relationship is that this aggressive species of bigger size could increase the risk of predation on the eggs and/or the young of neighbors (Thomas 1972, Skórka et al. 2014, Sládeček et al. 2014, Mel'nikov 2019. The results of a long-term study in more southern parts of Lake Baikal between 1972 and 2001 showed that the breeding success of marshland birds dropped to 20-40% probably as a result of increased nest predation by naturally dispersed big gulls (here common gull and Mongolian gull) which suffered from reduced production of fish waste (Meľnikov 2010). Thus, cryptic breeders avoided nesting near the nests of these gulls. Our episodic direct observations, images from camera-traps and collected remains of depredated eggs suggest that the carrion crows and Mongolian gulls were the dominant nest predators within the study area. We also confirmed the red fox as a nest predator in the area, but (according to the analysis of food remains at one fox den) it probably sought food mostly near the forest edge (Sládeček et al. 2014).
Our results suggest that cryptic breeders sought more protection from active defenders within a poorer community in 2013 with a greater predation risk. The loss of one dominant umbrella species, loss of water habitat with islets protected against mammalian predators and the settlement of a new aggressive avian predator in the community may have been the main drivers of this effort. Alternatively or additionally, all individuals from smaller populations in 2013 were better able to nest nearby active defenders compared with 1993.

Alternative explanations of spatial effects
There can be additional (and mutually non-exclusive) explanations for the non-specified spatial effects in addition to the effects of conspecific and heterospecific attraction. For example, individual variation in the propensity to aggregate or to avoid neighbourhoods mirroring individual differences in risk-taking behaviour (Steinhoff et al. 2020) may play a role in the choice of nest sites. First, related individuals may aggregate more (Caro 2005) and thus locally increase the magnitude of spatial effects. Second, a previous breeding experience may contribute to the decision making of individuals (Doligez et al. 2002, Sebastián-González et al. 2010b. We should also note that we analysed the data on an arbitrary scale of 100 × 100 m that does not allow us to detect the finer habitat parameters that could potentially play a role in some species (Dray et al. 2012). Although the dbMEM method employed in our study comprises all spatial scales present in our data, the results are still limited by the spatial resolution of the analyzed data.

Conclusions
Even though correlative analysis cannot reveal causality, in addition to the availability of suitable habitats, conspecific and heterospecific aggregation played a significant role in explaining spatial nest patterns of birds inhabiting a large pristine Siberian wetland. The magnitude of the aggregation effect throughout the community may fluctuate depending on the number and abundance of species that play a major role either as nest defenders or predators. The colonially breeding nest defenders with their protective umbrella effect may become key species of the community due to their attractiveness for many other species that may not be primarily colonial. Some of them can resolve the trade-off between breeding in heterospecific colonies even in suboptimal habitats and solitary breeding in optimal habitats. Their decision may vary with the presence and/or numbers of species that increase the risk of nest predation. Their presence or abundance may be just as important or even more motivating than the habitat quality itself and the disappearance of one or few key species from the community may initiate the community decline. Therefore, heterospecific aggregation can be a crucial attribute and indicator of long-term maintenance of prosperous and diverse communities of wetland birds. Integration of factors supporting the breeding of key species, such as small larids, may became key targets for comprehensive conservation measures in large wetlands.