Nitrogen loading effects on nitrification and denitrification with functional gene quantity/transcription analysis in biochar packed reactors at 5 °C

This study investigated the nitrogen transformation rates of different nitrogen-loading (20, 30, and 50 mg TN/L) biochar packed reactors (C:N:P = 100:5:1) within 125 days at 5 °C. The results showed that high nitrogen loading resulted in an NH4+ (TN) removal efficiency decline from 98% (57%) to 83% (29%), with biochar yielding a higher NH4+, TN and DON removal rate than conventional activated sludge. Moreover, all biochar packed reactors realized a quick start-up by dropping in temperature stage by stage, and the effluent dissolved organic nitrogen (DON) concentrations of R20, R30, and R50 were 0.44 ± 0.18, 0.85 ± 0.35, and 0.66 ± 0.26 mg/L, respectively. The nirS/amoA, nxrA/amoA, and amoA/(narG + napA) were deemed to be the markers of ammonium oxidation rate (SAOR), specific nitrite oxidation rate (SNOR), and specific nitrate reduction rate (SNRR), respectively. Compared with functional gene quantity data, transcription data (mRNA) introduced into stepwise regression analyses agreed well with nitrogen transformation rates. High nitrogen loading also resulted in the cell viability decreased in R50. Nitrogen loadings and operation time both led to a significant variation in cell membrane composition, and unsaturated fatty acids (UFAs) significantly increased in R30 (46.49%) and R50 (36.34%). High-throughput sequencing revealed that nitrogen loadings increased the abundance of nitrifying bacteria (e.g., Nitrospira) and reduced the abundance of denitrifying bacteria (e.g., Nakamurella, Thermomonas, and Zoogloea) through linear discriminant analysis (LDA).

Nitrogen-loading calculations are important for nitrogen-removal upgrades of WWTPs, and depend on hydraulic residence time or tank volume determination 13 . Nitrogen loading significantly influences NH 4 + and TN removal 14 , which determines wastewater processing capacity. To prevent low nitrogen removal efficiency, it is necessary to set the low nitrogen loading at temperatures lower than 15 °C 15 . Under low nitrogen loading conditions (<15 mg TN/L), functional genes (e.g., amoA, nirS, nxrA, napA, and narG) involved in nitrogen metabolism measured by the quantitative polymerase chain reaction (qPCR) were used to calculate and predict NH 4 + , NO 2 − , and NO 3 − transformation rates (p < 0.05) at 4 °C 16 . However, the rationality of this method for different nitrogen loadings deserves further study, which is one of the main purposes of this research.
In this study, we assessed the effects of nitrogen loading on nitrification and denitrification in biochar packed reactors at 5 °C. The main objectives of our investigation were: (i) to compare the performance of fundamental biochar-packed reactors with different nitrogen loadings; (ii) to assess biofilm performances at different total nitrogen loading, resulting in relative functional gene transcription, in terms of its capability of reducing ammonia nitrogen, nitrate, and nitrite; and (iii) to eventually investigate the association between bioreactor performance, nitrogen transformation rate, and the microbial community.

Results
Bioreactor performance. Figure 1 shows the fundamental effluent quality with different nitrogen loading (20,30, and 50 mg TN/L) within a 125-day operation period. In phase I, along with a temperature decline (from 20 °C to 5 °C), NO 2 − and NH 4 + showed a short-term accumulation in 50 mg TN/L nitrogen loading. NO 3 − initially accumulated with all nitrogen loadings, and effluent NO 3 − exhibited a greater decrease in R 20 and R 30 after the 15 th day. In phase II, the biofilm microorganisms adapted to 5 °C, and the fundamental effluent quality remained stable from day 30 to day 125. Effluent TN (NO 3 − ) concentrations of R 20 , R 30 and R 50 declined and stabilized to 8.49 ± 0.55 mg/L (7.81 ± 0.44 mg/L), 15.42 ± 1.13 mg/L (11.74 ± 0.72 mg/L) and 35.69 ± 1.42 mg/L (26.72 ± 0.95 mg/L), respectively, on the 30 th day (Fig. 1E). NH 4 + removal efficiency decreased from 89% (R 20 ) to 83% (R 50 ) as the influent nitrogen loading increased (Fig. 1B). Moreover, effluent DON concentration was found to be lowest in R 20 (0.44 ± 0.16 mg/L). In this study, the high dissolved oxygen (7.6-8.0 mg/L) severely limited toxic NO 2 − accumulation (<0.2 mg/L) during stabilization (Fig. 1C), even with 50 mg TN/L nitrogen loading. With the increase of influent COD loading (R 20 : 400 mg/L, R 30 : 600 mg/L and R 50 : 1,000 mg/L), effluent COD concentration ranged from 50 mg/L (R 20 ) to 125 mg/L (R 50 ) (Fig. 1A), and the COD removal rate reached its peak at 88.10% (R 30 ).
Integrity and composition of cell membrane. The cell membrane integrity rate indicated by the cell viability and fluorescence microscopy photographs were determined through BacLight Staining 7 (  Supplementary Table S4). Principal coordinate analysis (PCoA) was performed for the unweighted and weighted UniFrac distances to determine if PLFA construction shifted with different sampling times and nitrogen loadings. Significant variations existed in PLFA construction between three different nitrogen loadings (Fig. 3). This result confirms that nitrogen loading is an important factor in PLFA construction. Moreover, the PLFA construction diversity significantly changed with respect to operation time, particularly in 50 mg TN/L nitrogen loading. Unsaturated fatty acids (UFAs) (e.g., 18:2 cis 9,12, 20:4 cis 5, 8,11,14) significantly increased in R 30 (46.49%) and R 50 (36.34%) compared with R 20 (28.57%) on day 125 (Table S4).  MiSeq sequencer to fully understand the impacts of long-term biochar packed reactors operation at 5 °C (see Supplementary Tables S2 and S3). All samples of high-throughput Illumina sequencing were normalized to 9,648 OTUs, corresponding to the sample with the least number of reads. At the genus level, the OTUs, which showed up to 2% relative abundance, were selected to generate the heat map (Fig. 5A). After day 30 (phase II), when the reactors were at the peak of their nitrogen removal ability, the bacterial communities of the same reactor were similar on days 35, 75, and 125. Proteobacteria was the most abundant phylum in R 20 , R 30   Redundancy analysis (RDA) was used to evaluate the relationship between reactor environmental factors and bacterial population structure (>2‰) on days 35, 75, and 125 by relative abundance at the genus level (Fig. 5B). RDA1 and RDA2 explained 53.50% and 20.90% of the total variance, respectively. The distribution of bacterial community structures was clearly divided into three groups, R 20 , R 30 , and R 50 , because the angles between adjacent groups were nearly perpendicular, indicating that they were unrelated (Fig. 5B) Linear discriminant analysis (LDA) determines the features most likely to contrast change between classes by coupling standard tests for statistical significance with additional tests encoding biological consistency and effect relevance 17 . As shown in Fig. 4C, nitrogen loading increased the abundance of nitrifying bacteria (e.g., Nitrospira) and reduced the abundance of denitrifying bacteria (e.g., Nakamurella, Thermomonas, and Zoogloea).

Discussion
Low temperatures (5 °C) had little effect on ammonia nitrogen transformation of biofilm, and the high NH 4 + removal rate of biofilm (>98%) was usually reported at 4-12 °C 16,18,19 . The considerable effluent COD (>50 mg/L) of R 30 and R 50 implied an inadequate biofilm aerobic oxidation and denitrification ability at a low temperature 20 . Inadequate biofilm denitrification ability contributed to NO 3 − accumulation and the unsatisfactory TN removal rate of R 30 and R 50 , which was also reported in an activated sludge system (25% TN removal rate at 10-15 °C) at low temperatures 21,22 . Under the same temperature and nitrogen loading, R 20 represents a higher NH 4 + (98%) and TN (57%) removal rate than conventional activated sludge (NH 4 + removal rate: 92%; TN removal rate: 33%) 23 . Biochar addition promoted effluent NH 4 + and TN removal, which resulted from nitrifier/denitrificans immobilization below 8 °C and the increase of functional genes abundance (amoA, napA, nxrA, nosZ, narG, nirK, and nirS) 7,24 . Moreover, this might relate to many different functional groups on the surface of biochar, such as nitro, chloro, hydroxyl, amine, carbonyl and carboxylic 25 . Although more than 20 mg TN/L nitrogen loading facilitated a biofilm denitrification problem at 5 °C, a biochar packed reactor was helpful for DON removal with 20 mg TN/L nitrogen loading compared with traditional activated sludge reactors 23 . DON removal improvement from traditional activated sludge (3.45 ± 0.41 mg/L) to biochar (0.44 ± 0.16 mg/L) is essential in WWTPs because nitrogenous compounds are a potential threat to aquatic ecosystems 23 .
According to previous studies [26][27][28][29][30][31][32] , functional gene abundance can be calculated and effectively introduced into stepwise regression analyses with nitrification and denitrification transformations (e.g., SAOR, SNOR and SNNR). Specifically, for two consecutive steps of nitrification and denitrification, the amoA and nirS genes are thought to be NH 4 + to NO 2 − oxidation and NO 2 − to NO reduction markers, respectively 26,27 . Low NO 2 − concentration is beneficial to NH 4 + transformation, because NO 2 − has toxic effects on ammonia-oxidizing bacteria 28,29 . Similarly, nxrA is regarded as the NO 2 − to the NO 3 − oxidation marker involved in the nitrite oxidation process 30 . Both narG and napA are regarded as marker genes for NO 3 − into NO 2 − reduction, and are responsible for the first denitrification step. Generally, narG is dominant under anoxic and napA under oxic conditions 31,32 . Therefore, nirS/amoA, nxrA/amoA, and amoA/(narG + napA) were deemed to be the markers of SAOR, SNOR, and SNNR at 4 °C 16 .
According to the functional gene quantity and transcription data (see Supplementary Figs S2 andS3), we performed stepwise regression analyses between qPCR data and nitrogen transformation rates (SAOR, SNOR, and SNNR) on days 7, 15, 35, 75, and 125. As shown in Fig. 5 and Table 1, when functional gene quantity and transcription data (mRNA) were both introduced into stepwise regression analyses, the results agreed well with nitrogen transformation rates in 20 mg TN/L nitrogen loading (Pearson test, p < 0.05), whereas when only functional gene transcription data were introduced into stepwise regression analyses, the results did not agree well with nitrogen transformation rates in R 30 and R 50 (Pearson test, p > 0.05). The specific reasons are as follows. Extractive DNA originated from both live and dead cells, but extractive RNA only came from live cells. We hypothesize that DNA from dead cells contributed to the inaccuracy in the fitting between functional gene abundance and nitrogen transformation rates. We measured the cell viability in different nitrogen loadings with different sampling times, and the results showed obvious variations in the cell viability between different nitrogen loadings. Limited by cell viability variation, functional gene transcription data (mRNA) better fitted nitrogen transformation rates than functional gene quantity data, especially in R 30 and R 50 .
Low temperatures contribute to cell membrane transfer rate decreases 33 , and substrate mass transfer rate was another important factor for nitrogen transformation. For example, substrate (e.g., NH 4 + ) mass transfer depends on an ammonia transporter protein on PLFA. Meanwhile, nitrification and denitrification enzymes and electron transport of coenzymes (MQH 2 and QH 2 ) are located in the cell membrane 34 . The results showed that nitrogen loading and operation time both changed PLFA construction. From the observations shown in Fig. 4B, UFAs were significantly higher in R 30 (46.49%) and R 50 (36.34%) than in R 20 (28.57%), which enhanced the mass transfer efficiency of cells, and is necessary for efficient mass transfer at low temperature 35 . Moreover, Eliška, R. & Kateřina, H found that nitrogen loadings result in a PLFA content change 36 , which is often used to indicate the microbial community 23,33,35 . The cell viability decreased as the nitrogen loading increased, which might result from effluent metabolites, extracellular polymeric substances, and the microbial community structure. However, microbial community structure also plays an important role in nitrification and denitrification at low temperatures 37 .
Biochar promoted the adsorption of NH 4 + and NO 3 −12,38 , increased the abundance of nitrifying bacteria and the NO 3 − removal rate 39,40 , and changed microbial community construction for better nitrogen removal 41,42 . Pseudomonas, Janthinobacterium, Arthrobacter, and Flavobacterium occupied the highest abundance in R 20 , R 30 , and R 50 , which agreed with the activated sludge system at low temperatures 23,43 . Among them, Pseudomonas, Arthrobacter, and Flavobacterium were enriched in biochar with a high ammonification ability in the soil 44,45 . This bacterial genus enrichment might be a result of specific biochar surface characteristics, which concurrently partition microorganism spatial distribution, and prompt microorganism absolute abundance, as well as microbial community's structural stability 46 . Serratia only occupied high abundance in R 30 and R 50 , which could produce extracellular polymeric substances 47 (e.g. protein, polysaccharide, nucleic acid). Protein and nucleic acid are an important DON source after cell death from hydrolysis 22 . Curvibacter had positive correlation with NO 2 − , SNOR, and SNNR in Fig. 5B, which could use both nitrate and nitrite as an electron acceptor under both aerobic and anoxic conditions 48,49 . DON had no significant correlation with any bacterial genus, which was consistent with our previous study 22 , and no single bacterial genus could degrade DON effectively. In conclusion, we studied the impacts of nitrogen loading on nitrification and denitrification in biochar packed reactors at 5 °C. Three significant findings show the novelty of this research. (1) Biochar represents a high NH 4 + , (TN, and DON) removal rate at 5 °C. High nitrogen loading resulted in an NH 4 + (TN) removal efficiency decline from 98% (57%) to 83% (29%), but biochar represented a higher NH 4 + , TN and DON removal rate than conventional activated sludge. (2) Nitrogen loading has effect on cell viability, UFAs, and microbial community at 5 °C. The cell viability, UFAs and microbial community under different nitrogen loading is relatively less researched compared to other nitrogen-removal research at 5 °C. (3) Compared with functional gene quantity data, the transcription data (mRNA) introduced into stepwise regression analyses agreed agrees well with nitrogen transformation rates. This is an area that needs more research, but based on this research and the research that has preceded it 16,46 , it is well worth pursuing.

Materials and Methods
Bioreactor. Three SBRs (R 20 , R 30 , and R 50 ) fed with different total nitrogen loading (20,30, and 50 mg/L) were operated at 5 ± 0.1 °C (incubator: Sanyo Electric) with a working volume of 2 L. Seeding sludge was obtained from the aeration tank (suspended solids concentration about 4,400 mg/L) of a municipal wastewater treatment plant in Nanjing, China. Biochar was made from local rice husk charred at temperatures 600-800 °C for 8 h, and had a diameter of 0.5-3.0 mm and a specific surface area of 155.51 m 2 /g. This biochar was added to the experimental reactor (10 g/L). Effluent was discharged from the middle port of the reactor with a volumetric exchange ratio of 50%. The hydraulic retention time of the SBRs was 12 h and the feeding time, reaction time, settling time, and anoxic time were 0.5 h, 10 h, 0.5 h, and 1 h, respectively. The reactors were generally operated with fixed mixed liquor suspended solids (MLSS) with a concentration of 4,000 mg/L and an oxygen concentration of 7.6-8.0 mg/L. Synthetic wastewater. The synthetic wastewater took glucose, NH 4 Cl, and KH 2 PO 4 as carbon, nitrogen, and phosphorus sources, respectively, and the mass ratio of C:N:P was 100:5:1 (Table 3). Furthermore, the synthetic wastewater was set with different TN loadings (20, 30 and 50 mg/L). One liter of synthetic wastewater contained 0.6 mL (R 20 ), 0.9 mL (R 30 ), and 1.5 mL (R 50 ) of metals solution, and one liter of the metals solution contained 11 g CaCl 2 ·2H 2 O, 1.5 g FeCl 3 ·6H 2 O, 0.15 g H 3 BO 3 , 0.03 g CuSO 4 ·5H 2 O, 0.18 g KI, 0.12 g MnCl 2 ·4H 2 O, 0.06 g Na 2 MoO 4 ·2H 2 O, 0.12 g ZnSO 4 ·7H 2 O, 0.15 g CoCl 2 ·6H 2 O, and 10 g EDTA. Furthermore, NaHCO 3 was added to keep the alkalinity in the SBRs in the range of pH = 7.0 ± 0.5. Fluorescent staining. The AS samples were obtained from running SBRs (reaction time). Immediately, the samples were stained using a 2 μL LIVE/DEAD BacLight bacterial viability kit, and were then observed with an epifluorescence microscope after the samples had been diluted 10 times. The BacLight bacterial viability kit consists of two stains, propidium iodide (PI) and SYTO9, which were both used to stain the nucleic acids. Green fluorescing SYTO9 was used to permeate all cells, which enabled us to acquire total cells counts, whereas red fluorescing PI only enters cells with damaged cytoplasmic membranes. The obtained microscopic images were visualized with Image J software (a Java-based image processing program, National Institutes of Health, USA). After taking an average of 20 samples, we counted the live/dead cells rate 7 .
Microorganism molecular analysis. DNA and RNA extraction. DNA extraction of 0.5 g from each sample was performed using a FastDNA ® SPIN Kit for soil (MP Biomedicals, OH, USA) following the manufacturer's protocol.
An E.Z.N.A. TM Cycle-Pure Kit (Omega Bio-tek Inc., USA) was used to extract and purify total genomic DNA. RNA extraction (0.5 g suspended solids per sample) was performed using a TransZol Up Plus RNA Kit (TransGen Biotech, China). RNA samples were purified with MagicPure TM RNA Beads (TransGen Biotech, China), and reverse transcription to cDNA through an EasyScript ® One-Step gDNA Removal and cDNA Synthesis SuperMix (TransGen Biotech, China). The amount and purity of DNA, RNA, and cDNA were determined using a NanoDrop ® Spectrophotometer ND-1000 (Thermo Fisher Scientific, MA, USA) based on an absorbency of A260 and a ratio of A260:A280. Extracted DNA, RNA, and cDNA were checked by 2% agarose gel electrophoresis and stored at −80 °C.
High-throughput sequencing of 16 S rRNA gene. PCR amplification was carried out using the forward primer (5′-ACT CCT ACG GRA GGC AGC AG-3′) and reverse primer (5′-TCT CAN VGG GTA TCT AAT CC-3′) for the V3-V4 hypervariable region of the 16 S rRNA gene 51 . PCR was performed in a Veriti ® 96-Well Thermal Cycler (Applied Biosystems, USA). The quadruplicate PCR reactions for each sample preparation were purified again using the E.Z.N.A. TM Cycle-Pure Kit (Omega Bio-tek Inc., USA) and stored at −80 °C. About 500 ng of the purified PCR product for each sample was mixed and sequenced on an Illumina MiSeq sequencer. After sequencing, Python scripts were written to perform quality filtering of the raw reads with the Sickle and Mothur programs to remove low-quality sequences and reduce noise. The filtered sequences were assigned to taxa by the RDP classifier. The detailed data analysis is available in the Supplementary Methods and Results.
Other methods. Effluent  , and TN were determined as per the Standard Methods 52 . The differences between TN and total inorganic nitrogen (NH 4 + , NO 2 − , and NO 3 − ) was determined to be DON 53 . The determination of SAOR, SNOR, and SNRR measurements followed Wang et al. 54 . The PLFA extraction and analysis method was according to our previous study (2018).