Biofilm thickness controls the relative importance of stochastic and deterministic processes in microbial community assembly in moving bed biofilm reactors

The study addresses how biofilm thickness influences the balance between deterministic (niche-based selection) and stochastic (drift, dispersal, diversification) processes during microbial community assembly in engineered systems. Prior research shows both process types can be important, but their relative contributions vary across systems and conditions. Biofilms are central to many engineered and natural environments and can enhance retention of slow growers and system stability. However, prior results on assembly in biofilms are mixed, with stream biofilms often showing deterministic assembly while microbial electrolysis cell biofilms often show stochastic assembly. The authors hypothesized that both processes contribute to biofilm assembly, that selection operates across all thicknesses but differs by thickness, and that thin biofilms (more homogeneous environment and competition for space) would experience stronger homogeneous selection and potential competitive exclusion, while thicker biofilms (greater diffusion limitation and microenvironmental heterogeneity) would exhibit more variable selection and/or stochasticity. They tested these hypotheses by comparing established biofilms across five controlled thicknesses (50–500 µm) using neutral modeling and null-model phylogenetic analyses.

- Engineered microbial environments frequently exhibit combined stochastic and deterministic assembly, but relative contributions vary with selection strength and dispersal rates. - Studies report divergent conclusions: mature stream biofilms tend to assemble deterministically, whereas biofilms in microbial electrolysis cells often assemble stochastically. Confounding selective regimes in many studies obscure the specific contribution of biofilm formation per se. - Biofilm thickness is an emergent property influenced by hydrodynamics, loading, and community composition; it affects gradients and microhabitats and thus selection pressures. Despite extensive modeling of thickness effects on treatment performance, experiments that control thickness are fewer, though such control can be leveraged to manage communities and functions. - Prior work with two thicknesses pointed to deterministic differences but did not explicitly assess stochasticity. The authors build on this by explicitly incorporating neutral modeling, considering a defined source community, and expanding thickness resolution.

- Experimental system: Two continuous-flow nitrifying moving bed biofilm reactors (MBBRs) operated for ~350 days with controlled maximum biofilm thickness using AnoxK Z-carriers (grid-covered saddles limiting biofilm to grid height). Five carrier types: Z50, Z200, Z300, Z400, Z500 (50–500 µm). Reactor 1 (3 L) contained Z200, Z300, Z400, Z500 (200 carriers each); Reactor 2 (1.5 L) contained Z50 (293 carriers). Both reactors fed unfiltered wastewater treatment plant effluent (Källby, Sweden) supplemented with 50 mg/L NH4-N and 0.5 mg/L PO4-P; nominal loading 2 g N m−2 d−1. Observed average loadings: 1.84 and 2.23 g m−2 d−1 with 88% and 92% N removal for R1 and R2, respectively. Operating conditions: 20°C, pH 7.5, DO 4.5 mg/L, HRT 2 h (to minimize homogenizing dispersal). - Sampling: Carriers sampled in duplicate at T1 (day 168 R1; day 123 R2) and T2 (day 275 R1; day 230 R2); single carriers at T3 (day 386 R1; day 341 R2) for Z50, Z200, Z400, Z500. Influent sampled five times per reactor from two batches. - Biomass processing and DNA methods: Biofilm biomass detached with sterile brush; suspended biomass collected by centrifugation. DNA extracted with MP Biomedicals FastDNA Spin Kit. Total cell numbers estimated by qPCR of 16S rRNA gene (primers 1055F/1392R) with SYBR Green chemistry; per-sample 16S rRNA gene copy number corrected to cell numbers using CaRcone and rrnDB-based average copy number prediction. - Amplicon sequencing: V3–V4 16S rRNA region amplified (PRK341F/PRK806R). Libraries sequenced on Illumina MiSeq. - Bioinformatics: DADA2 (v1.6) for QC, merging, and ASV inference; taxonomy via SILVA v128. Samples rarefied to 37,378 reads. Analyses in R using phyloseq, phangorn, vegan, Hmisc, picante, MicEco, ade4, adegraphics, ggplot2. - Neutral community modeling: Sloan’s neutral model used to quantify stochastic assembly from influent (source) to carriers (targets). Model run: all biofilms (ten runs of five carriers each; mean reported) and separately per thickness. N_T m (community size × migration rate) estimated; with qPCR-derived N_T, migration rate m was calculated using sampling correction (Sloan et al.). - Deterministic processes/phylogenetic structure: Compared phylogenetic alpha diversity (Faith’s PD) and beta-diversity metrics βMNTD/βMPD against richness-preserving null models (999 randomizations) to derive βNTI and βNRI. Negative values indicate clustering; values beyond |2| indicate significant deviation from null. - Community ordination: Double principal coordinate analysis (DPCoA) and constrained DPCoA; weighted UniFrac in phyloseq. - Data availability: Sequence data deposited at GenBank PRJNA322602.

- Neutral model fit varied by thickness: Biofilms ≥200 µm showed better neutral model fits than 50 µm biofilms (Spearman ρ: Z50 0.298; Z200 0.561; Z300 0.546; Z400 0.476; Z500 0.456; all biofilms combined 0.495). - Fraction of community consistent with neutral assembly: In thicker biofilms (≥200 µm), 86–94% of ASVs, comprising 70–75% of relative abundance, fit the neutral model. In 50 µm biofilms, only 74% of ASVs (42% relative abundance) fit, with many taxa under selection. - Selection in thin biofilms: Z50 had 87 ASVs selected for (22.7% relative abundance). In thicker biofilms, selected-for ASVs were fewer and lower in abundance (Z200: 48 ASVs, 2.4%; Z300: 36, 1.4%; Z400: 55, 3.9%; Z500: 62, 5.5%). - Migration rates (neutral model with qPCR-based N_T): Very low and similar across thicknesses, on the order of 10^−7 per birth event (approximately 7.4×10^−8 to 1.7×10^−7), indicating infrequent establishment from the influent and limited dispersal contribution in mature biofilms. - Habitat filtering across all biofilms: Phylogenetic analyses showed clustering relative to the influent (negative standardized effect sizes; average −6.0 to −7.3), with the strongest clustering trends at short phylogenetic distances and most pronounced in 50 µm biofilms. - Within-thickness phylogenetic structure: Thin biofilms (50 µm) were more similar to each other than expected by chance at short phylogenetic distances (homogeneous selection), whereas thicker biofilms tended toward phylogenetic overdispersion (βNRI) and greater phylogenetic beta-diversity, consistent with variable selection or stochastic drift effects among replicates. - Taxa driving differences: Biofilm communities, irrespective of thickness, were enriched in Nitrospira (Nitrospirae), indicating strong habitat filtering favoring nitrite-oxidizers; Nitrosomonas remained abundant across biofilms. - Functional implications (from prior linked work): Thin biofilms exhibited high nitrification efficiency; thicker biofilms showed lower nitrification efficiency but enhanced micropollutant transformation.

The findings demonstrate that biofilm thickness modulates the balance between stochastic and deterministic assembly processes in MBBRs. Thin (50 µm) biofilms experienced stronger deterministic selection, producing communities that were more phylogenetically clustered and more similar to each other than expected by chance. Likely drivers include stronger hydrodynamic shear at the biofilm surface, higher surface-area-to-volume ratios, and intense competition for limited space and resources in a relatively homogeneous environment, promoting selection for closely related, fast-growing, shear-tolerant taxa. In contrast, thicker biofilms (≥200 µm) exhibited better fits to the neutral model, indicating a larger role for stochastic processes in shaping community composition. These thicker biofilms also showed phylogenetic overdispersion at broader evolutionary scales and greater variability among replicates, patterns consistent with variable selection arising from heterogeneous microenvironments (e.g., diffusion gradients) and/or with drift and historical contingency given low migration rates. Across all thicknesses, habitat filtering was evident through enrichment of Nitrospira, reflecting selection imposed by nitrifying conditions. Very low migration rates derived from neutral modeling suggest that dispersal from the influent plays a minor role in mature biofilms, thereby enhancing the importance of drift in thicker biofilms once established. Despite functional steady state in nitrogen removal, community compositions varied over time within thickness groups, implying ongoing drift and/or shifts in selection rather than strict compositional steady states. Collectively, the results provide a mechanistic link between physical biofilm architecture (thickness) and community assembly pathways, with implications for managing biofilm-mediated treatment functions.

Assembly of MBBR biofilms involves both stochastic and deterministic processes whose relative contributions depend on biofilm thickness. Thin biofilms are dominated by homogeneous selection, yielding closely related and similar communities across replicates and high nitrification efficiency. Thicker biofilms show stronger signatures of stochastic assembly (notably drift) and/or variable selection, greater phylogenetic beta-diversity, and enhanced micropollutant degradation. Dispersal (migration from influent) appears to play a minor role in mature biofilms and does not vary with thickness, suggesting surface area and hydrodynamic conditions are more relevant drivers of colonization than thickness per se. These insights suggest that controlling biofilm thickness can be a lever to manage microbial community structure and function in engineered systems. Future work should quantify effective community sizes to refine migration estimates, track assembly dynamics over longer timescales to test for compositional steady states, and experimentally disentangle variable selection from drift in thick biofilms, ideally by testing neutral models and null-model phylogenetic metrics using well-defined source communities.

- Experimental design: The thinnest biofilms (Z50) were grown in a separate reactor and initiated 45 days later than the thicker biofilms, potentially introducing reactor effects despite identical operating conditions and stable influent communities. - Inference limits: Phylogenetic overdispersion and increased beta-diversity in thicker biofilms could arise from either variable selection or stochastic drift; the study cannot decisively distinguish between these mechanisms. - Migration rate estimation: Neutral-model-derived migration rates depend on accurate estimates of local community size (N_T). While qPCR-informed N_T was used with sampling corrections, uncertainties remain; effective community size concepts may better capture ecologically relevant N_T. - Temporal scope: Sampling focused on mature biofilms after functional steady state; earlier colonization phases were not captured, potentially missing periods of higher migration and different assembly dynamics. - Generalizability: The nitrifying selection imposed by ammonia supplementation and the specific hydrodynamic and loading regimes may limit transferability to other systems or selective contexts.