Jul 16, 2021 Long-term operation assessment of a full-scale membrane-aerated biofilm reactor

• Data from a full-scale hybrid MABR op- erating under Nordic conditions were analyzed. • OTRispositivelycorrelatedwithammo-nium concentrations. • A nitrifying bio ﬁ lm could be established within 3 weeks. • Multivariate analysis revealed different operational patterns of behavior. • MABR shows lower energy consumption and footprint than conventional systems. Membrane-aerated bio ﬁ lm reactor (MABR) technology is an exciting alternative to conventional activated sludge, with promising results in bench and pilot-scale systems. Nevertheless, there is still a lack of long-term and full-scale data under different operational conditions. This study aims to report the performance of a full- scale hybrid MABR located in the North of Europe. In ﬂ uent, ef ﬂ uent, and exhaust data were collected for 1 year(September2019to September 2020) using online sensors/gas-analyzersand off-linelaboratoryanalysis. Next,oxygentransferrate(OTR),oxygentransferef ﬁ ciency(OTE),andnitri ﬁ cationrates(NR)werequanti ﬁ edas process indicators. Finally, multivariate methods were used to ﬁ nd patterns among monitored variables. Obser- vationsrevealedthat lowerair ﬂ ows achievedhigherOTE atthesamevalues of OTR and OTR wasstronglycorre-latedtoammonia/um concentration intheMABR tank (NH x,eff ).The dynamics betweenoxygen concentration in the exhaust (O 2,exh ) and NH x,eff indicated that a nitrifying bio ﬁ lm was established within 3 weeks. Average NR were calculated using four different methods and ranged between 1 and 2 g N m − 2 d − 1 . Principal component analysis (PCA) explained 81.4% of the sample variance with the ﬁ rst three components and cluster analysis (CA)dividedtheyearlydatainto ﬁ vedistinctiveperiods.Hence,itwaspossibletoidentifytypicalNordicepisodes with high frequency of heavy rain, low temperature, and high variations in pollution load. The study concludes that nitri ﬁ cation capacity obtained with MABR is robust during cold weather conditions, and its volumetric value is comparable to other well-established bio ﬁ lm-based technologies. Moreover, the aeration ef ﬁ ciency (AE) obtained in this study, 5.8 kg O 2 kW h − 1 , would suppose an average reduction in energy consumption of 55% compared to ﬁ ne pore diffused aeration and 74% to the existing surface aeration at the facility.

Wastewater treatment plants located in areas with large seasonal variations in temperature, such as the Nordic region, face challenges related to low water temperatures (T) and substantial influent dilution during the colder months to successfully achieve nutrient removal. Moreover, there is an ongoing trend toward stricter effluent requirements and more sustainable operation, putting pressure on finding technological solutions that can simultaneously accommodate those demands.
Nitrification robustness is essential for facilities with stringent nitrogen effluent requirements during cold weather conditions. Nitrification has a high temperature sensitivity; for every 10°C drop in temperature, the maximum specific growth rate for nitrifiers halves, and therefore the minimum solids retention time doubles (Henze et al., 2008). Operating at longer solids retention times can increase the biological nutrient removal performance during cold weather periods. However, operating at longer solids retention times -15-30 days for oxidation ditches, common in Northern Europe-requires more extensive infrastructure and often results in higher energy consumption and poorly settling sludge (Rittman and McCarty, 2012). Biofilm-based technologies can help utilities maximize the use of their existing infrastructure. Indeed, the addition of media to support biofilm growth in activated sludge can be traced back to the 1930s (Stensel and Makinia, 2014). Nonetheless, conventional co-diffusional biofilms typically suffer from mass transfer limitations, resulting in higher dissolved oxygen concentration requirements and higher energy requirements (Metcalf and Eddy, 2014).
Membrane-aerated biofilm reactor (MABR) technology has gained considerable interest recently, emerging as an alternative to conventional biological nutrient removal with the potential to combine the benefits of biofilm-based technologies and low-energy aeration supply (Martin and Nerenberg, 2012). The MABR is based on the use of counter-diffusional biofilms, in which air is supplied through hollow-fiber or flat-sheets membranes, which at the same time function as biofilm support. Since microbial activity on the surface of the membrane acts as a catalyzer for oxygen transfer , the oxygen transfer efficiency (OTE) in the field can be higher than in clean water, with the potential of achieving significant reductions in energy requirements compared to diffused aeration (Castrillo et al., 2019).
Besides the well-known advantages of biofilms for wastewater treatment, such as the possibility of working at high biomass concentrations and the ease of separation between biofilm and product (Van Loosdrecht and Heijnen, 1993), the counter-diffusional nature of MABR offers additional benefits. In counter-diffusional biofilms, an aerobic and an anoxic/anaerobic zone can coexist within the biofilm, providing the environmental conditions to simultaneously sustain both nitrification and denitrification (Timberlake et al., 1988). Furthermore, this allows the MABR to be retrofitted into existing anoxic activated sludge zones to provide additional nitrification capacity. This type of MABR configuration is referred to as "hybrid", in which anoxic activated sludge surrounds the MABR membranes and can result in complete total nitrogen removal, i.e., achieving removals approaching 100% at high organic loadings (Downing and Nerenberg, 2008).
The development of MABR technology started in the 1970-80s first as a means to achieve efficient oxygen supply (Yeh and Jenkins, 1978;Wilderer et al., 1985;Côté et al., 1989), to later combine the function of the membrane as biofilm support (Timberlake et al., 1988;Abdel-Warith et al., 1990). After decades of research to understand the fundamentals and deal with the challenges associated with biofilm control (Brindle et al., 1998;Casey et al., 1999), following studies focused on taking advantage of its counter-diffusional nature, for example, to achieve total nitrogen removal (Downing and Nerenberg, 2008), to reduce nitrous oxide emissions (Kinh et al., 2017), or to establish alternative mainstrem processes such as nitrite shunt (Mehrabi et al., 2020), partial nitrification and anammox (Bunse et al., 2020) and simultaneous nitrification-denitrification and anammox . Some of the most recent developments include improved mathematical modeling (Acevedo Alonso and Lackner, 2019), studying the influence of protozoa predation (Kim et al., 2020) and the importance of membrane material selection (Wu et al., 2019). Despite these recent advancements, operational performance results from full-scale or pilot-size installations are still very scarce (Nerenberg, 2016), and the effect of seasonal variations on MABR performance under full-scale operation in cold areas has yet to be addressed.
This study presents results from one year of operation of a full-scale hybrid MABR in an existing enhanced biological phosphorus removal (EBPR) type of configuration located in the Nordic region (Denmark). The main objectives of the study are: (a) assess the overall performance with regards to nitrification capacity, oxygen transfer, and required startup time, among others (b) understand how several operational conditions, including seasonal variations, impact the performance, and (c) assess the methods employed to monitor performance in full-scale installations. The results presented in this study, one of the first of its kind, will provide important insight into the performance and monitoring of full-scale MABRs as well as shed some light into the factors that govern its behavior. Volumetric air loss between inlet and outlet w A i r flow rate, kg s −1 Xo 2,in/out Mol fraction of oxygen in atmospheric air and in the MABR exhaust ρ o2 Oxygen density under normal conditions, kg m −3

Ejby Mølle WRRF
The Ejby Mølle Water and Resource Recovery Facility (EMWRRF) in Odense, Denmark, has a 410,000 population equivalents treatment capacity. The main liquid treatment train (Fig. 1, units 2-9) is comprised of grit removal and screening (6 mm), chemically enhanced primary treatment with the addition of polymers and ferric sulfate, and Bio-Denipho™ nutrient removal with phased isolation ditches including anaerobic zones for biological phosphorus removal and tertiary sand filters. There is also a secondary line treating wet weather wastewater (Fig. 1, units 10-15), including storm tanks and first and second step trickling filters with intermediate secondary clarifiers. Sludge processing (Fig. 1, consists of waste activated sludge thickening, anaerobic digestion with a combined heat and power biogas cogeneration system to produce electricity and heat, and digested sludge dewatering before hauling sludge offsite for composting. Dewatering centrate is treated in a granular sidestream deammonification system. Fig. 1 shows a simplified process flow diagram of the EMWRRF and the MABR pilot location and setup.

MABR demonstration setup, operational conditions & process characteristics
The MABR reactor consisted of a sidestream 25 m 3 circular reactor adjacent to the EBPR bioreactor facility's anaerobic zone at the EMWRRF. A full-scale hollow-fiber MABR cassette with a total volume of 11.3 m 3 and a total membrane surface area of 1920 m 2 was installed inside the pilot reactor in 2019 (see Fig. 1A). The reactor was set up as a continuous stirred-tank reactor fed with mixed liquor from the fullscale anaerobic zones (i.e., primary effluent with return activated sludge). The feed was pumped from the anaerobic zone and introduced to the tank near the bottom using a distribution grid. A recirculation pump was utilized to ensure proper mixing of the tank. The effluent was located at the top of the tank (see Fig. 1B, C, D). The main influent/operational characteristics are summarized in Table 1.
Low-pressure air was supplied to the MABR unit for process oxygen supply and mixing. Mixing airflow was maintained at 12 m 3 h −1 throughout the study period due to concerns regarding the feed's high suspended solids concentrations (See Table 1). Process and mixing airflows to the MABR were adjusted manually.
External aeration (EA) was provided with fine bubble diffusers on one side of the tank to counteract the feed's extremely low redox conditions (often less than −300 mV), which had created problems in previous phases of the study. The EA was operated intermittently and controlled using two different strategies: time-based control during the first 3 months and using an ORP-based feedback control during the study's remainder (see Fig. 1 E, F, G). The ORP-based feedback control also included managing the feed flow to the tank: Decreasing the flow when the ORP reached a pre-defined setpoint.
Pressure before and after the units was measured using a pressure transmitter connected to the SCADA system. Airflow before the MABR was measured using an analog rotameter; therefore, the process airflow values were recorded manually in the SCADA system. Different probes were used to monitor continuously liquid measurements both in the feed and inside the MABR reactor. NH x concentrations and T were measured using AmmoLyt® plus, and ORP was measured using a SensoLyt® ORP, both from WTW (see Fig. 1H, I, J). A sample from the exhaust gas after the MABR unit was taken semi-continuously to a gas-monitoring GASloq 1200 from ABB. The system contained a gas analyzer Uras 26 Easyline and was designed to operate as a multi scan and measuring point analysis system with a carbon dioxide (CO 2,exh ) and oxygen (O 2,exh ) measurement. Both probes and gas analyzer were connected to the SCADA system, and measurements were recorded every 10 s (see Fig. 1K).  Finally, Fig. 1 also illustrates the MABR's main process characteristics. The hollow-fiber membrane plays two primary roles: act as a biofilm carrier and transfer gas from the lumen to the biofilm through the membrane gas-permeable wall. Substrates are supplied to the biofilm from the biofilm-liquid boundary layer, and the concentration of the substrates decreases with the biofilm thickness toward the biofilm-membrane interface. Oxygen diffuses through the membrane into the biofilm, and the oxygen concentration in the biofilm decreases with the biofilm thickness toward the biofilm-liquid boundary layer. This counter-diffusional nature enables the co-existence of oxygenrich and oxygen-depleted regions within the biofilm, resulting in its stratification into possibly aerobic, anoxic, and anaerobic zones.

Analytical methods
Time-proportional composite samples (24 h) from the feed and the pilot reactor's interior were collected regularly and analyzed throughout the study. The concentrations of NH x , nitrite (NO 2 ), nitrate (NO 3) , phosphate (PO 4 ), and total and dissolved chemical oxygen demand (COD) were measured using NANOCOLOR standard tests and a spectrophotometer NANOCOLOR UV/VIS II (Macherey -Nagel). Mixed liquor suspended solids (MLSS) was determined using standard procedures (SM 20th 2540D).

Oxygen transfers rate (OTR), oxygen transfer efficiency (OTE) and aeration efficiency (AE)
The oxygen transfer rate (OTR) is the measurement of the flux of oxygen gas that diffuses from the lumen's interior into the biofilm over time. The calculation of OTR (g O 2 m −2 d −1 ) was based on the exhaust oxygen model reported in Houweling and Daigger, 2019 (Eq. 1, (2)): where Q air,in is airflow (N m 3 d −1 ), x o2,in is the mole fraction of oxygen in atmospheric air, V loss is the volumetric air loss between inlet and outlet, x o2,out is the mole fraction of oxygen in the exhaust, ρ o2 is the oxygen density under normal conditions (kg m −3 ), and A is the membrane surface area (m 2 ). OTE refers to the ratio of the oxygen transferred, and the oxygen supplied through the membranes. OTE was calculated based on the principles stated in Houweling and Daigger, 2019 (see Eq. (3)): The aeration efficiency (AE) of a given system refers to the mass of oxygen delivered per unit of power required. The AE (kg O 2 kW −1 h −1 ) was calculated according to Metcalf and Eddy (2014) (see Eqs. (4) and (5) P ¼ w RT 28:97 n e p 2 p 1 where P (kW h) is the estimated power requirement to deliver Q air,in , w (kg s −1 ) is the airflow rate, R (−) is the universal gas constant, T is the temperature (C) of the influent gas to the blower, n is (k-1)/k where k is the specific heat ratio, e is the blower efficiency and p 1 and p 2 represent the pressure at the inlet and outlet of the blower, respectively.

Nitrification rates
Nitrification rates (NR) represent the quantity of NH x oxidized to NO x . NR were calculated using four different methods (see Eqs. (6), (7),(8)): NR using NH x concentrations from NH x sensor (NR 1 ) and composite samples (NR 2 ) were described in Eq. (6) where NH x,inf and NH x,eff (g m −3 ) represented the concentration of NH x obtained from laboratory composite samples analysis and online signals and Q inf was the influent flow rate (m 3 h −1 ).
NR based on oxygen utilization (NR 3 ) assumed that 100% of the oxygen transfer through the membranes (g O 2 m −2 d −1 ) was used for the conversion of NH x into NO 3 and it was calculated as described in Eq. (7), where 4.57 is the theoretical oxygen demand for nitrification (g O 2 /g NH x -N) (Henze et al., 2008).
NR based on batch tests (NR 4 ) was calculated from batch tests performed in the pilot reactor as follows: feed flow and EA were stopped while the MABR unit continued in operation. The concentration of NH x in the tank was continuously monitored through the online sensor signal. The test was stopped when the NH x concentration reached 5 g N m −3 . The rate was estimated from the slope of the linear decrease in the NH x concentration (ΔNH x,eff ) over time (Δt) after the feed flow to the reactor was turned off. In Eq. (5), V was the volume of the reactor (m 3 ). More details are provided in Eq. (8).
2.6. Statistical analysis All analyses were performed using R software. A Shapiro-Wilk statistical test was applied to check the assumption of normality in the distribution of the data. When the test showed that normality could be assumed, parametric tests such as the Pearson correlation test and t-test were used. On the other hand, if the sample showed a nonnormal distribution, then non-parametric tests such as Kendall rank correlation tests and Wilcoxon signed-rank tests were selected instead. Logarithmic transformations were not explored. Correlation factors were considered "strong", "moderate" and "weak" corresponding to absolute loading values >0.5, between 0.5 and 0.3 and < 0.3. Results were considered "significant" when the p-value was <0.05 and "highly significant" when the p-value <0.001.

Preliminary data analysis
The preliminary data analysis included the treatment of missing values and removal of data when the system was temporarily out of service for different reasons: maintenance of equipment or emptying of the reactor tank. The data was aggregated into daily averages before applying the following multivariate techniques.

Multivariate techniques
A multivariate analysis of the sensors/gas analyzers data was performed using principal component analysis (PCA) and cluster analysis (CA). PCA extracted the eigenvalues and eigenvectors from the covariance matrix of the scaled and centered variables. The PCs are the uncorrelated (orthogonal) variables obtained by multiplying the original correlated variables with the eigenvectors. Each eigenvector consists of a vector of coefficients (loadings). PCA allowed reducing the dimensionality of the original data set with a minimum loss of information. The selection of the PCs was based on eigenvalues >1. A Promax rotation was performed on the PC eigenvalues to achieve a simpler, easier to interpret structure. CA is an unsupervised pattern recognition technique that classifies the objects or a system into categories or clusters based on the nearness or similarity of data points. In this paper, kmeans clustering based on Euclidian distance was performed on the data set. We computed 30 different methods to choose the appropriate number of clusters and chose the one suggested by most methods.

Oxygen transfer
The MABR unit was operated at 8 m 3 h −1 process airflow during most of the study period, with shorter periods at 6, 10, and 12 m 3 h −1 . The outlet gas was manually set to 0.25 bar gauge pressure to allow for condensate removal since water must overcome the hydrostatic pressure to leave the system and to prevent the intrusion of water in the event of leaks. During the study period, the average discharge pressure was 0.26 ± 0.04 bar, while the average inlet pressure was 0.45 ± 0.04 bar. The pressure loss across the system can vary depending on the airflow applied, condensate build-up, or occurrence of leaks. Several manual air mass flow measurements before and after the membranes were carried out to confirm Eq. (2). The average volume loss observed was 3%, while the average volume loss calculated as in Eq. (2) was 5%. Fig. 2A shows OTR and OTE under the described operational conditions, calculated using Eqs. (1), (2) & (3). The entire period's average values were 9.54 ± 1.79 g O 2 m −2 d −1 (Eq. (1)) and 25.7 ± 5.75% (Eq. (3)), respectively. The AE averaged 5.79 ± 1.42 kg O 2 kW −1 h −1 (see Eq. (4)). Fig. 2A also reveals that different process airflows changed the relationship between OTR and OTE. More specifically, lower airflows achieved higher OTE (up to 47%) at the same values of OTR. This means a higher quantity of oxygen is delivered into the biofilm per unit of volumetric airflow, which has a strong implication from an energy consumption perspective (see Eq. (5)). On the other hand, at higher airflows, the OTE was lower comparatively and energy consumption higher (see Eq. (4)).
Finally, Fig. 2B reveals that OTR has a strong correlation to the NH x concentration in the MABR tank (Kendall's rank correlation τ = 0.56, p-value <2•10 −16 ). At higher NH x concentrations, an increased OTR was achieved (and vice versa). Fig. 2B also shows that OTR values can be accurately predicted from NH x concentrations using a linear regression model (R 2 = 0.55, p-value <2•10 −16 ). This relationship is a unique property of counter-diffusional biofilms, where microbial activity in the biofilm acts as a catalyzer for oxygen transfer. Therefore, at higher NH x concentrations, more nitrification activity occurs in the biofilm, and more oxygen is transferred.
OTR values were not correlated to process airflow rates (Kendall's tank correlation τ = 0.006, p-value <2•10 −16 ). However, when we combine both NH x,eff , process airflows, OTE and pressure using a multiple linear regression model (R 2 = 0.95, p-value <2•10 −16 ), we observe that for a given pressure and NH x concentration, an increase in a unit (m 3 h −1 ) of process airflow will result in an increase in OTR of 0.29 g O 2 m −2 d −1 (p-value = 3•10 −4 ). And when predicting OTE, for a given NH x,eff , pressure and OTR (R 2 = 0.96, p-value <2•10 −16 ), decreasing process airflows in one unit results in an increase in OTE of 3% (p-value <2•10 −16 ). Therefore, operating at higher NH x , eff and low process airflows will result in both high OTR and OTE, increasing the overall energy efficiency of MABRs, but higher absolute OTRs will be achieved at higher process airflows at a given NH x concentration.

Startup
During the initial 6 weeks of operation, the startup phase, feed flow into the reactor was slowly increased starting from 7 m 3 h −1 to 15 m 3 h −1 , corresponding to a NH x,l,inf of 1.1-4.88 g N m −2 d −1 , depending on NH x,inf concentrations. The feed to the demonstration reactor, activated sludge from EMWRRF, is expected to contain nitrifying organisms, and their attachment to the biofilm could help speedup the development of a nitrifying biofilm (Martin and Nerenberg, 2012). The influent and effluent NH x sensors wrongly indicated that nitrification was taking place in the reactor from the beginning, due to sensor drift, which was detected when the first grab sample was taken, 2 weeks into operation. The first composite sample analysis, carried out approximately 6 weeks after startup, showed a nitrification rate of 1.6 g N m −2 d −1 .
It is generally well-known that microbial activity in MABR biofilms increases oxygen diffusion across the membrane lumen Nerenberg, 2016) and that bulk NH x concentrations are generally the limiting substrate for nitrifying counter-diffusional biofilms treating domestic wastewater (Houweling and Daigger, 2019). Moreover, the higher oxygen utilization to growth ratio of nitrifying versus heterotrophic organisms (19:1 and 0.5:1), indicates that even low fractions of nitrifying organisms can be responsible for high rates of oxygen utilization in the biofilm (Houweling and Daigger, 2019). In this study, we use the relationship between oxygen transferred and NH x concentrations (see Fig. 2b) as an indicator for nitrifying activity in the biofilm. Fig. 3 shows average hourly data from the O 2,exh, and the NH x sensors in the MABR reactor for 6 days during the first 6 weeks of operation. Measured data shows a lack of response from the MABR to changing NH x concentrations during the first 2 weeks after installation, as indicated by the graph's flat slopes. At week 3, a change in the slope could be seen, indicating the presence of nitrifying organisms in the biofilm. From week four on, the slope is more accentuated, indicating a strong presence of nitrifying organisms in the biofilm, which drives the higher consumption of oxygen at higher NH x concentrations. From the generated results, we conclude that OTR in the MABR unit varies along with NH x concentrations, and it took approximately 3 weeks to establish a nitrifying biofilm. The first point was already discussed in Section 3.1. Fig. 4 depicts box plots of the NR values quantified using the previously described four different methods: online sensors (NR 1 ), composite samples (NR 2 ), oxygen transfer (NR 3 ), and batch mode (NR 4 ). The results show that NR 1 , 3 had a larger number of samples than NR 2,4 ; while NR 3 , 4 had a lower variability, with median absolute deviations of 0.39 and 0.25 respectively, while NR 1 , 2 had the highest variabilities (median absolute deviations are 0.49 and 0.46). Calculations reveal significant differences between the estimated NR values among NR 3 and NR 1 , 2 , 4 (all three p-values <0.001). This can be attributed to an overestimation of NR 3 based on the assumptions of Eq. (5), where all the oxygen transferred is used only for nitrification. Hence, the differences with NR 1 , 2 , 4 methods indicate that a significant amount of oxygen is used for other purposes, i.e., removal of COD demanding compounds such as organics and sulfides. Indeed, in systems where there is a high amount of (aerobic) COD removal within the biofilm, the use of NR 3 should be revised. Another factor affecting the other three NR methods (NR 1 , 2 , 4 ) is hydrolysis in the mixed liquor surrounding the membranes. The reactor is   Table 1); it is expected that some of the particulate compounds present will be solubilized despite the low hydraulic residence time of the reactor. These soluble organic compounds can contribute to the biofilm's oxygen consumption by creating an additional oxygen requirement not captured by the NH x sensors or the influent/effluent analysis. Therefore, the real NR value lies in between NR 1 , 2 , 4 and NR 3 .

Multivariate analysis of sensor data
PCA revealed that the first component was able to explain 49.9% of the variability of the data, and the second and third 19.8 and 11.7%, respectively. An analysis of the eigenvalues (Fig. S1) confirmed that the appropriate number of PCs was three. Both Table 2 and Fig. 5A show that PC1 has a strong correlation to NH x loading (NH x,l,inf ) and moderate loadings on O 2,exh , CO 2,exh , and OTR, all negative except O 2,exh . In PC2, the EA and Q inf are strongly correlated in different directions, with moderate contributions from T (positive) and NH x,l,inf, and ORP inf (negative). In the last PC (PC3), both NR 1 and EA were strongly correlated, and ORP inf , ORP eff had moderate loadings, all in the same direction.
During CA, we computed 30 different methods using Euclidian distances and a minimum of two centroids and a maximum of six. We found that 14 of the methods suggested 5 were the appropriate number of clusters to adequately represent the data (Fig. S2). The PCA biplot with the projected scores and the loadings, grouped into the five clusters obtained with the CA k-means algorithm, can be seen in Fig. 5A and B. Cluster IV and V mostly overlap chronologically (Fig. 5B), indicating they represent a period of time with episodes characterized by different conditions in the variables that govern PC2: EA, Q inf , T and NH x,l,inf (see Table 3 and Fig. 5A). Cluster III and V were on opposite sides of PC1, despite corresponding to periods close to each other in time, showing significant differences in the variables that govern PC1, such as NH x,l,inf, and the response from the exhaust gas side: OTR, CO 2,ex and O 2,ex (see Table 3 and Fig. 5A).
A Wilcoxon Signed-Rank test was included in determining each variable's importance within the cluster with respect to the total sample. Table 3 contains the calculated median values for the selected variables used during the multivariate analysis for each of the identified clusters, including a symbol indicating whether the median of the variable for each cluster was significantly different from the entire sample's median. Cluster numbers II, IV and V represent the more extreme conditions in the data set, containing respectively 10, 11 and 9 out of 11 variables that differ significantly from the total data set. On the other hand, clusters I and III correspond to periods of time closer to average conditions with 6 and 5 significantly different variables.
Cluster I is characterized by having the highest EA contribution from all periods considered (see median values in Table 3). In this case, the increased use of fine bubble aeration was not a response to low ORP because during this period, the EA was controlled using a time-based control, independent from ORP (see Section 2.2). Therefore, the ORP in the tank was significantly higher than for the rest of the study (see Fig. 6E and F). However, it is worth noting that the ORP eff was still within an anoxic range, not aerobic (<0 mV). The MABR during this time achieved high OTR and the highest NR, with an OTR/NR ratio of 6.18 (see Fig. 6G and H), the closest to the theoretical value of 4.57.
Cluster II had the lowest T, NH x,eff and NH x,l,inf and relatively high Q inf , (see average values in Table 3). The ORP for both the influent and the effluent was the highest. Correspondingly, the use of EA was the lowest. This period can be considered representative of winter conditions in the Nordic region in Europe. During this period, the response of the MABR was stable nitrification (1.18 g N m −2 d −1 ) (not significantly different from the median of the sample) and relatively low OTR. The low OTR results are in line with the observation made in Section 3.1, where we showed OTR values were correlated to NH x,eff (see Fig. 5). Lastly, a high percentage of the oxygen transferred was used for nitrification; the OTR/NR ratio was 6.38 (see Fig. 6E and F).
In cluster III, T was low, and there was almost no use of EA, but unlike cluster II, it also had the highest Q inf and high NH x,l,inf into the MABR reactor (see average values in Table 3). This period would correspond to typical wet-weather conditions in the spring, where treatment plants  can receive high hydraulic loadings (see Fig. 6A). NR was very similar to those in cluster II (1.20 g N m −2 d −1 ), while OTR was higher, with an OTR/NR ratio of 8.12. This is in accordance with what the analysis of PC1 indicated, i.e., OTR was strongly correlated to NH x,eff, and what was previously shown in Section 2.1 (see Fig. 2B). Clusters IV and V are characterized by the highest T, high NH x,eff , and the lowest values for ORP both in the influent and the effluent, which caused a high use of EA (ORP-feedback control) (see average values in Table 3). This period represents summer conditions in Northern Europe (see Fig. 6A,C & F). However, the extremely low ORP inf values (<−300 mV) are specific to this facility and the operating conditions of the EBPR zones. Clusters IV and V mainly differ in the NH x,l,inf , which was almost double in cluster IV than in cluster V and in the ORP inf , which was lowest in cluster V (see Table 3 and Fig. 6D, F). Fig. 6B shows how the hydraulic retention time (HRT) in cluster IV and V exhibits a different behavior than in the rest of the clusters. This is caused by the ORP-based feedback control, which not only increased EA but also reduced Q inf when the ORP reached a low setpoint. Cluster V's extremely low ORP inf caused the reduction in NH x,l,inf . NR 1 in both clusters was the lowest of all periods, although higher in cluster V than in cluster IV, suggesting a further increase in HRT could have improved NR. There was also a very high OTR/NR ratio, indicating a very significant decoupling of these two indicators. The high concentration of CO 2,ex during these periods could be an indicator for increased heterotrophic activity in the biofilm: while nitrifiers in the biofilm consume CO 2 , heterotrophic organisms produce it, a higher concentration in the exhaust indicates a shift in the balance between those two populations of microorganisms, toward a higher heterotrophic activity.

Analysis of composite samples
Over the study period, 16 composite samples were collected and analyzed. There was a data gap during late winter 2020 and spring 2020 due to restrictions related to COVID19 that impeded the acquisition of samples (see Fig. S3). As mentioned in Section 2.2, the reactor was fed with sludge from the EBPR zones at EMWRRF, where biological phosphorus release is expected to occur, and the main characteristics are summarized in Table 1. Effluent COD, NH x and PO 4 concentrations were on average: 36.0 ± 7.8 g COD m −3 , 8.77 ± 4.98 g N m −3 and 6.55 ± 6.22 g P m −3 . NO 2 and NO 3 concentrations in the reactor averaged 0.02 ± 0.01 and 0.36 ± 0.19 g N m −3 , respectively. More information about the dynamic profiles can be found in the Supplementary Information section.
A correlation matrix (see Fig. 7) shows that NR 2 had a strong positive correlation with NO 2 and NO 3 concentrations and a negative correlation with PO 4 effluent concentrations. Neither total nor dissolved COD concentrations in the influent showed a significant correlation to NR 2 . The positive correlation between NH x oxidized and NO X produced was expected. However, the negative correlation between NR 2 and PO 4 concentrations was not. This could be attributed to the fact that the higher PO 4 concentration in an EBPR zone would be an indication of less readily available COD, which should mean less competition for oxygen between nitrifiers and ordinary heterotrophic organisms. In this case, PO 4 concentrations were positively correlated to dissolved COD concentrations. This could be explained by the fact that the EMWRRF not only relies on EBPR to remove phosphorus but also uses chemical precipitation in the primary settling tanks (see Fig. 1). At higher T, dissolved COD produced by hydrolysis in the EBPR zones or upstream could be present in excess, related to the amount of PO 4 available after chemical precipitation. MLSS concentrations were negatively correlated to dissolved COD, NH x,inf/eff and T; because at higher T, the EMWRRF is operated at longer solids retention times, increasing the MLSS concentration in the EBPR zone. The detailed mechanisms involved in the observed denitrification and phosphorus removal in the MABR reactor are currently under investigation and will be reported at a later stage.

Performance results
Although small pilot and lab-scale studies have shown much higher OTE (up to 69% (Li and Zhang, 2018)) and OTR (up to 15 g O 2 m −2 d −1 (Côté et al., 2015)); the results from this study (OTE = 25.7% and OTR = 9.54 g O 2 m −2 d −1 ), are in line with recently reported experiences with full-scale MABR units, with OTR values ranging between 10 and 12 in MABR with surface areas comprising 2200 and 1920 m 2 (Koch et al., 2019;Houweling et al., 2017). Process airflow rates showed no correlation to OTR when analyzed independently, but showed a significant correlation when combined with NH x concentrations. Recently reported results from laboratory-scale studies (Bunse et al., 2020), showed airflow had a significant impact on OTR especially at low operating pressures. In our regression model, pressure was not statistically signifcant, this could be caused by the higher pressure required when operating full-scale MABR, in order to prevent water intrusion in case of leaks and/or improve condensate removal.
Regarding AE, results up to 19.6 kg O 2 kW −1 h −1 have been estimated for laboratory-scale MABRs, which meant a reduction in energy consumption of 83.7% compared to conventional fine bubble diffusers (3.19 kg O 2 kW −1 h −1 at 27.89% OTE) (Castrillo et al., 2019). The AE obtained in this study, 5.8 kg O 2 kW −1 h −1 , compared to a fine pore diffusion system (up to 2.6 kg O 2 kW −1 h −1 (Wagner and Stenstrom, 2014)), would be 55% more efficient. Moreover, the use of MABR would suppose an average reduction of the energy used for aeration of 74% compared to the existing surface aerators at EMWRRF (up to 1.5 kg O 2 kW −1 h −1 ) (Andreasen, 2010). The use of more advanced aeration strategies, such as periodic venting or intermittent aeration, could further improve the oxygen transfer performance (Perez-Calleja et al., 2017) and even promote shortcut nitrogen removal (Mehrabi et al., 2020).
The MLSS concentrations in the feed and the reactor, which averaged 5892 ± 938 g SS m −3 , are the highest reported MLSS for a hybrid pilot or full-scale MABR. Hydrodynamics can negatively affect the performance of the MABR, influencing the thickness of the concentration boundary layer at the biofilm-liquid interface, among others (Syron and Casey, 2008). Mixing was maintained as high as the dedicated blower would allow, in this case, 12 m 3 h −1 out of concern that MLSS concentrations above 5000 mg/L would affect the hydrodynamics of the system. Table 3 Summary of performance data grouped by cluster number. Statistical significance of the distributions per cluster compared to the overall distribution using a Wilcoxon signed-rank test are indicated in bold (p-value <0.001) and italics (p-value <0.05). Statistical significance of the distributions per cluster compared to the overall distribution using a Wilcoxon signed-rank test are indicated in bold (p-value < 0.001) and italics (p-value <0.05).
Under conventional application, process airflows are typically recycled and used for mixing/scouring. Under these conditions, the need for mixing/scouring might be the determining factor for choosing a process airflow value. When mixing/scouring and process aeration are decoupled (like in this case study), operating at lower airflows could benefit energy consumption while achieving high OTR, as shown in Fig. 2A. NR values in this study (1.03-2.14 g N m −2 d −1 ) were in the same order of magnitude of previously published studies. Downing and Nerenberg, 2008 reported a NR value of 1 g N m −2 d −1 in a laboratory-scale hybrid MABR fed with synthetic wastewater, while Koch et al., 2019 observed up to 1.3 g N m −2 d −1 in a full-scale MABR pilot fed with primary effluent. Houweling et al., 2017 reported values from four pilot studies that achieved up to 2.6 g N m −2 d −1 , using the method stated in Eq. (7), close to the median value found in this study for NR 3 of 2.14 g N m −2 d −1 . As expressed before in section 3.3, since hydrolysis was not accounted for, causing the underestimation of NR 1,2,4 , and N 3 did not take into account the use of oxygen for purposes other than nitrification, this study's real NR value lies somewhere in between the results obtained by NR 1,2,4 and N 3 (1.03-2.14 g N m −2 d −1 ). Therefore, special consideration should be taken when estimating NR using OTR as in NR 3 in systems with high organic loading; and using NH x concentrations as in NR 1,2,4 where hydrolysis is expected to occur. It is important to highlight that the volumetric NR calculated in this study correspond to 175-364 g N m −3 d −1 which are comparable to conventional intensification biofilm-based technologies well established in the market, such as the moving bed biofilm reactor (MBBR) (100-400 g N m −3 d −1 (Metcalf and Eddy, 2014)). However, the upperend NR in MBBRs are obtained at high dissolved oxygen concentrations in the bulk (up to 6 g O 2 m −3 ), while the NR found in this study occurred in anoxic or even anaerobic bulk conditions due to the counter-diffusional nature of the biofilm (see Fig. 1). This allowed the occurrence of simultaneous nitrification and denitrification in this study, with NO 3 concentrations in the reactor below 1 g N m −3 .
The central findings obtained in this study, i.e., energy-saving and footprint reduction should be analyzed in-depth considering other factors such as CAPEX and OPEX when intensifying existing activated sludge systems. For that purpose, process simulation could be an excellent choice to run a conceptual design assessment, and there are already available mechanistic models of an MABR that could be used for that purpose (Lackner et al., 2008;Ni et al., 2013;Acevedo Alonso and Lackner, 2019). In this way, it would be possible to assess aspects such as membrane area (and cost) versus degree of nutrient removal and energy consumption. Moreover, the relationship between OTR and NH x concentrations and how that would affect a system's design, including several MABR, could be further explored. This is a fascinating investigation but out of the scope of this study.

Factors affecting performance
From the multivariate analysis of the long-term sensor data (see Table 3, Fig. 5), we can conclude that NR in this study was not significantly affected by winter or spring conditions: low T, high hydraulic loading, and low NH x concentrations. This positions MABR as a promising process intensification technology to be considered under Nordic conditions. In contrast, Houweling et al., 2017 found that rain events or snowmelt may have a significant negative response on the oxygen transferred as measured by OTR and NR in a hybrid MABR. Even during winter and spring, the high MLSS and low ORP from this study could have acted as a "buffer", preventing aerobic conditions in the bulk. These differences could be attributed to the different characteristics of the sludge surrounding the MABR in both studies (one anaerobic and the other anoxic) or to the biofilm thickness (See Fig. 1). T sensitivity in conventional co-diffusional biofilm-based technologies, such as MBBR, has been extensively reported (Regmi et al., 2011), but T did not impact NR in this study significantly. Moreover, NH x,l,inf was only correlated to the gas side's responsehigher NH x loadings lead to higher OTR and vice versa-but not to NR.
The fact that NR 2 showed no correlation to COD concentrations in this study is in disagreement with previously reported results: LaPara et al., 2006 found heterotrophs that outcompeted nitrifying autotrophic organisms in a laboratory-scale MABR fed with sodium acetate and ammonium chloride at high COD concentrations. Houweling et al., 2017 found NR were affected by NH x load and C/N ratio in a pilot-scale MABR. However, Downing and Nerenberg, 2008 showed that in a hybrid configurationlike the one in this study-, a 1920 cm 2 MABR fed with synthetic wastewater, NR were insensitive to BOD loadings. Moreover, the fact that NR 2 was strongly and negatively correlated to PO 4 concentrations in the influent is also in disagreement with the general understanding that the higher the PO 4 concentration in EBPR zones, the lower the soluble COD available could compete for oxygen with nitrification. This phenomenon is currently under investigation.
The difference between NR 3 and the other three methods presented, NR 1,2,4 can be explained by the high OTR/NR ratio found in this study. As we could see in Fig. 1, ideally, the biofilm's aerobic layer would be entirely comprised of autotrophic organisms oxidizing NH x ; however, this will not always be the case in practice. Kunetz et al., 2016 also observed in a hybrid MABR that at higher NR, a higher fraction of O 2 was used for nitrification. The influence of the biofilm thickness and different layer composition on the MABR was out of this work's scope but will be further explored in the future.

Methods for assessment of performance
The volume loss across the membranes, and therefore OTR, OTE, and AE was calculated using Eq. 2, which only considers oxygen diffusion in the volume loss across the membranes. However, this method proved to be accurate enough in this study (average 5% versus 3% obtained with manual mass flow measurements). Therefore, to estimate OTR, it seems unnecessary to use mass flow measurements before and after the MABR units; airflow measurements in the inlet and O 2,exh concentration readings suffice, significantly reducing the instrumentation required.
The relationship between OTR and NH x concentration reported in this study had been previously observed Houweling and Daigger, 2019). Pellicer-Nàcher et al., 2013 observed ammonia concentrations increased biofilm activity and therefore OTR, and suggested that oxygen transfer could be further exploited in MABRs when operating under excess NH x concentration and limiting O 2 . Houweling et al., 2020 showed this relationship, referred to as Natural Ammonia-Based Aeration Control, could be exploited to use gas measurements of oxygen as a surrogate for NH x concentrations in the bulk liquid. Results from this study, presented in Fig. 2B, support that statement. Fig. 3 shows it can also be used to determine startup times or even system upsets affecting nitrifying populations, especially when NH x are unreliable or unavailable.
This study presents CO 2,exh concentrations in the exhaust of an MABR reactor for the first time. From the analysis of PC1 (Fig. 5 and Table 2), we can derive that there is a strong correlation between NH x,l,inf and the response from the exhaust gas side (O 2,exh , CO 2,exh ,) and OTR. Indeed, at higher OTR, the O 2,exh is lower, and CO 2,exh is higher. We can also see in Table 3 that the CO 2,exh concentration was highest when the OTR/NR ratio was highest in cluster IV and V. This forms a strong indication that CO 2,exh measurements in the exhaust, could provide useful additional information on whether the balance between autotrophic or heterotrophic activity in the biofilm shifts over time.
Besides, the accumulation of CO 2,exh in the exhaust of the MABR, together with other gases such as nitrous oxide (not presented in this study), opens up possibilities for downstream processing of these airflows, as it is commonly done with other streams of exhaust emissions (Sharif et al., 2021). The results from NR 4 confirm that NR in the reactor can be entirely attributed to oxygen transferred from the MABR since batch tests were performed without EA. Batch tests proved to be a simple way to accurately monitor NR using one single sensor, reducing the accumulated error caused by two signals (influent and effluent) and the equipment and staff requirements to perform composite samples analysis. This could prove useful in full-scale installations where it is possible to operate in batch mode for a few hours every week, such as installations with reactors operated in parallel or during low load situations (i.e., night time), and it could be easily automated and integrated in the SCADA system.
The use of multivariate techniques during this study allowed the assessment of how the many different variables involved in the process affected the MABR reactor's performance. Some of these variables were manipulated: Q inf or EA, while others could not be: ORP inf, or T. The use of these advanced statistical techniques, while still not widespread, has already proven useful when investigating the complex processes governing wastewater treatment. Vasilaki et al., 2018 andBellandi et al., 2020 used multivariate analysis to identify the dependencies and underlying patterns between nitrous oxide emissions and operational variables.
Last, microbial analysis techniques could help elucidate some of the mechanisms behind the observations made in this study. Techniques such as FISH (Schramm et al., 2000;Terada et al., 2003;Matsumoto et al., 2007;Downing and Nerenberg, 2008) and PCR (Terada et al., 2010;Taşkan et al., 2019;Zhao et al., 2021) have been successfully used to identify key organisms in laboratory-scale MABR biofilms. However, it is worth noting the difficulties associated with biofilm sampling in full-scale installations, where a crane is required on-site to lift-up the MABR from the tank to acquire the samples. Very few samples were collected during the study period and are currently under investigation.

Conclusions
The main findings of the presented study can be summarized in the following points: • NH x concentrations governed the oxygen transfer processes. No relevant correlation was found between process aeration and OTR (τ = 0.06, p < 2•10 −16 ), and OTR was strongly correlated to NH x concentration up to the point that could be accurately predicted using a linear model (R 2 = 0.55, p < 2•10 −16 ). Also, lower airflow values achieved higher OTE (up to 47%) at the same range of OTR. • The relationship between oxygen transferred and NH x concentration was used to identify the startup time necessary to establish a nitrifying biofilm, which was about 3 weeks in the case of this study. • NR calculations revealed significant differences between oxygenbased methods (NR 3 ) and methods based on NH x measurements (NR 1 , 2 , 4 ) (all three p-values <0.001). It is possible that under high organic loading, N 3 will underestimate NR, and when hydrolysis of the MLSS surrounding the membranes occurs, N 1,2,4 underestimate NR. Special consideration should be taken under these conditions. • Three PCs were able to explain 81% of the data variability, and the five clusters obtained by CA allowed to identify distinct periods of operation from online data (sensors/analyzers). Besides, interesting process dependencies could be found such as NR remained stable under low-temperature periods corresponding to spring/ winter Nordic conditions. Unexpectedly, nitrification rates were lowest during the summer, possibly due to the low ORP in the feed (< −300 mV). • The Correlation matrix with Pearsons' coefficients for composite samples (laboratory analysis) revealed some expected correlations, for example, NR and NO 2 and NO 3 concentrations. Nevertheless, it was difficult to explain the trade-offs between nitrification rates and effluent phosphate, and additional investigations will be necessary. • Significant energy savings could be obtained with MABR, where the average AE obtained was 5.8 kg O 2 kW −1 h −1 , which would mean an average reduction in energy consumption of 74% compared to the existing aeration at the facility. With respect to footprint, the volumetric NR values reported in this study (175-364 g N m −3 d −1 ) are comparable to other biofilm-based technologies in the market. • This study demonstrates MABR is a suitable technology to aid facilities in the Nordic region (or similar weather conditions) to intensify existing activated sludge-based biological nutrient removal achieving significant energy savings simultaneously. The mechanisms behind some of the results presented and how the biofilm thickness and microbial stratification played a role in these observations are not yet fully understood and further experiments are required.
Software availability