Microbial community overview of size-fractionated PN/A granules
The pilot-scale reactor had been in operation continuously for more than 500 days before the sampling. The influent was rich in ammonium (1724 ± 846 mg N L-1) and contained residual organic carbon (336 ± 116 mg COD L-1) (Table S1). The reactor was oxygen limited via continuous aeration resulting in dissolved oxygen DO values of 0.5-1.0 mg L-1. Stable PN/A performance was achieved in terms of an average total nitrogen removal rate of 0.82 kg N m-3 d-1 and efficiency of 87%.
The microbial community structure of the PN/A granule samples was characterized using the combination of 16S rRNA gene sequencing and metagenomic sequencing based on Illumina technology. Taxonomic analyses based on 16S rRNA gene sequencing showed that, microbial communities in granules with different sizes were dominated by members of Proteobacteria (24~46%), Chlorobi (14~34%), Bacteroidetes (11~20%), Chloroflexi (8~11%) and Planctomycetes (0.5~10%) (Fig. 1). These bacterial phyla are known to be relatively abundant and ubiquitous in PN/A granules [20, 22, 23] and these five phyla accounted for at least 87% of all sequencing reads in granules of all size-fractions in this study. In addition, it was observed that the relative abundance of these members changed with increasing of granule size. Proteobacteria, accounted for 46% of all sequences in <0.2 mm size fraction, while decreasing to 27% in >1.0 mm size-fraction. Similarly, Bacteriodetes were consistently enriched in the smaller size-fractions, where they constituted of 20% all sequences in 0.2-0.5mm size-fraction (versus 11% in >1.0mm size-fraction). In contrast, clades of Planctomycetes and Chloroflexi were consistently enriched in the larger size-fractions. Sequences matching WS6 and Acidobateria represented relatively minor components of the total amplicon sequences (1.6±1% and 1.3±1%, respectively), but showed more abundances in large size-fractions. A previous study has predicted a fermentative and denitrifying role in the anaerobic granule core for WS6 and Acidobacteria, respectively . It was surprising that each size-fraction also contained different bacterial communities, with a range of 2.3-7% of unique OTUs in each size-fraction (Fig. S2), supporting the fact that they harbor distinct habitats. For example, Zoogloea was only present in the <0.2 mm size-fraction, whereas Opitutus was only observed in the >1.0 mm size-fraction.
Although similar taxonomic distributions were observed in metagenomic datasets, some differences were found in terms of species abundance, such as a slightly higher representation of Chloroflexi, Planctomycetes and Acidobacteria in metagenomic sequencing compared to 16S rRNA gene sequencing. Also, some phyla were less represented (e.g., Bacteriodetes and Chlorobi) in metagenome than in 16S rRNA gene sequencing datasets. Such differences could be partly explained by biases derived from the PCR process in amplicon sequencing . Additionally, 16S rRNA gene recovered from the metagenomics data can span the entire length of the gene, whereas the PCR-based amplicon approach only targets specific region. The differences in accuracy of taxonomic assignments with various regions of the 16S rRNA gene may lead to the observed discrepancies with the two approaches .
Metagenome co-assembly and binning yielded 22 high-quality (> 80% completeness and <4% contamination rates) MAGs (Table 2). The phylogenetic inference (Fig. 2) showed that these 22 MAGs spanned 8 phyla, including Proteobacteria, Planctomycetes, Bacteriodetes, Chloroflexi, Verrucomicrobia, Acidobacteria, Chlorobi, and Omnitrophica. The recovered genomes in the PN/A granules mostly resemble those of previously detected organisms in other PN/A system . Among the recovered genomes were characterized taxa known to be involved in the canonical PN/A process: Nitrosomonas (AOB) encoded the key enzyme for ammonia oxidation (amo) (Fig. S3); and Candidatus Kuenenia (AMX) contained the core genes for hydrazine metabolism (hzs and hdh) and the gene cluster for nitrite reduction (nir), all essential for the anammox process. In addition, a few of organisms (e.g., PTB1-4, and CFX1-3) were capable of nitrate/nitrite respiration, which was evident from the genes encoding nitrate/nitrite reductase (e.g., nar and nir). Among them, only PTB3 and PTB4 encoded a complete denitrification pathway from nitrate to dinitrogen gas (Fig. S3), suggesting their potential roles in removing organic carbon and nitrate.
Larger size fractions harbor more diverse communities and sustain more redundant functions
In order to investigate the difference of the species diversity among size-fractions, we calculated alpha-diversity indices (including ACE and Chao1 for richness, and Shannon H and Simpson D for diversity) based on the 16S rRNA gene sequencing results (Table 1). Granule size explained the increasing number of OTUs (species) in larger size-fraction, ranging from 356 (<0.2 mm) to 434 (>1.0 mm). Similarly, the microbial richness, as indicated by the estimators ACE and Chao1, generally showing an increasing trend along with the granule size (except for 0.2-0.5mm size-fraction being the lowest). In addition to granule size, microbial richness could also be affected by biotic factors such as competition and facilitation, and abiotic factors such as nutrient availability and disturbance frequency . In terms of microbial diversity, Shannon H rose from 3.19 in <0.2 mm size-fraction to 3.83 in >1.0 mm size-fraction; whereas Simpson D showed the lowest value of 0.0735 for 0.8-1.0 mm size-fraction and the highest of 0.1347 for <0.2 mm size-fraction (overall showing a downward trend). Both the increasing Shannon H and decreasing Simpson D in larger granules indicate that larger size-fraction contain higher microbial diversity in this PN/A system.
According to the theory of island biogeography, there is a strong species-area relationship at macro-scale. As shown in Fig. 3A, taxonomic diversity contained within each size-fractions increased with increasing granule size, showing an approximately linear growth. The granules with >1.0 mm size-fraction present the most diverse community, whereas the <0.2 mm size-fraction had the lowest diversity. Considering the granule have a three-dimensional habitat distribution, which means microorganisms can grow on the inner layer of granules, rather than only growing on the surface of an island. Thus, we test whether the ecological distribution of microorganisms in granules follows a species-volume relationship (Fig. 3B). A sharp increase of the taxonomic diversity was observed in granules with volume smaller than 0.38 mm3 (<1.0 mm size-fractions), while there was a slow increase in taxonomic diversity in granules larger than 0.38 mm3 (>1.0 mm size-fraction). The elevated diversity in larger granule has been shown to coincide with previous studies of PN/A granules in lab-scale reactors [19, 20]. Taxa are usually adapted to a particular habitat, and large granules contain more distinct layered habitats due to the mass transport limitation of nutrient [30-33]. Therefore, large granules could support the survival of more species.
Functional diversity has been used to reveal the diversity and distribution of functional traits (or functional genes) across communities [34, 35]. Here we used Shannon’s diversity index to indicate functional diversity , estimated from the diversity of the protein-coding genes categories identified from annotated metagenomic reads. Similar to the distribution of taxonomic diversity among size-fractions, the communities found in small size-fractions typically had lower level of functional diversity (Fig. 3C and D). This is likely a result of that large granule having more anaerobic habitats suitable for the growth of anaerobes, therefore expressing a broader array of genes associated with anaerobic pathways (e.g., methanophenazine hydrogenation and dipicolinate synthesis). With granule size larger than 1.0 mm, interestingly, no increase in the functional diversity was observed despite increasing species richness along with the increasing granule size. It is widely assumed that community diversity could be decoupled with functioning [36, 37]. With a certain microbial diversity, additional species would not bring significant additional functional benefits . Multiple microbial taxon can carry out the same process within an ecosystem, thus resulting in the functional redundancy, where the system can perform just as well with fewer species. Considering both functional diversity and the taxonomic diversity, we found that the overall diversity of functional gene categories found in a given sample was, to some degree, predictable from taxonomic diversity of the microbial communities. It should be noted that the samples collected in this study were from one steady-state PN/A system. Future investigations on more granule-based systems might be required to confirm our observations of this work and would be better for a more comprehensive comparison.
Nitrogen conversion microorganisms and metabolic pathways in different size-fractionated granules
A PN/A bioreactor is a nitrogen-driven ecosystem, and its stable performance relies on achieving an appropriate balance between the activities of the principal microbial protagonists therein, including AMX and AOB, while suppressing or washing out nitrite-oxidizing bacteria (NOB). Here, we specifically analyzed the populations involved in nitrogen cycle of the five size fractions based on both 16S rRNA gene and metagenomic sequences. Anammox bacteria were proportionally enriched in larger size-fractions (Table S2). Sequencing matching anammox Candidatus Kuenenia with a similarity of 99.6% (classified as Planctomycetes in Fig. 1) peaked in relative abundance (9.4% of total) in >1.0 mm size-fraction, while decreased to 5.3% in 0.5-0.8 mm size-fraction and were nearly negligible in <0.2 mm size fraction (0.3%). Their elevated abundances in the larger size fractions are understandable since Candidatus Kuenenia is inhibited by oxygen , which favors more anaerobic space in larger granules due to a longer oxygen permeation depth [32, 33]. Conversely, the relative abundance of the sequences matching a typical AOB, i.e. Nitrosomonas (classified as Proteobacteria in Fig. 1) showed an opposite trend, peaking at 33.8% of the total <0.2mm size-fraction amplicons, compared with 5.7-6.4% in >0.8mm fractions. Since smaller granules have higher surface/volume ratio, which provide larger aerobic volume fraction , the enrichment of Nitrosomonas in small granules is consistent with these bacteria being aerobes adapted for the efficient use of oxygen. The presence of Nitrospira (one of the typical NOB), which could negatively affect nitrogen removal efficiency by competing nitrite with anammox bacteria, was negligible with an abundance of 0.01-0.3% across all size-fractions. Additionally, a typical heterotrophic denitrifier, Denitratisoma (classified as Proteobacteria in Fig. 1), was also observed in all the granules, which can utilize the organic carbon for nitrate/nitrite removal. The relative abundance of this organic carbon-scavenger increased from 3.8% (<0.2 mm-granules) to 8.6% (0.8-1.0 mm-granules), because larger granules harbour more anaerobic zones supporting their growth by taking up nitrite/nitrate as the electron acceptors.
Functional screening of the metagenomic and metatranscriptomic datasets revealed genes and transcripts involved in nitrogen metabolisms across size-fractions, which were correlated with the presence of the functional microorganisms (Fig. 4 and Table S2). In both the gene-abundance and gene-expression levels, the key nitrogen-related genes, were detected in all granules with different sizes. Specifically, the core genes and transcripts involved in hydrazine metabolism (hzs and hdh), which are all essential for anammox process, were expressed in larger granules, suggesting that larger size fractions predominantly performed the anammox process. In contrast, nitrification is the predominant nitrogen metabolism in the smaller granules, where the key genes and transcripts encoded for ammonia oxidation enzymes (amo and hao) were more abundant. Additionally, genes that encode catabolic enzyme specific for denitrification, including napA, nar, nirK, nor and nosZ, were low in abundance and expression compared with those involved in anammox and nitrification process. The counts of genes and transcripts generally show an increasing trend along with granule size, whereas 0.8-1.0 mm granules have the highest transcript hits.
Based on the established mathematical model, the microbial stratification inside the granules was also investigated (Fig. S4). Overall, anammox bacteria were shown to dominate the larger granules (0.6-1.2 mm) while AOB were predicted to be present in the smaller granules (0.1-0.3 mm. The model showed that the NOB were out selected in all granule sizes in the competition with anammox and heterotrophic denitrifiers. These model predictions are in accordance with metagenomic and metatranscriptomic results. The good match between the experimental and simulation results supported the relationship between microbial stratification and granule size. Specifically, in 0.1-0.3 mm granules, AOB occupy the granule surface where oxygen and ammonium are sufficiently present, while heterotrophs dominate the inner layer. In the granules with radius of 0.6-1.2 mm, anammox bacteria gradually develop at the inner side of the granule, with the relative distribution of AOB and heterotrophic bacteria unchanged at outer layer. Anammox bacteria are outcompeted by heterotrophs at outer shell of the granule mainly due to the inhibition of oxygen at the granule surface. The predicted microbial stratification inside biofilm (Fig. S4) remains to be experimentally verified by using advanced methods (e.g., 3D-FISH and microsensor).
Overall, the fraction-specific distribution was evident for key nitrogen-involving taxa, genes and transcripts, as well as the microbial stratification inside biofilm matrix. These patterns suggest that the small granule plays a dominant role in ammonium oxidation that is critically dependent on access to sufficient oxygen, and the removal of the produced nitrite will mostly rely on the anammox and heterotrophic denitrification in large granule. However, as large granules provide habitats for both aerobes and anaerobes, it is capable of complete nitrogen removal, but the availability of nitrite will restrict the overall nitrogen removal rate. Additionally, the stratified biofilm was determined by the availability and inhibition effect of specific substrates (e.g., nitrite and oxygen). Variable granule size distribution may drive heterogeneity in nitrogen removal performance in PN/A process.
Granular aggregate as a microscopic model system for island biogeography
Island biogeography theory is one of the most influential paradigms in ecology, which is central to several disciplines (e.g., conservation ecology and metapopulation theory) and continues to provide insights into research of island-like systems [41-44]. For example, Ali, Wang, Salam, Hari, Pronk, Van Loosdrecht, Saikaly  observed that the species sorting was significantly higher in flocs and small granules than large granules in aerobic granular sludge system, while putative functional populations were enriched with an increase in microbial aggregates size. Here, we adopted the theory of island biogeography as a conceptual framework for understanding the microbial ecology of microbes within a granule-based PN/A bioreactor. Granules, as insular habitats, provide favorable microscopic habitats for microorganisms. Viewing granules as islands in a bioreactor of potential colonists, increasing taxonomic diversity in larger size-fractions was shown to be congruent with the theory of island biogeography as the microbial richness was largely dependent on the size of the island  and larger islands house more taxa  (Fig. 3).
The observed increase in microbial diversity can be partially explained by the increased micro-habitats diversity and substrate-transmission complexity that associated with a larger granule. Firstly, larger granule provides better habitat differentiation and the enhanced heterogeneity, which foster more functional populations. Specifically, the large granules do have an anaerobic layer at the inner core that accommodates populations of other physiological types of bacteria. Secondly, each of the granule houses a micro-ecosystem that derives its nutrients and energy from wastewater. The availability of nutrients/substrate that may favor different bacterial species to colonize and perform differently. In PN/A reactor, the availability of oxygen, to a large extent, differentiates the microbial communities of small granules from large granules. The substrate (e.g., oxygen) drawdown from the outer layer to the inner core of granules significantly restructures the microbial ecosystem in granules, resulting in the differences in the abundances and the activities of aerophilic or anaerobic microbes. Large granules offer more favourable habitat for anaerobic bacteria (e.g., heterotrophic bacteria) that are capable of anaerobically decomposing organic matter within the anoxic core of the granules. Aerobic COD degradation in small granules and on the surface of large granules, together with anaerobic COD degradation in the core of large granules contributed to the COD removal in this PN/A system (as indicated by the COD removal efficiency in Table S1).
It is worth to note that PN/A granules are recirculated in the bioreactor, which make the scenario of the interaction and overview of microbial ecology in granules more complicated. Taking the bioreactors as a whole, where a group of granules closely scattered in the bioreactor, numerous of granules form an archipelago and there were reciprocal interactions between these islands. Under such circumstance, a higher immigration rate of species is expected, which might reduce the species differences between size-fractions.
Overall, granules function as islands for housing microorganisms, in which larger granules not only harbor higher microbial diversity, but also support more functional diversity. The theory of island biogeography is instructive to understand the microbial ecology associated with granule. A better understanding of the microbial ecology in PN/A reactor could guide us to better design and operate bioreactors to foster the development of specific microbial communities that can support the desired functional processes .
Integrating nitrogen cycling microorganism and activities in granules
As indicated by the different expression of aerobic and anaerobic ammonium oxidation genes in size-fractionated anammox granules (Fig. 4), an ecological model based on the granule PN/A system with size distribution was delineated in Fig. 5. The small-size granules are dominated by AOB and heterotrophs, which catalyze aerobic oxidation of ammonium to nitrite and aerobic oxidation of organic carbon to CO2. On the other hand, both aerobic and anaerobic microorganisms (AOB, anammox and heterotrophs) are observed in large-size granules, in which ammonium, produced nitrite and nitrate are removed to dinitrogen gas. It is also worth to note that NOB are rarely detected in all granules, thus the nitrate involved is mainly produced from anammox reaction and partially consumed by heterotrophic denitrification .
In small granules, AOB contribute to ammonium oxidation while heterotrophs contribute to COD oxidation. As the oxygen permeates along the entire small granules (0.1-0.3 mm) (Fig. S4), the anaerobic processes (anammox and heterotrophic denitrification) are both inhibited in small granules. Thus, the produced nitrite will escape from small granules and transfer to the larger granules. In large granules (0.6-1.2 mm), the oxygen permeation from surface ceases after 300-400 μm, and the anaerobic layer can expand. In the anaerobic layer, both anammox and heterotrophic denitrification work together to remove nitrite either transferred from small granules or from the AOB growing at outer layer of the large granule itself. It is worth to note that the kinetics of anammox and heterotrophic denitrifiers selected were slightly higher than the average values . Therefore, NOB are out-selected during the competition for nitrite with AMX and heterotrophic denitrifiers in the present model. In the real reactor, the out-selection of NOB could be due to the high free nitrous acids (FNA) and free ammonia (FA) concentrations under sidestream condition .
To quantify the nitrogen flows among granules, mass balance was analyzed based on the simulation results and the calculated biomass volume of each granule size. The granules with radius of 0.1 mm, only contributed 4.2% to nitrite production which can be removed inside the same granules while the majority (95.8%) will be transferred to the larger granules with radius of 0.6 mm (14.8%), 0.9 mm (46.5%) and 1.2 mm (34.5%). Similarly, about 14.4%, 45.3% and 33.7% of the produced nitrite in granules with radius of 0.3 mm will be transferred to 0.6 mm, 0.9 mm and 1.2 mm granules, respectively (Fig. 5). This cross-feeding of nitrite from small to large granules indicates the mutualism and competition in PN/A system occur not only inside granules, but also among granules with different sizes . However, this cross-feeding of substrate is unidirectional (Fig. 5). The granules with radius of 0.6-1.2 mm could offer both aerobic and anaerobic conditions, which can foster AOB, anammox bacteria and heterotrophic bacteria at different micro-habitats (Fig. S4). By utilizing the ammonium, oxygen and organic carbon from the feeding, and the nitrite transferred from the smaller granules, these large granules can develop a stratified microbial community to utilize all these inputs by themselves with a balance. Thus, there is no substrate transferred from the large granules to the small granules. However, these results do not suggest that larger granules are better than smaller granules. Indeed, it was observed that the N2O production increased significantly in larger granules due to incomplete denitrification . The optimal granule size in practice is determined based on the balance between nitrogen removal performance, greenhouse gas emission (N2O) and also the operational cost (maintaining large granules requires more energy consumption), which warrants further studies.
Although the microbial ecology of anammox granules was studied in previous studies [22, 23, 25, 49, 50], the size-fractioned anammox granule system proposed here delineates the interactions among granules with different sizes from micro-scale view for the first time. In a recent study, the cross-feeding among microorganisms was investigated based on metagenomics analysis . However, the present study touched another dimension to demonstrate the interactions between granules based on more comprehensive data (metagenomics and metatranscriptomics). The cross-feeding model proposed here offers an explanation for the granule size distribution impacts on PN/A performance.