Biodiversity and structure of livestock-associated bacteriophage communities in Brazil and Spain

Livestock waste represents a complex microbial ecosystem in which bacteriophages can play essential roles in shaping bacterial populations and mediating gene exchange. However, comparative analyses of livestock-associated viral communities across geographic regions remain scarce. We conducted a unified metagenomic analysis of 23 Brazilian and 29 Spanish livestock-waste metagenomes, encompassing broiler litter and swine manure samples. Viral operational taxonomic units (vOTUs) were identified and used to compare phage community composition across production matrices, waste-handling stages, geographic regions, and countries. More than 100,000 vOTUs were identified, of which over 98% were classified as bacteriophages. Brazilian samples showed distinct phage-community profiles between broiler and swine production systems, whereas Spanish swine manure samples exhibited modest regional variation at finer taxonomic resolution. Cross-country comparisons indicated broad convergence at the phage class level but substantial compositional differences at finer taxonomic resolution. These findings indicate that livestock-associated phage communities are shaped by ecological context, production matrix, and waste-management practices. This study expands current knowledge of viral biodiversity in agricultural environments and contributes to a One Health understanding of microbial dynamics in livestock-waste systems.

Livestock is a significant contributor to the global economy and food supply, yet intensive farming practices often lead to environmental challenges, including the production of animal waste (Beltrame et al., 2022; Andretta et al., 2021). Among livestock-associated pathogens, bacteria and viruses affect swine production, requiring a comprehensive understanding of their ecological roles, particularly within waste (Ramesh et al., 2021; Kim et al., 2024). Swine manure and broiler litter are examples of these byproducts and constitute complex microbial ecosystems that integrate microorganisms from animals, forming diverse microbiomes composed, for example, of bacteria and viruses. Among these viruses, bacteriophages are abundant particles that replicate within bacterial cells, play key roles in bacterial population dynamics, and may influence horizontal gene transfer (Elois et al., 2023; Sher et al., 2025; Park and Kim, 2020).

Despite their ubiquity and ecological significance, the viral composition of agricultural and livestock waste, such as swine manure, remains poorly characterized compared with other environments (Cook et al., 2021a,b). Few studies have analyzed complex matrices such as broiler litter and swine manure, and these studies highlight the diversity of viral communities, including phages (Ramesh et al., 2021; Cook et al., 2021a,b; Kwok et al., 2020). Phages can also carry auxiliary metabolic genes related to virulence factors, antibiotic resistance, and other functions (Cook et al., 2021a,b). Direct comparison of livestock-associated viromes across countries and production systems remains challenging because biological variation is often related to context and methodological heterogeneity. At the biological level, sample type, production system, antimicrobial use, and farm biosafety practices can structure metagenomic profiles in swine production environments, indicating that differences among regions may reflect management-associated ecology as much as geography itself (Mencía-Ares et al., 2020; He et al., 2021). At the technical level, virome studies are also influenced by how viral fractions are obtained and analyzed, since metagenomes and enriched viromes can yield different ecological signals depending on sample matrix and sequencing strategy (Kosmopoulos et al., 2024). For this reason, cross-country interpretation of livestock viromes requires caution, because apparent geographic contrasts may otherwise combine ecological divergence with differences in sampling context, production practices, and analytical approach. Geographical variations in production systems and antimicrobial usage have been shown to significantly influence the resistome and mobilome profiles of swine production environments, suggesting that regional factors may similarly structure associated phage communities (Mencía-Ares et al., 2020).

Within this perspective, the One Health concept, which recognizes animal, environmental, and human health as interconnected domains, provides a relevant conceptual basis for investigating viral communities in agricultural environments. Despite the recent advances in viral genome catalogs, the reference space available for viral metagenomics still requires improvement, and this limitation is particularly relevant for livestock-associated systems. Metagenomic sequence discovery continues to outpace the number of experimentally characterized viral isolates, making sequence-based inference essential but often without adequate taxonomy characterization (Simmonds et al., 2023). Although resources such as IMG/VR v4 now contain more than 15 million viral genomes and genome fragments, they still provide only a partial view of the unresolved virosphere (Camargo et al., 2023). These findings indicate that livestock viromes remain only partially represented in current databases, which helps explain the limited taxonomic depth typically recovered in comparative studies of agricultural environments.

In this context, the present study performs a comparative metagenomic analysis of bacteriophage communities based on 23 Brazilian metagenomes, including paired swine biodigester influent and effluent samples as well as broiler litter, thereby extending the investigation across distinct livestock production systems. These data were compared with previously available Spanish swine manure metagenomes, retaining only raw manure samples to ensure conceptual comparability with Brazilian swine influent prior to biodigestion. Brazil and Spain were selected as representative livestock-production settings with relevance to swine farming. Given the relevance of both countries in livestock production, this study examines how viral communities are structured in animal production systems within each national context and how they differ across geographic livestock settings contributing to the advancement of knowledge in the context of One Health in agri-food environments.

Animal metagenomes

Brazilian samples comprised eight poultry litter metagenomes collected in 2019 from broiler production farms in the western region of Santa Catarina State, Brazil. These included seven samples from conventional broiler flocks and one sample from a farm without antimicrobial use. Swine farm samples (n = 15) were obtained from the biodigester sampling framework described by Fongaro et al. (2014) and comprised seven paired biodigester influent and effluent samples from conventional swine farms collected annually between 2013 and 2019, plus one additional biodigester influent sample collected in 2019 from a swine farm without antimicrobial use (BioProject PRJEB63971). Spanish swine manure metagenomes were obtained from an available dataset (Mencía-Ares et al., 2020). Raw FASTQ files were downloaded from the NCBI Sequence Read Archive (BioProject PRJNA628671), which includes 105 metagenomes from 38 farms. Only manure samples were retained (n = 29), distributed across 10 regions (Castilla León 7; Andalucía 6; Murcia 4; Galicia 3; Castilla-La Mancha 2; Extremadura 2; Valencian Community 2; Aragón 1; Asturias 1; Catalonia 1). The full per-sample metadata and group sizes are provided in Supplementary Table S1b. Supplementary Table S1 is organized into Supplementary Table S1a (biodigester influent-effluent taxon turnover), Supplementary Table S1b (per-sample sequencing statistics), Supplementary Table S1c (phage vOTU taxonomy), Supplementary Table S1d (all-vOTU taxonomy before phage-only filtering), Supplementary Table S1e (per-sample phage abundance matrices), Supplementary Table S1f (non-phage vOTUs removed by filtering), and Supplementary Table S1g (complete statistical test outputs, raw p and BH-FDR q values for every test). In this study, downloaded reads were reprocessed using the same quality-control, decontamination, and assembly workflow applied to the Brazilian dataset to ensure analytical comparability across countries.

Bioinformatic analyses

Quality control and preprocessing

To minimize batch effects, Brazilian and Spanish raw reads were processed under a harmonized pipeline, and all downstream analyses were conducted using the same filtering and taxonomic criteria. Briefly, sequencing adapters and low-quality bases were removed using Cutadapt (Martin, 2011) (v2.8) and Trim Galore (v0.6.5-1), respectively. To eliminate non-microbial background sequences associated with animal hosts, feed, and bedding materials, reads were mapped against reference genomes using Bowtie2 (Langmead et al., 2019) (v2.4.1) with default parameters. The reference genomes used for host and environmental contamination removal included Gallus gallus domesticus (GenBank accession GCA_016699485.1), Sus scrofa (RefSeq accession GCF_000003025.6), and plant genomes commonly associated with animal production systems, including Eucalyptus grandis (GCF_016545825.1), Zea mays (GCF_902167145.1), Glycine max (GCF_000004515.6), Triticum aestivum (GCF_018294505.1), and Sorghum bicolor (GCF_000003195.3). All reference genomes were downloaded in November 2025.

vOTU abundance

Cleaned paired-end reads were mapped against multiple viral reference databases [IMG/VR high-confidence genomes (Simmonds et al., 2023), INPHARED (Cook et al., 2021b), RVDB (Chin et al., 2025), and RefSeq Viral (O'Leary et al., 2016) representative genomes] using minimap2 (short-read preset, -ax sr). Alignments with mapping quality ≥10 were retained, while unmapped, secondary, and supplementary alignments were excluded (samtools -F 2308; minimap2 –secondary=no). Viral abundance was quantified as mapped read counts per reference sequence using samtools idxstats (Danecek et al., 2021).

The deduplication of identical reference sequences across databases was applied using exact-sequence MD5 hashing (after sequence normalization), generating canonical reference identifiers. Only canonical references observed in at least one sample (≥3 mapped read) were retained for downstream analyses. vOTUs were then defined by clustering the observed canonical sequences at ≥95% nucleotide identity and ≥85% coverage using MMseqs2 (Mirdita et al., 2021).

Phage-specific taxonomic retention in vOTUs

The vOTU taxonomy was obtained using a conservative hierarchical pipeline. Taxonomy was assigned sequentially to IMG/VR, INPHARED, BLASTn fallback, and an optional gap-filling step based on NCBI Taxonomy. For each vOTU, member canonical sequences were reconstructed from the canonical-to-vOTU mapping, and a representative sequence was defined for downstream similarity-based analyses. Primary taxonomic annotation was obtained by linking vOTU members to IMG/VR identifiers (UVIGs) through a curated reference mapping and extracting viral taxonomic lineages from the high-confidence IMG/VR metadata. When multiple IMG/VR references were associated with a given member, the annotation with the greatest taxonomic depth was retained.

vOTU-level taxonomy was then inferred using a conservative rank consensus approach, whereby taxonomic ranks were assigned only when supported by at least two concordant annotations, and assignment proceeded hierarchically from higher to lower ranks. vOTUs remaining unclassified after IMG/VR annotation were queried against the INPHARED database using Mash v2.3 (Huang et al., 2024), restricting analyses to representative sequences to reduce redundancy; only best matches meeting conservative distance thresholds were retained, and taxonomy was transferred from curated INPHARED metadata using the same consensus criteria. vOTUs still unresolved were subjected to a BLASTn-based fallback against IMG/VR high-confidence genomes, applying filters of 95% nucleotide identity and 85% query coverage and assigning taxonomy through a conservative lowest-common-ancestor consensus of top-scoring hits, with explicit rank caps to prevent overassignment. Finally, an additional post-processing step queried NCBI Taxonomy to fill higher-rank gaps only when lower-rank assignments were resolved. Using taxonomic classification, sequences assigned to bacteriophage-related classes were identified and curated separately from the overall viral dataset, allowing the extraction of phage-specific vOTUs. This approach was applied consistently to both the Brazil and Spain datasets.

Diversity metrics and statistical analyses

Sample-by-vOTU count matrices were filtered by dataset-specific minimum phage-assigned library size prior to ecological analyses (Brazil phage-only ≥1,000 reads; Spain phage-only ≥50 reads). These threshold were set as a permissive lower bound for sample inclusion; in practice all retained libraries greatly exceeded them (minimum 933,324 phage-assigned reads in Brazil and 283,306 in Spain), so no samples were excluded on this criterion. Feature prevalence filtering was applied only to beta-diversity matrices, while alpha diversity was calculated from unfiltered count tables after sample filtering, to avoid biasing within-sample diversity estimates by excluding rare taxa. Alpha diversity was estimated using Shannon and Gini–Simpson indices. Rarefaction depth was set to 80% of the minimum library size among retained samples (capped at 200,000 reads), with 50 iterations to obtain mean estimates and 95% confidence intervals. Beta diversity was assessed using complementary dissimilarity frameworks: Jaccard (presence/absence), Bray–Curtis (TSS-normalized abundances), Hellinger-Euclidean, and Aitchison (CLR-transformed). Ordinations were generated using principal coordinate analysis (PCoA). Group differences were tested using PERMANOVA with 9,999 permutations, and homogeneity of dispersion was assessed using PERMDISP.

Groups with fewer than two samples were excluded from inferential testing, and results for singleton groups were reported descriptively. P -values were adjusted using the Benjamini–Hochberg false discovery rate (BH-FDR). Among the Brazilian samples, those with antibiotic use will be referred to as W/A, and those without antibiotic use will be referred to as No/A. The sample-type PERMANOVA (broiler litter vs. biodigester influent vs. effluent) was restricted to conventionally managed (W/A) samples, excluding two organically managed (No/A) samples, one broiler litter and one biodigester-influent, reported descriptively not included in inferential testing. The resulting comparison included n = 21 samples (broiler litter = 7, biodigester influent = 7, biodigester effluent = 7). The host contrast (broiler litter vs. biodigester effluent) and the influent–effluent contrast used the same W/A-only subsets (n = 14 each, 7 per group). The Spanish regional test (Results 3.6) was restricted to regions with ≥2 samples (seven regions; n = 26 of 29) and conducted at the phage vOTU level. The pooled cross-country test (Section 3.7) compared Brazil (n = 15, biodigester influent NoA + influent + effluent) against Spain (n = 29) at the phage class level. All PERMANOVA tests used Jaccard, Bray-Curtis, Hellinger-Euclidean, and Aitchison CLR-Euclidean distances with 9,999 permutations. Where appropriate, additional permutation-based approaches were applied to evaluate within vs. between-group distances and paired biodigester influent–effluent differences. For the Spain total-virome regional analysis, alpha diversity differences were additionally evaluated using non-parametric tests (Kruskal–Wallis and pairwise Mann–Whitney U), with BH-FDR correction.

To further visualize dissimilarity patterns between countries, we additionally calculated the median class level Bray–Curtis dissimilarity between each Spanish region and three Brazilian reference groups, biodigester influent, biodigester effluent, and the single organic swine raw manure sample. This comparison is descriptive since the organic raw manure group contains only one sample and cannot support statistical testing, consistent with the treatment of singleton groups described above. The goal was to show how consistently each Spanish region relates to different types of Brazilian swine samples, complementing the pooled two group test reported for the cross country comparison.

Taxonomic differential abundance and taxon-shift analyses

Phage taxonomic groups were analyzed using compositional approaches. Between-group differences were tested on CLR-transformed counts using two-sided Mann–Whitney U tests. Effect size was reported as median CLR difference (ΔCLR) and rank-biserial correlation. P -values were adjusted using BH-FDR (α = 0.05). Taxa were retained for visualization when they met significance thresholds (BH-FDR <0.05) and showed large effect sizes (|ΔCLR| ≥ 2). Paired biodigester turnover (influent vs. effluent) was evaluated by calculating per-pair log 2 fold-change in TSS-normalized abundances and testing median shifts using Wilcoxon signed-rank tests (minimum three valid pairs), followed by BH-FDR correction.

vOTU recovery and taxonomic overview

Across matrix-wide profiling after abundance/prevalence filtering (reads ≥ 3; prevalence ≥ 2), Brazil included 64,528 detected vOTUs (n = 23 samples), whereas Spain included 36,978 detected vOTUs (n = 29 samples). After phage-focused retention, Brazil retained 63,348 phage vOTUs (98.17% of the total detected set), whereas Spain retained 36,772 phage vOTUs (99.44%).

Phage biodiversity in Brazil dataset

Taxonomic composition of resolved phage vOTUs was assessed at class, family, and genus levels across three Brazilian sample groups (Figures 1A–C). Across all groups, 96.5% of phage vOTUs (61,109/63,348) were resolved to class level, while family and genus-level resolution remained limited (3.2% and 3.6%, respectively), reflecting the novelty of the recovered phage diversity. At the class level (Figure 1A), Caudoviricetes dominated across all three groups, comprising 99.1% of class-resolved vOTUs in Broiler litter, 98.8% in Biodigester influent, and 99.0% in Biodigester effluent samples. Minor contributions were observed from Faserviricetes (~0.5%−0.6%), Malgrandaviricetes (~0.2%−0.5%), and Leviviricetes (~0.1%−0.2%) across all groups.

At the family level (Figure 1B), distinct compositional profiles emerged among groups. Broiler litter was enriched in Autographiviridae (16.6%), Herelleviridae (14.1%), and Straboviridae (13.1%), alongside a contribution from Inoviridae (12.7%). Biodigester influent samples showed a more even distribution, with Microviridae (15.6%), Straboviridae (14.9%), Autographiviridae (14.8%), and Inoviridae (14.7%) as the dominant families. Biodigester effluent samples were characterized by Straboviridae (18.0%), Autographiviridae (14.9%), and Inoviridae (14.6%), with a moderate contribution of Microviridae (10.7%). Overall, Peduoviridae was exclusively identified in Broiler (7.6%), while Herelleviridae was higher in broiler (14.1%) than in biodigester effluent (6.6%). At the genus level (Figure 1C), Felixounavirus and Tequatrovirus were the most abundant genera in Broiler litter (5.1% and 4.3%, respectively), followed by Lederbergvirus (3.4%), Kayvirus (2.6%), and Punavirus (2.4%). In Biodigester influent, Tequatrovirus was the dominant genus (7.2%), followed by Felixounavirus (4.9%) and Microvirus (4.7%). Biodigester effluent samples were similarly dominated by Tequatrovirus (8.3%), with lower contributions from Microvirus (3.0%) and Felixounavirus (2.8%). The genera collectively accounting for the remaining ~50% of resolved genus-level vOTUs are described in Supplementary Figure S1. This low-resolution came from two factors: the genomic novelty of livestock-waste phages and gaps in reference databases below the class level. Family and genus level results should therefore be read as exploratory hypotheses, not definitive conclusions.

Phage diversity in Spain dataset

Across Spanish regions, 98.2% of phage vOTUs (36,124/36,772) were resolved to class level, while family and genus-level assignments were recovered for only 1.6% and 1.6% of vOTUs, respectively. At the class level (Figure 2A), Caudoviricetes was the dominant class in all ten Spanish regions, ranging from 98.8% (Asturias) to 99.6% (Aragon, Castilla-La Mancha, Catalonia, Extremadura, Murcia, Valencian Community). Faserviricetes, Malgrandaviricetes, and Leviviricetes were consistently present as minor components (<1%) across all regions, with no region showing marked enrichment of any minor class. At the family level (Figure 2B), regional variation was more pronounced. Inoviridae was dominant in Andalucia (22.1%) and Asturias (27.4%), while Herelleviridae was the primary family in Aragon (19.9%), Castilla-La Mancha (17.9%), and Catalonia (12.6%). Schitoviridae was consistently prominent across Catalonia (16.5%), Extremadura (17.4%), Galicia (15.3%), and Murcia (14.6%). Autographiviridae and Microviridae were distributed across all regions, each contributing 8%−15% depending on the region. Microviridae was abundant in the Valencian Community (19.6%) and Murcia (15.4%). At the genus level (Figure 2C), Punavirus emerged as the most consistently recovered genus across all Spanish regions, ranging from 1.8% (Murcia) to 5.8% (Asturias). Region-specific patterns were also observed, where Felixounavirus was the leading genus in Extremadura (4.1%) and Castilla-La Mancha (3.1%); Moineauvirus and Skunavirus co-dominated in Murcia (both at 2.2%); and Samwavirus was most prominent in Asturias (3.2%), Catalonia (2.8%), and Galicia (1.9%). Tequatrovirus, a dominant genus in the Brazilian Biodigester groups, was present in Spain only at low levels in Castilla León (1.7%) and Valencian Community (2.3%). As with Brazil, genera collectively representing the remaining ~50% of genus-level resolved vOTUs for Spanish phage dataset are detailed in Supplementary Figure S2. As in the Brazilian dataset, this reflects the genomic novelty of environmental phages, family and genus-level results should likewise be treated as exploratory.

Alpha diversity and beta diversity

Alpha and beta-diversity patterns for the Brazilian phage dataset are summarized in Figures 3A–F. At the phage vOTU level, alpha diversity differed among the three sample types for Observed richness (Kruskal–Wallis: H = 7.24, p = 0.027) and Gini-Simpson (H = 6.88, p = 0.032), but not for Shannon diversity (H = 4.47, p = 0.107). Median Shannon values were 5.213 (broiler litter), 5.677 (biodigester influent), and 4.144 (biodigester effluent); median Gini–Simpson values were 0.961, 0.956, and 0.908, respectively. After Bonferroni correction (three pairwise contrasts), only the broiler litter vs. biodigester effluent contrast remained significant for Gini–Simpson (Mann–Whitney U = 46, q = 0.012), whereas no other pairwise contrast reached significance. Paired within-digester comparisons were likewise non-significant for both Shannon and Gini-Simpson (Wilcoxon V = 7, p = 0.297).

Beta diversity of taxonomically aggregated phage class profiles, assessed using four distance metrics (Jaccard, Bray–Curtis, Hellinger-Euclidean, and Aitchison CLR-Euclidean), revealed robust compositional differences among sample types (broiler litter, biodigester influent, and biodigester effluent; n = 21, 7 per group), as indicated by significant PERMANOVA results for all metrics (Jaccard F = 4.168, R 2 = 0.317, p = 0.0294, q = 0.029; Bray–Curtis F = 4.960, R 2 = 0.355, p = 0.0041, q = 0.006; Hellinger F = 7.774, R 2 = 0.463, p = 0.0001, q <0.001; Aitchison F = 4.993, R 2 = 0.357, p = 0.0045, q = 0.006). PERMDISP tests were non-significant across all metrics (q > 0.1), indicating that these differences reflect true shifts in community composition, not differences in within-group dispersion.

The host contrast was evaluated separately for broiler litter vs. biodigester influent and broiler litter vs. biodigester effluent. Broiler litter differed significantly from biodigester influent for Bray–Curtis (q = 0.028), Hellinger (q <0.001), and Aitchison (q = 0.010), but not for Jaccard (q = 0.155). Broiler litter differed significantly from biodigester effluent across all four metrics (Jaccard q = 0.016, Bray–Curtis q = 0.005, Hellinger q <0.001, Aitchison q = 0.0014). PERMDISP was non-significant for all metrics in both comparisons, with q above 0.1 throughout, indicating that these differences reflect genuine compositional separation. Together, these two comparisons show a consistent host effect. Broiler litter differs from the swine biodigester system whether compared to the influent or the effluent stage, reinforcing that the observed separation reflects the host system rather than the digestion stage. By contrast, no significant compositional difference was detected between biodigester influent and effluent by PERMANOVA in any metric (Jaccard q = 0.506, Bray–Curtis q = 0.313, Hellinger q = 0.316, Aitchison q = 0.680), and PERMDISP was likewise non-significant for all metrics. Overall, phage community structure, interpreted through taxonomically annotated phage vOTU profiles, was associated with the host system (broiler vs. biodigester), with no detectable effect of biodigestion on either community composition or within-group dispersion.

Taxon shift after biodigestion

Paired influent–effluent comparisons (n = 7 digesters) showed no statistically significant taxa shifts after multiple-testing correction at class level (all BH-FDR q > 0.05), and sign-flip sensitivity tests were likewise non-significant. Overall, the results do not support statistically significant class-level turnover during biodigestion. Exploratory results are demonstrated in Supplementary Table S1a.

Phage vOTUs from manure in Spain

Alpha diversity of bacteriophage communities at the vOTU level was comparable across Spanish regions, with median Shannon values of 4.535 (Andalucia), 5.053 (Castilla León), 5.629 (Murcia), 5.632 (Extremadura), 5.274 (Castilla-La Mancha), 5.616 (Valencian Community), 5.256 (Catalonia), 5.668 (Galicia), 3.621 (Asturias), and 5.649 (Aragon); median Gini-Simpson ranged from 0.872 (Asturias) to 0.957 (Aragon; Figures 4A, B). Consistent with these patterns, Kruskal–Wallis testing among replicated regions (n ≥ 2) did not detect significant differences in Shannon or Gini–Simpson diversity.

Beta diversity at the phage vOTU level was evaluated using Jaccard, Bray–Curtis, Hellinger-Euclidean, and Aitchison CLR-Euclidean distance matrices, respectively shown in Figures 4C–F. PERMANOVA detected significant variation among the better-represented regions for Bray-Curtis (F = 1.856, R 2 = 0.369, q = 0.003), Hellinger (F = 1.609, R 2 = 0.337, q = 0.003), and Aitchison (F = 1.302, R 2 = 0.291, q = 0.031) after FDR correction, whereas the Jaccard presence-absence metric was not significant (F = 1.191, R 2 = 0.273, q = 0.076). PERMDISP was non-significant across all metrics (Jaccard: F = 5.412, p = 0.094; Bray–Curtis: F = 2.635, p = 0.476; Hellinger: F = 2.834, p = 0.425; Aitchison: F = 1.723, p = 0.664), supporting interpretation of observed differences as centroid shifts. Together, these results indicate modest but significant regional structuring in the abundance structure of phage vOTUs across Spanish manure samples, without corresponding shifts in coarse class-level composition.

Comparison of phage vOTU profiles between Spain and Brazil

Cross-country class-level ordination and similarity patterns are summarized in Figure 5A. The Spanish manure samples represent raw, non-digested swine manure and are conceptually comparable to Brazilian Biodigester influent samples. At the phage class level, both countries were overwhelmingly dominated by Caudoviricetes (~99% of resolved phage vOTUs), indicating broad conservation of coarse taxonomic structure. Even so, the pooled Spain-vs.-Brazil swine comparison was significant by PERMANOVA (F = 5.56, R 2 = 0.117, p = 0.014). However, PERMDISP was also significant (F = 15.5, p <0.001), indicating that the observed separation reflects both centroid differences and unequal within-group dispersion and should therefore be interpreted with caution. By contrast, the comparison between Brazilian biodigester influent and biodigester effluent remained non-significant (PERMANOVA: F = 1.396, R 2 = 0.104, p = 0.313; PERMDISP p = 0.278). Figure 5B shows that each Spanish region was about equally dissimilar to all three Brazilian groups, biodigester influent, effluent, and biodigester influent without antibiotic use, this third group (n = 1) was included descriptively only, as a single-sample group cannot support formal beta-diversity testing.

Despite their ecological importance, structured cross-country comparisons of livestock viromes remain scarce. In this study, we analyzed livestock-associated metagenomes from Brazil and Spain, including swine biodigester influent and effluent compartments, to enable a structured comparison of phage communities based on vOTUs. By reprocessing both Brazilian and Spanish datasets under a unified analytical framework, we aimed to improve comparability and address limitations of previous studies, such as Ramesh et al. (2021), which focused primarily on pathogen discovery. The pattern observed was consistent with the high bacterial biomass sustained by manure matrices, which supports abundant bacteriophage populations (Sharma et al., 2023; Huang et al., 2024).

It is important to interpret the present dataset as a conservative representation of the recoverable and classifiable viral fraction, rather than a complete inventory of the environmental virome. This limitation arises from methodological constraints such as nucleic-acid recovery efficiency, assembly performance, abundance filtering, and the incompleteness of environmental virus databases.

Within Brazil, alpha diversity of taxonomically annotated phage profiles differed among sample groups for Observed richness and Gini-Simpson, whereas Shannon diversity did not (as reported in Results). These differences among diversity indices reflect their distinct sensitivities to taxon abundance. Observed richness captures rare vOTUs, the Gini-Simpson index is driven by dominant lineages, and the Shannon index gives intermediate weight across the full abundance range. The significant differences in richness and Gini–Simpson, but not in Shannon, suggest that group differences occur mainly at both ends of the abundance distribution, among rare taxa and among dominant phage lineages, while Shannon, by integrating across all abundances, does not reach significance. Pairwise comparisons between biodigester influent and effluent were not significant in either unpaired or paired within-digester tests, indicating that the overall Kruskal–Wallis signal is driven by the host-system contrast and not by the biodigestion process itself. Notably, broiler litter showed higher diversity than swine-associated samples, consistent with a more even distribution of phage lineages across groups. Together, these observations support the idea that phage community profiles in livestock waste are shaped primarily by production matrix and host context.

Across all Brazilian sample types, Caudoviricetes accounted for approximately 99% of classified phage vOTUs, consistent with the known dominance of tailed dsDNA phages in anaerobic livestock waste systems (Rossi et al., 2022; Zhang et al., 2017). Because only ~3% of Brazilian phage vOTUs could be resolved to family level (and even fewer to genus), the family and genus taxonomic the findings described below should be regarded as exploratory and interpreted as hypotheses. At the family level, broiler litter was enriched in Autographiviridae, Herelleviridae, Straboviridae, and Inoviridae, whereas biodigester samples showed prominent contributions from Straboviridae, Autographiviridae, Inoviridae, and Microviridae, a pattern partially concordant with pig gut phageome surveys (Hu et al., 2024; Chen et al., 2025). Autographiviridae, which encompasses strictly virulent T7-like podophages targeting Gram-negative hosts such as E. coli, were among the top families in broiler litter and swine influent (Boeckman et al., 2022; Walker et al., 2020). Their prominence is consistent with the presence of Gram-negative enteric hosts, including Proteobacteria-associated taxa, in poultry litter and swine waste environments (Rossi et al., 2022). Herelleviridae, which comprises large-genome contractile myophages infecting Firmicutes including Bacillus, Staphylococcus, Lactobacillus, and Enterococcus, was recovered across all Brazilian groups, reflecting the Firmicutes-dominated nature of both broiler gut and swine manure microbiomes. Straboviridae (T4-like myophages) and Schitoviridae (N4-like podophages) also contributed across groups, both targeting Gram-negative bacteria and confirmed as core members of pig gut phage communities (Hu et al., 2024).

Within the subset of vOTUs resolved to family level, Microviridae accounted for 15.6% of family-resolved vOTUs in biodigester influent and 10.7% in biodigester effluent. Microviridae did not rank among the dominant families in broiler litter. This pattern suggests that the higher representation of Microviridae in swine-associated samples reflects specifically host context, instead of effect of anaerobic digestion, consistent with the high E. coli density of concentrated swine waste and the reported abundance of Microviridae in pig fecal viromes (Rossi et al., 2022; Wang et al., 2019), though this remains exploratory given the limited family-level resolution.

Spanish swine manure samples were also dominated by Caudoviricetes across all 10 regions, with mean sample-level proportions ranging from 99.61% (Andalucia) to 99.94% (Murcia and Aragon). Faserviricetes remained a consistently minor component (<0.5%) in every region. Thus, both Brazilian and Spanish phage datasets converge on a strongly Caudoviricetes -dominated class profile, and the Brazil–Spain contrast does not involve replacement of the dominant phage class. Ecological differences between these systems appear to be represented at finer taxonomic resolution or in vOTU-level abundance structure.

Because three of these regions, Aragón, Asturias, and Catalonia, are represented by a single farm each, and several others by only two or three, these regional patterns at the family level describe individual farms or very small groups. Given this, and that only a small fraction of Spanish vOTUs could be assigned to family or genus, these results should nevertheless be regarded as exploratory. Because Inoviridae establish chronic productive infections in Gram-negative bacteria, whereas Herelleviridae infect Firmicutes and Schitoviridae represent N4-like phages frequently associated with Gram-negative hosts, these regional family-level differences may reflect variation in the underlying bacterial communities and farm-management context (He et al., 2021; Sukhum et al., 2021; Kumari et al., 2015). Genus-level profiles further support this interpretation. Tequatrovirus was the dominant resolved genus in Brazilian biodigester influent and effluent but was detected only at low relative abundance in Spain, whereas Punavirus was the most consistently recovered genus across Spanish regions. These contrasting genus-level profiles indicate that ecologically meaningful differences between the two countries emerge more clearly at finer taxonomic resolution than at the class level, where both datasets remain dominated by Caudoviricetes.

Within Brazil, taxonomically aggregated phage class profiles showed significant community structuring among broiler litter, biodigester influent, and biodigester effluent across all four beta-diversity metrics, including Jaccard, Bray-Curtis, Hellinger, and Aitchison, with no evidence of heterogeneous dispersion by PERMDISP at the global level. These results indicate that differences among Brazilian sample groups are detectable both in coarse taxonomic composition and in relative abundance structure across production matrices, although the pairwise broiler-vs.-swine contrast showed significant dispersion heterogeneity for Aitchison with PERMDISP non-significant for all metrics, indicating that these differences reflect genuine compositional shifts, not heterogeneous within-group dispersion. The Brazil–Spain contrast reflects differences in waste-processing systems, with Brazilian swine samples originating from anaerobic biodigesters, whereas Spanish samples represented raw manure stored in outdoor slurry pits. These systems differ in oxygen exposure, retention time, pH, temperature, and microbial host composition, all of which are expected to influence viral community structure (Rossi et al., 2022; Fan et al., 2023; Shi et al., 2024). However, both countries remained dominated by Caudoviricetes at the phage class level, indicating that the contrast between systems does not involve replacement of the dominant phage class. In the pooled class-level comparison, Spain and Brazil differed significantly by PERMANOVA, but PERMDISP was also significant, indicating that this contrast should be interpreted cautiously because it reflects both centroid separation and unequal dispersion. Because country and waste-processing type are linked in this dataset, it is not possible to separate their individual effects. The observed differences may reflect manure management, but this cannot be tested independently here.

We did not detect a difference between biodigester influent and effluent in alpha diversity, in class level composition, or in paired taxon turnover. This suggests that the phage community entering the digester is similar to the community leaving it. A likely explanation is that a stable bacteriophage community lives inside the reactor, kept in place by the high bacterial biomass and constant host availability, and is not reorganized by the digestion process. Anaerobic digestion can create environmental stress that triggers prophage induction (Rossi et al., 2022), and studies of other full scale biogas plants have found real differences in viral diversity between the waste going in and the digestate coming out (Shi et al., 2024), so a lack of change is not guaranteed across all systems. Swine manure digesters specifically can also host active, inducible temperate phages that respond to changing reactor conditions, as shown in a recent study on mineral amendments and ammonia stress (Peng et al., 2025). It is therefore possible that biodigestion shifts phage abundances or triggers specific prophages without replacing the dominant Caudoviricetes community. We cannot currently tell true community stability apart from changes that our method is not able to detect. This gap is mainly a limitation of viral reference databases and of taxonomy-based classification itself. As viral databases become more complete over time, family and genus level classification should improve, and host range prediction methods, which do not depend on taxonomy, could help resolve this question sooner.

Among Spanish regions, alpha diversity at the phage vOTU level was comparable, with no significant Kruskal–Wallis signal for Shannon or Gini–Simpson across replicated regions. In contrast, vOTU-level beta diversity revealed significant regional structuring for abundance-sensitive metrics (Bray–Curtis, Hellinger, Aitchison; q ≤ 0.031), whereas the presence-absence Jaccard metric was not significant. By contrast, no robust regional structuring was detected at the phage class level across any of the four metrics. Together, these results indicate that regional differences in Spanish slurry are detectable in the abundance structure of phage vOTUs, but do not translate into major shifts in general class-level composition. Extensive Iberian systems were characterized by outdoor farming, lower animal density, later slaughter age, and lower antimicrobial use, whereas intensive farms showed distinct slurry-associated structuring. Notably, nine extensive farms lacked slurry pits, indicating that the slurry subset already represented heterogeneous manure-management conditions.

Because the Spanish farms spanned heterogeneous management systems (intensive and extensive, differing in antimicrobial use and housing), the weak class-level regional structuring is unlikely to be simply a consequence of uniform management. Instead, coarse composition appears broadly conserved across these varied conditions, with differences emerging mainly in vOTU-level abundance. These patterns should nonetheless be interpreted as exploratory. Our interpretation is supported by both the present data and the literature. He et al. (2021) and collaborators showed that virome composition was more associated with farm management and biosafety practices than with geographic distance within a country (He et al., 2021), consistent with the possibility that waste-processing regime is an important source of variation in these systems. The biodigester and slurry contrast may represent an ecological gradient in livestock waste systems associated with the large physicochemical differences between these environments. Taken together, these results suggest that matrix and production context are associated with phage community profiles, potentially alongside geography. They also highlight that communities can maintain comparable local richness and evenness while detectably reorganizing their taxonomic composition. These patterns suggest that the mode of ecological differentiation varies with the level of taxonomic resolution examined. Within Brazil, PERMANOVA significance extended to the presence-absence Jaccard metric, indicating shifts in phage class membership across production matrices.

In contrast, for Spanish vOTU-level regional comparisons, Jaccard was not significant while abundance-weighted metrics were, indicating abundance reallocation within the shared vOTU pool. Together, these observations suggest that environmental filtering operates across multiple layers: production matrix and host context restructure coarse phage profiles, whereas regional heterogeneity within a single management system modulates the proportional representation of shared members without replacing them. This is consistent with metacommunity theory (Leibold et al., 2004), which predicts that regional species pools are often shared among local sites, with local assembly rules determining realized composition.

In livestock waste systems, the available regional phage pool may be similar across farms, while local ecological conditions, substrate composition, bacterial host community, and also waste-processing regime determine which members proliferate and at what abundance. Future studies combining host-range prediction with longitudinal sampling across biodigester compartments could test whether the compositional stability observed here reflects true ecological equilibrium or a transitional state influenced by retention time and organic loading. The taxonomic depth achieved here also underscores a current limitation in the field, even with the best available reference resources, many bacteriophage sequences remain unclassifiable at fine resolution, and expanding curated viral databases will be essential for accurate ecological interpretation of livestock-associated viromes.

By comparing the biodiversity of viruses in livestock metagenomes from Brazil and Spain under a unified workflow, this study shows that the taxonomically resolved read-based viral fraction was strongly enriched in bacteriophage-associated sequences, with high taxonomic recovery at least at one rank (>98%) and high phage-retained fractions in both datasets (98.17% in Brazil and 99.44% in Spain).

In Brazil, community structure was associated with the host production matrix by the host production matrix, with clear separation of broiler litter samples and no significant influent-effluent shifts in swine biodigesters. In Spain, alpha diversity was comparable across regions, whereas vOTU-level beta diversity revealed modest regional structuring without major shifts in coarse class-level composition. Cross-country comparisons suggested broad convergence at the phage class level, whereas remaining between-country differences are more likely concentrated at finer taxonomic resolution. Collectively, these findings support a model in which livestock-associated virome composition is associated with production matrix, waste-handling regime, and local management context, including differences among intensive and extensive production systems, while biogeographic patterns contribute an additional layer of between-country differentiation.