Agricultural intensification is one of the major factors driving biodiversity loss. However, most studies in human-dominated landscapes have used taxonomic diversity in their analysis, ignoring evolutionary relationships. Consequently, the relationship between landscape structure and phylogenetic diversity is not well understood. Here, we tested the hypothesis that landscape heterogeneity is positively related to bird phylogenetic indexes of diversity and structure, leading to over-dispersed phylogenies in very heterogeneous landscapes. We analyzed phylogenetic responses in interfaces between forest edges and anthropogenic matrices (forest-pasture and forest-eucalyptus) using generalized linear mixed models. We also compared these indexes between land covers to assess which one best preserves the phylogenetic history of communities. We used both traditional phylogenetic indexes and those corrected for species richness. Our results showed that phylogenetic diversity varied significantly between land cover types, but this did not occur when we removed effects associated with species richness, suggesting that all land covers preserve similar levels of evolutionary history. Additionally, our best models showed a positive relationship between landscape heterogeneity and bird phylogenetic indexes of diversity and structure, but the strength of these relationships may be land-cover dependent. In summary, our work highlights the influence of landscape heterogeneity on the phylogenetic diversity and structure of bird communities, reinforcing the need for its incorporation into conservation-based studies.

1. Study sites and landscape metrics

We created land-use maps for this study area using satellite images of high resolution (ArcGIS 10.3 base map imagery, Digital Globe satellites 2010–2011, 1:4,000; Figure S2). We made this map with manual digitization based on visual interpretation of patch differences in color, texture and shapes. We considered 14 different land covers or human land uses: old-growth forest, pasture, eucalyptus plantations, second-growth, wetland, cropland (mainly maze), sugarcane, water bodies, urban areas, rural homesteads, urban or suburban homesteads, paved roads, buildings, and bare soil (for more details, see Barros et al., 2019). Based on this map, we chose 34 sample landscapes to represent the regional gradient of compositional landscape heterogeneity. Of these, 16 had their points located across forest-eucalyptus interfaces and 16 across forest-pasture interfaces. The eucalyptus plantations were, in general, homogenous, lacking a native vegetation understory, while the pastures were not exposed to intensive grazing. In addition, we used two control landscapes in continuous forest areas so as to include forest specialist species in the phylogenetic tree. This allowed the inclusion of forest species that originally inhabited our study area but which are now only seen in large forest blocks.

Landscape sites were defined using 1.2-km buffers around the centroid of the two sampling sites (one at forest edge, the other in pasture or at eucalyptus edge). We chose this spatial scale based on previous evidence from multiscale analysis of bird responses to landscape structure (Barros et al., 2019). We also defined a 2-km minimum distance between sampling landscapes to avoid recounting the same individuals in different landscapes.

For each landscape, we calculated the compositional landscape heterogeneity (Fahrig et al., 2010) via a Shannon diversity index using Fragstat v.4 software (Mcgarigal et al., 2012). This index is based on the number of physical land cover types or human land uses and their evenness within the landscape. When the value of the Shannon index is zero, the landscape contains only one patch (i.e., homogeneous landscape). Landscape heterogeneity increases as the number of different patch types increases and/or the proportional distribution of area among patch types becomes more equitable.    

2. Bird sampling

We used point counts, with a 50-m survey radius, to sample bird communities (Sutherland et al., 2006). We selected paired forest-matrix ecotones as sample sites in each landscape, one in forest edge and other in the adjacent pasture or eucalyptus plantation. The paired sampling sites were located around 70-100 m from the edge (140-200 m from each other, Figure S2), while control sampling sites consisted of a single point at each location, situated in the interior of continuous forests, and at least 1 km from the edge. We sampled each site for ten minutes, three times, on different days during the first three hours after sunrise, during two consecutive breeding seasons, totaling 30 minutes per site. We sampled half of sampling sites from September 2014 to January 2015, and the remainder from October to December 2015. Only birds obviously using the habitats were recorded (e.g. birds flying overhead were not included). For each sampling site, data recorded from the three replicates were combined into a single community database, except for the abundance data. We set the abundance data of each sampling site as the highest value recorded in one single day. This is a conservative value adopted to avoid overestimated data.

3. Bird phylogenetic trees

We considered all species recorded in the 34 landscapes as the regional pool of species. Then, to construct the phylogenetic tree, we used the phylogeny database of bird species built by Jetz et al. (2014), available for download at http://birdtree.org/, using the Hackett All Species option as the source of trees. From the 1,000 trees pruned for our full set of species, we created a consensus tree using the Mesquite 2.75 program (Maddison & Maddison, 2010).

4. Phylogenetic metrics

To evaluate the phylogenetic metrics for each community, we used the “picante” package (Kembel et al., 2010) in R statistical software (R Core Team 2018). We chose six phylogenetic metrics to represent the phylogenetic diversity and structure of sampled bird communities:

1. Phylogenetic diversity (PD): described as the total sum of phylogenetic history, this is measured through the total branch length of a phylogeny representing the species in a community (Faith, 1992);

2. Standard effect size of PD (SES.PD): PD of communities may be correlated with their species richness (SR, or number of species; Swenson, 2014). Thus, the effect of SR can be removed by comparing the PD values of studied communities with that of communities of equal species richness generated by null models drawn randomly from the regional species pool;

3. Mean nearest taxon distance (MNTD) abundance-weighted: MNTD is weighted by abundance, and represents the average phylogenetic distance between an individual and the most closely related non-conspecific individual (Webb, 2000; Webb, 2002); high MNTD indexes indicate the co-occurrence of distantly related species within communities (phylogenetic over-dispersion), while low levels indicate the co-occurrence of closely related species (phylogenetic clustering);

4. Standard effect size of MNTD (SES.MNTD) abundance-weighted: MNTD may also be correlated with SR. Communities with higher than expected MNTD values for a given SR indicate the co-occurrence of distantly related species (phylogenetic over-dispersion or SES.MNTD > +1.5), while low values indicate the co-occurrence of closely related ones (phylogenetic clustering or SES.MNTD <-1.5). The SES.MNTD index is a better metric for limiting similarity relationships (phylogenetic over-dispersion);

5. Mean pairwise distance (MPD): represents the average phylogenetic distance between all pairwise species combinations present in a community, and it is influenced by relationships in deep evolutionary time (Webb, 2002). High MPD indexes indicate the co-occurrence of distantly related species within communities (phylogenetic over-dispersion), while low levels indicate the co-occurrence of closely related species (phylogenetic clustering);

6. Standard effect size of MPD (SES.MPD): MPD may be correlated with SR. Communities with higher than expected MPD values for a given SR indicate co-occurrence of distantly related species (phylogenetic over-dispersion or SES.MPD > +1.5), while low values indicate co-occurrence of closely related ones (phylogenetic clustering or SES.MPD <-1.5). The SES.MPD index is a better indicator of environmental filter relationships (phylogenetic clustering; Webb, 2000; Kraft, Cornwell, Webb, & Ackerly, 2007);

To determine whether SES values differed from the community expected by chance, we compared observed values between individuals to expected SES values for 999 communities using an independent swap algorithm (Gotelli, 2000). Additionally, to assess whether such traits were conserved across the phylogeny (Data S1), we calculated the phylogenetic signal (Blomberg, Garland, & Ives, 2003) for six bird functional traits. For this, we conducted a phylogenetically independent contrast analysis using the “aotf” function in Phylocom software (Webb, Ackerly, & Kembel, 2008).  

5. Statistical analysis

To assess whether the phylogenetic diversity indexes and species richness differed among all land covers, we used one-way ANOVA, and to compare differences among matrices and forest edges, we used paired-sample t-tests (Rezende & Diniz-Filho, 2012). To compare only matrices (pasture and eucalyptus plantations), we performed simple t-tests. Results were interpreted via P-values where P < 0.05 indicates significant differences between land covers. Statistical analyses were performed using R Statistical Software (R Core Team, 2018). 

We used generalized linear mixed models (Table 1) to analyze the relationship between phylogenetic and landscape metrics (Glmm; function “lmer”, package “lme4”; Bates, Mächler, Bolker, & Walker, 2015; Bolker et al., 2009).  For each phylogenetic metric of diversity (PD and SES.PD) and structure (MPD, MNTD, SES.MPD and SES.MNTD), we produced simple and additive models using spatial heterogeneity and land cover as predictive variables. We fitted these models using the landscape identification (a code for each landscape) as a random effect to account for spatial dependence of the sampling points present in the same landscape. In some models, we also tested, as random effects, the landscape heterogeneity and land cover. Next, we ranked these models and, using the Akaike Information Criterion corrected for small samples (AICc), estimated which ones best predicted the phylogenetic diversity and structure of bird communities. We considered best models to be those with the lowest AICc values. Models with ΔAICc differences less than two were considered as equally plausible to explain observed patterns (Martensen, Ribeiro, Bankes-Leite, Prado, & Metzger, 2012). For each model, we also calculated the Akaike weights wi (range from 0 to 1; the highest values were considered the most plausible models). We calculated all these indexes using the AICctab function of the “bbmle” package (Bolker, 2017).

After selecting the most plausible models (ΔAICc < 2, wi > 0.1), we tested for spatial autocorrelation in the residual distribution using Moran's I test (“DHARMa” package in R). Because of spatial dependence, we used spatial regression models that included the spatial effect (function “fitme”, package “spam”; Rousset & Ferdy, 2014; Table 2). In these models, we excluded landscape identification as a random effect as these regression models already took into account the spatial pattern present in the data. Finally, we visually checked the model fits and residual distributions and selected the best ones that controlled the spatial effect (Figure S3). All these analyses were performed using the R statistical software package (R Core Team 2018).