Multi‐marker DNA metabarcoding reveals spatial and sexual variation in the diet of a scarce woodland bird

Abstract Avian diet can be affected by site‐specific variables, such as habitat, as well as intrinsic factors such as sex. This can lead to dietary niche separation, which reduces competition between individuals, as well as impacting how well avian species can adapt to environmental variation. Estimating dietary niche separation is challenging, due largely to difficulties in accurately identifying food taxa consumed. Consequently, there is limited knowledge of the diets of woodland bird species, many of which are undergoing serious population declines. Here, we show the effectiveness of multi‐marker fecal metabarcoding to provide in‐depth dietary analysis of a declining passerine in the UK, the Hawfinch (Coccothraustes coccothraustes). We collected fecal samples from (n = 262) UK Hawfinches prior to, and during, the breeding seasons in 2016–2019. We detected 49 and 90 plant and invertebrate taxa, respectively. We found Hawfinch diet varied spatially, as well as between sexes, indicating broad dietary plasticity and the ability of Hawfinches to utilize multiple resources within their foraging environments.

understanding of the complex interactions which birds have within their environment . In turn, this provides essential information for the conservation and management of avian species and their associated habitats (O'Donnell et al., 2012;Ontiveros et al., 2005). Intraspecific dietary niche separation is important for ecological dynamics (Cloyed & Eason, 2017) and has a number of drivers. Environmental factors such as habitat type and quality may influence diet due to potential changes in species composition between habitat types (β-diversity; Shutt et al., 2019). Within forests, it has been found that invertebrate species richness can differ between tree taxa (Murakami et al., 2008;Shutt et al., 2019). Intrinsic factors, such as sex can also be a source of dietary variation, driven by morphological and behavioral differences, as well as differing nutritional requirements (da Silva et al., 2020;Mata et al., 2016). Environmental factors can also interact with intrinsic factors to influence dietary variation between demographic groups, such as dietary variation between sexes only being apparent during the breeding season, when males and females have differing reproductive demands (da Silva et al., 2020;Durell et al., 1993).
The Hawfinch (Coccothraustes coccothraustes) is one of many woodland specialists to have shown major declines in the UK (Kirby et al., 2018). Hawfinch breed across the Palearctic, with Britain its western range limit (Kirby et al., 2015). Hawfinch ecology is poorly understood, and in Britain, Hawfinch are now too rare to have regular status assessments by national annual monitoring schemes (Kirby et al., 2015). Instead, population change is inferred from distribution data compiled from bird atlas surveys (Balmer et al., 2013). These atlas data indicate a 76% reduction in the number of 10 km squares occupied between 1968 and 2011 (Kirby et al., 2015(Kirby et al., , 2018. This decline was previously evidenced by Langston et al. (2002), who estimated a 40% population decline between the mid-1980s and late-1990s. Hawfinch show a very localized distribution within the UK, with population strongholds exhibiting a strong westerly bias (Kirby et al., 2018), and the factors causing their decline are unknown. Landscape modification, decreased invertebrate abundance, and changes in woodland management have been suggested (Fuller et al., 2005;Kirby et al., 2018). Further potential contributory factors may include under-planting of ancient woodland with conifers in the 1970s, and a storm in 1987 which caused the loss of many cherry trees, which are thought to be an important food resource for Hawfinch (Kirby et al., 2018;Spencer & Kirby, 1992).
Hawfinch dietary studies are limited, with all previous information obtained through visual observation (Mountford, 1957;Newton, 1967). Hawfinch are thought to be dietary specialists adapted to utilize large seeded tree species due to their large and powerful beak (Mountford, 1957). During the breeding season (typically from April to June), Hawfinch have been observed to feed on buds and flowers of cherry (Prunus sp.), hornbeam (Carpinus betulus), beech (Fagus sylvatica), and Wych elm (Ulmus glabra; Mountford, 1957). Hawfinch have also been observed to incorporate invertebrates into their diet during the breeding season, consuming Lepidoptera, Coleoptera, Hemiptera, Annelida, Gastropoda, and Araneae (Mountford, 1957).
Metabarcoding requires minimal a priori knowledge of the target organism's diet (Alberdi et al., 2017;Valentini et al., 2009), and a wide range of ingested taxa can be identified to fine taxonomic levels . Morphology-based methods often record ingested taxa at a coarse taxonomic resolution, missing subtle differences in the taxa consumed (Alberdi et al., 2017). This limits the opportunities to make fine scale inferences relating to species' ecology (da Silva et al., 2020;Mata et al., 2016). The accurate identification of components within an omnivorous diet is, however, still considered challenging (da Silva et al., 2019;De Barba et al., 2014;Tercel et al., 2021). Due to the costly, laborious, and taxonomically demanding nature of exploring omnivorous diet, studies attempting to elucidate all dietary aspects are rare (Pompanon et al., 2012;Robeson et al., 2018;Tercel et al., 2022).
Metabarcoding studies have been used to successfully identify details of the diet of farmland birds (Cabodevilla et al., 2021) and insectivorous species (Evens et al., 2020;McClenaghan et al., 2019;Mitchell et al., 2021). Metabarcoding therefore has the potential to improve our knowledge of the diet of omnivorous woodland passerines. Despite this, studies focusing on omnivorous birds are rare (but see da Silva et al., 2020;Spence et al., 2022;Tang et al., 2022 (Clements, 2013;Kirby et al., 2018). All feed sites were baited with supplementary sunflower (Helianthus sp.) seeds, provided ad libitum from December to July. Study sites were broadly typical of British mixed broadleaved woodland, with sites in the Wye Valley and north Wales dominated by beech, oak, and ash (Fraxinus excelsior). The study site located in East Anglia was a mixed woodland consisting of lime (Tilia sp.), ash, and maples (Acer sp.). The New Forest site was dominated by oak, with an understorey flora comprising of Holly (Ilex sp.) and bramble (Rubus sp.). All site locations are approximate for anonymity.
Hawfinches were caught by trained bird ringers, operating under British Trust for Ornithology (BTO) approved licenses using either mist or whoosh nets. Sex was determined through secondary wing feather coloration following Svensson (1992). All birds were aged as adults or juveniles, with birds aged through width of the outer web of the sixth primary feather and sharpness of the greater alula feather, following Fornasari et al. (1994). Juveniles were defined as young of that year which had already fledged. Hawfinch show strong sexual plumage color dimorphism of secondary feathers, enabling the sexing of individuals in both the non-breeding and breeding season (Svensson, 1992). We placed Hawfinches within individual clean, paper bags which were then placed inside a cotton bag and left for 10-20 min until the bird defecated. To avoid excessive stress, if birds had not defecated within this time frame they were processed and released. We removed each fecal sample using a new, clean, plastic toothpick to minimize contamination. Samples were frozen to −20°C after collection.

| Dietary analysis
We undertook DNA extractions in a dedicated pre-PCR laboratory. We extracted DNA from a total of 365 fecal samples using the TATAAACTTC-3′, (Simon et al., 1992), following selection and modification by Stockdale (2018) for amplification of a 306-bp fragment of the COI region (Davies et al., 2022). Stockdale (2018) validated primer sets to ensure DNA amplification of the expected range of taxa. The PCR process involved amplification of the target regions with the addition of a unique combination of 10-bp molecular identifier tags (MID-tags), with samples having a unique pairing of forward and reverse tags for subsequent sample identification. This was undertaken separately for each marker, so all samples could be uniquely identified for both markers. Within each PCR 96-well plate, 12 negative (extraction and PCR) and two positive controls were included following Taberlet et al. (2018). We categorized all products from each individual PCR plate based on band brightness after gel electrophoresis (very faint, faint, medium, bright). We quantified the DNA concentration from a minimum of three representative PCR products per plate from each brightness category using a high sensitivity assay with a Qubit Flourometer (Thermo Fisher Scientific) to confirm whether estimating relative DNA concentration by eye from a gel photo was accurate. For each marker, we pooled each PCR plate according to concentrations determined by the Qubit Fluorometer to ensure equimolar concentration of all PCR products in each pool.
We multiplexed the amplicons into five pools, each containing between 63 and 186 samples. We undertook library preparation for Illumina sequencing using a NEXTflex Rapid DNA-Seq kit (Bioo Scientific), with a unique adapter added to each pool for subsequent bioinformatic identification. We diluted pools to 4 nM and quantified them using Qubit dsDNA High-sensitivity Assay Kits. Finally, we combined the diluted pools equimolarly and sequenced the pools on a MiSeq desktop sequencer via a v2 chip with 2 × 250 bp paired end reads (expected capacity 24-30,000,000 reads). Due to the unbalanced nature of the amplicon libraries, we added a 15% PhiX buffer to the sequencing run to improve crosstalk and phasing calculations.

| Identification of dietary taxa
The Illumina runs generated 6,328,388 and 12,307,560 ITS2 and COI reads, respectively. The bioinformatics pipeline was undertaken following Drake et al. (2021) and Davies et al. (2022). We used FastP (Chen et al., 2018) to trim, align, and quality check sequencing reads. Tagged reads were labeled by unique sample ID using Mothur (Schloss et al., 2009) and subsequently demultiplexed. We implemented Unoise 3 within Usearch (Edgar, 2016) to remove chimeras and cluster reads to generate zero-radius Operational Taxonomic Units (zOTUs) using a 100% clustering threshold, as well as creating a read abundance matrix for zOTUs produced. The Blastn algorithm in Blast+ (Camacho et al., 2009) matched reference sequences contained within GenBank to the sequences generated. We used Megan6 (Huson et al., 2016) to identify unique dietary items using the top hit for each zOTU (based on bit-score). As described in negative controls for each respective zOTU were removed prior to statistical analysis. For example, if a negative control had a read count of 150 sequences for a known dietary zOTU, we removed 150 reads from all samples the zOTU was found in. We opted not to use a standard threshold for retaining zOTUs (such as ten reads per sample, as commonly used within dietary metabarcoding studies). Instead, we determined thresholds for read removal using the percentage at which known artifacts assigned to known positive controls were removed Davies et al., 2022;Drake et al., 2021). To be retained, the read count of a zOTU had to be greater than 3% and 1% of the maximum read count for all detections of that zOTU across all samples for ITS2 and COI, respectively (Davies et al., 2022; Appendices S1-S4). This further cleaning step was performed to account for over-represented taxa tag jumping or "leaking" from samples with extremely high read counts across multiple samples. Data from respective ITS2 and COI libraries were aggregated together to form a single taxon list for each marker, and non-dietary species such as fungi and gastrotrichs were removed.

| Statistical analysis
For all statistical analyses, the presence/absence of each taxonomic unit within a sample was used instead of read count. The use of relative read abundance (RRA) was not deemed suitable, due to the large number of reads detected from the high densities of sunflower seeds provided at all sites, which may skew ecological conclusions generated from RRA data. Furthermore, there are inherent biases present throughout the HTS workflow, including differential DNA extraction success and PCR amplification rates between taxa detected within the diet (Lamb et al., 2019). Additionally, count-based inferences are not advised if little a priori knowledge of the communities analyzed exists (Lamb et al., 2019). Control samples were excluded from the analyses. All statistical analyses were carried out in R version 4.2.2 (R Core Team, 2020) unless otherwise stated. All analysis was undertaken at the genus level to standardize taxonomic resolution, as not all taxa could be identified to species.
To evaluate the most prevalent plant and invertebrate taxa within Hawfinch diet, we calculated frequency of occurrence (FOO). We did this by totalling the number of instances that a given taxon occurred across all hawfinch samples. We then calculated this as a percentage of the total number of samples (% FOO), by dividing the frequency of occurrence by the total number of hawfinch fecal samples collected and multiplying by 100. We also aimed to quantify dietary diversity and reveal if the number of samples collected in this study was sufficient to represent Hawfinch diet within the UK. To achieve this, we used coverage-based rarefaction and extrapolation. We used Hill diversity metrics to estimate species diversity (Roswell et al., 2021;Tercel et al., 2022) as opposed to species-accumulation curves, as these have been shown to poorly represent true community diversity (Roswell et al., 2021). Coverage, as defined by Roswell et al. (2021) "describes to what extent a sample captures the true diversity of the whole community, including species that have not yet been detected". For example, a coverage of 0.80 means that 20% of the individuals in the community being sampled belong to species that have not been found (Roswell et al., 2021). We produced coverage-based rarefaction and extrapolation curves and estimates of species richness (Hill-richness), Shannon's entropy (Hill-Shannon), and Simpson's index (Hill-Simpson) in R using the iNEXT package (Hsieh et al., 2016).
To test whether Hawfinch diet differed between populations and sexes, we applied a single multivariate generalized linear model (MGLM) containing the predictor variables region and sex, using the function manyglm within the package mvabund (Wang et al., 2012).
Due to the small numbers of juveniles (Table 1), they were removed from all analyses. We broadly categorized regions into the following: Wye Valley, north Wales, north Cardiff, New Forest, and East Anglia, as some regions contained multiple catching sites. If an individual was sampled more than once, data was used from the first capture only to avoid pseudo-replication.
Year of sampling was included as a covariate in all models to control for any variation through time.
To assess if any specific plant or invertebrate taxa were responsible for dietary differences regionally or between sexes, p-values from univariate tests were extracted with the function anova.manyglm using the p.uni = "adjusted" argument. This adjusts p-values using Holm's step down resampling algorithm (Westfall & Young, 1993).
We applied parametric bootstrap (Monte Carlo) resampling, ensuring inferences took into account correlation between variables, as this is recommended for presence-absence multivariate data (Wang et al., 2012). When necessary, we performed pairwise comparisons using the pairwise. comp function of anova.manyglm. For all models, we checked quantile-quantile (Q-Q) diagnostic plots for normality, with homoscedasticity checked by plotting Dunn-Smyth residuals against fitted linear predicted values (Bates et al., 2015;Wang et al., 2012). We visualized dietary differences using non-metric multidimensional scaling analysis (NMDS) via the function metaMDS in the vegan package (Oksanen & Wagner, 2019). The nMDS was performed with Jaccard dissimilarities in two-dimensional space (k = 2).
We produced spider plots using nMDS results via ordispider and plotted through ggplot2 (Wickham, 2016). Singletons and outliers were removed to visually improve the ordination. commonly eaten, or evenly consuming dietary taxa, as Hill-richness is most sensitive to rare species (Roswell et al., 2021;Tercel et al., 2022).

| Overview of Hawfinch diet composition
We successfully amplified plant DNA from 262 fecal samples.
We identified 49 plant zOTUs belonging to 17 families, of which 86% could be identified to species and 100% to genus (

| Variation in diet composition between sites and sex
MGLMs revealed that diet differed significantly between sampling regions (MGLM: Wald = 416.5, p = <.001; Figure 4)    Dietary composition differed significantly between sexes (MGLM: Wald = 148.4, p = <.001; Figure 5); however, no univariate tests showed significant differences (Appendix S5). Many taxa were consumed at different frequencies between sexes. Females were found to be feeding more frequently on invertebrates than males (52% of all invertebrate taxa detected were consumed more frequently by females) however, differences for all taxa were smaller than 10% (

| DISCUSS ION
Detailed dietary information is vital for understanding species' ecology and implementing effective conservation management plans (Chua et al., 2021). This study used a metabarcoding approach Hawfinch diet composition not previously recorded using nonmolecular methods. Although this study documents over 100 taxa consumed, it also confirms that previously recorded common food resources such as beech, cherry, and Lepidoptera (Mountford, 1957;Newton, 1967) (Schweiger et al., 2015).
Buds of ash, maple, and beech, as well as Lepidoptera, became important food resources during spring and summer (Mountford, 1957). The importance of beech as a food resource was confirmed in this study, being the most prevalent plant taxon (detected in 38.5% of samples). It is well understood that birds must balance food handling times with net energy intake, and a resource is deemed more profitable if it has a higher energy reward per unit handling time (Molokwu et al., 2011). It is known that Hawfinch  commonly feed on beech nuts during autumn and winter months (Mountford, 1957) due to the high fat and carbohydrate levels compensating for energy losses during winter (Renner et al., 2013). The onset of the breeding season can drive changes in feeding preferences as nutritional needs diversify (Lima, 2009). As the sampling in this study began during the pre-breeding season and continued to the end of summer, Hawfinch may have been gaining a high energy reward from feeding on any remaining available beech nuts, but also obtaining similar nutritional benefits from the increased availability of beech buds in the spring. Beech buds have been shown to contain >15% fat (Lebl et al., 2010), and this may be an important energetic requirement for Hawfinch to boost condition before and during the breeding season.

Percentage of samples testing positive for invertebrate taxa (%FOO)
Our results also revealed that Hawfinch consumed, albeit in low frequencies, herbaceous plant taxa such as dandelions (Taraxacum) and Rubus sp. This would indicate that Hawfinch are consuming these taxa while foraging on the ground or in the shrub layer during the breeding season. Hawfinch were previously thought to be purely canopy feeders during the breeding season (Perea et al., 2014); however, our results suggest that Hawfinch may be able to exploit seasonal food resources which occur in the lower canopy and on the ground.
Hawfinch egg laying begins around mid-April (Kirby et al., 2019), and the presence of invertebrates within the diet during the breeding season is similar to other omnivorous passerine dietary studies conducted over similar temporal periods (Shutt et al., 2020). As invertebrate taxa were found more frequently in female diet than male, it can be presumed that invertebrates may help to provide specific nutrients beneficial to breeding physiology, such as egg production, as well as providing high protein food for chicks (Marshall et al., 2016). The dietary patterns found herein are commonly observed in other passerine species such as Chaffinch (Fringilla coelebs, Linnaeus; Holland et al., 2006). Lepidoptera, Coleoptera, Hemiptera, Annelida, Gastropoda, and Araneae were all observed as prey at the order level (Mountford, 1957), and all (excluding Annelida) were detected within this study. The high prevalence of winter moth within Hawfinch diet is not unexpected, as this larva is an important food resource for other woodland passerine species, such as nestling tits (Perrins, 1991). In contrast, the high prevalence of tree slug within the diet was unexpected, as it was previously thought only snails were consumed (Mountford, 1957). This may be explained by the availability of algae and lichens within woodland, which are the main components of tree slug diet (Kappes, 2006). During wet weather, tree slugs feed on algae growing on tree trunks, but remain under the bark of dead timber during unsuitable weather (Kappes, 2006).
Thus, tree slugs may be taken during periods of high rainfall when foraging efficiency for defoliating Lepidoptera is reduced ( We revealed sexual differences in Hawfinch diet, likely due to behavioral and nutritional differences between males and females.
Behavioral differences between sexes may be a factor due to

TA B L E 3 (Continued)
differences in reproductive roles (da Silva et al., 2020;Freeman, 2014). An important aspect to consider within any DNA metabarcoding study is prey detection biases, which can impact the results and subsequent ecological interpretation of metabarcoding studies . The choice of primers is considered to be one of, if not the most important steps for reducing biases . The choice of primers can impact amplification efficiency and taxonomic classification of subsequent amplicon sequences (Brandon-Mong et al., 2015). As a result, primer choice could influence the understanding of prey composition within the diet and thus the interpretations of foraging ecology (Alberdi et al., 2018;Forsman et al., 2022). Primer biases are considered particularly problematic when undetected taxa are ones which TA B L E 4 Genus level frequency of occurrence and the difference between sexes (Δ).

TA B L E 4 (Continued)
contribute substantially to the foraging ecology of the study species . Hawfinch have previously been observed to feed mainly on Lepidoptera (Mountford, 1957), which was well represented at high frequency of occurrences across sites. The primer pair used in this study were originally used to characterize the diet of blackbirds (Turdus merula) and song thrushes (Turdus philomelos) and subsequently were designed to amplify a broad range of invertebrate taxa (Stockdale, 2018). It is important to note however that no primer pair can provide a completely unbiased and comprehensive account of species' diet due to highly degraded DNA failing to amplify in PCRs, primer biases and differences in mitochondrial copy number per cell (reviewed in Clare, 2014). A one-locus-severalprimer approach should be used more readily within DNA metabarcoding studies in order to maximize taxonomic coverage and minimize false negatives (Corse et al., 2019).
It is also important to acknowledge the possibility of secondary consumption via lepidopteran taxa within the diet, which may result in indirect species associations . Secondary consumption can result in falsely inflated detection of plant taxa through co-amplification of plant DNA within the guts of lepidopteran taxa consumed by Hawfinch. Ecologically, it is known that Hawfinch feed primarily within the canopy (Mountford, 1957;Perea et al., 2014).
This suggests that most invertebrate taxa were obtained from the vegetation or bark within the tree canopy, resulting in possible accidental ingestion of plant taxa when gleaning prey items from trees.
Due to metabarcoding methods being unable to determine which plant tissue is being consumed, in conjunction with Hawfinch also feeding on the same plant taxa as their prey at similar times of the year, differentiating what is "true" secondary predation in this study was extremely challenging.
It is important to note that these conclusions are based upon a small number of samples collected from north Cardiff (n = 7) and East  (Driessen et al., 2013).
In conclusion, this study has provided the first molecular insight into the generalist diet of Hawfinch, at a finer resolution than previous work. We demonstrate that the diet of Hawfinch, as predicted, varies both spatially and between sexes. This dietary variation suggests Hawfinch can respond to changing resource availability by showing dietary plasticity. The results of this study were only possible due to the high taxonomic resolution available through metabarcoding methods. As metabarcoding is becoming more prevalent within ecological research, it becomes increasingly important to understand how taxonomic resolution can impact ecological studies, although species-level identification may not always be necessary, depending on hypotheses studied (Brown et al., 2014;Renaud et al., 2020). Finally, this study demonstrates how the utilization of DNA metabarcoding can increase ecological understanding and improve insights into fine scale ecological patterns. Conceptualization (equal); funding acquisition (equal); investigation (equal); methodology (equal); supervision (lead); writing -original draft (equal); writing -review and editing (equal). Will Kirby:

CO N FLI C T O F I NTE R E S T S TATE M E NT
The authors declare no conflict of interest.

O PE N R E S E A RCH BA D G E S
This article has earned an Open Data badge for making publicly available the digitally-shareable data necessary to reproduce the reported results. The data is available at https://doi.org/10.5061/ dryad.0p2ng f25d.

DATA AVA I L A B I L I T Y S TAT E M E N T
All data supporting the work presented in this study are openly available on Dryad at https://doi.org/10.5061/dryad.0p2ng f25d.