The study focused on six public wastewater treatment plants in Gyeongsangbuk-do, Republic of Korea. Monthly sampling was conducted from June to December 2021, with integrated samples collected once a month from each plant. The process characteristics are provided in Table 1. Samples were collected in the following sequence: influent (including associated treatment water), after primary sedimentation, anaerobic stage, after secondary sedimentation, and effluent. The collected samples were stored in ice-filled coolers and transported to the laboratory. Field investigation parameters included temperature, pH, electrical conductivity, and salinity, measured using a YSI Professional Plus device. Four-liter samples were collected at each point for physicochemical analysis.

The following water quality parameters were analyzed according to the Water Pollution Standards Method [22]. Chemical Oxygen Demand (COD_Mn), Total Organic Carbon (TOC, NPOC method), Biological Oxygen Demand (BOD₅), Total Nitrogen (TN), Total Phosphorus (TP), Mixed Liquor Suspended Solids (MLSS), Ammonia Nitrogen (NH₄⁺-N), Nitrate Nitrogen (NO₃⁻-N), and Nitrite Nitrogen (NO₂⁻-N). The analysis equipment used included the JP/TOC-LCPH (Shimadzu) for TOC, the SAN++ (Skalar) for TN and TP, and the CH/850 Professional IC (Metrohm) for ion parameters (NH₄⁺-N, NO₃⁻N, and NO₂⁻N).

Microbial samples for analysis were collected in July from six WWTPs. For each sample, 100 mL was taken from the aerobic stage, and samples were transported to the laboratory in 4°C. Upon arrival, each sample was divided into 25 mL aliquots in four conical tubes. After immediate centrifugation, the supernatant was discarded, and the pellets were preserved at −20°C until further analysis. All WWTPs followed the same sampling protocol. For the AD facility, additional sampling was conducted on the media. Microbial material was scraped from the spherical plastic media (20 cm in diameter) containing plastic sponges (35 cm ’ 13 cm ’ 5 cm) and was washed with phosphate buffered saline solution (Thermo Fisher Scientific) to removal residual plastic. After centrifugation at 3000 rpm, the supernatant was discarded, and the remaining microbial pellets was stored for further analysis.

The analysis of microbial composition was conducted on the samples using 16S rRNA gene sequencing. Genomic DNA was extracted using the FastDNA^TM SPIN kit for soil (MP Biomedicals). Sequencing libraries were prepared following the Illumina 16S Metagenomic Sequencing Library protocols to amplify the V3 and V4 region. Genomic DNA (2 ng) was subjected to PCR amplification using Herculase II fusion DNA polymerase Nextera XT Index Kit V2 (Agilent Technologies). The target range for the V3–V4 region was amplified using the primer set 341F (CCTACGGGNGGCWGCAG) and 805R (GACTACHVGGGTATCTAATCC) [10, 23]. The final product was quantified using qPCR and qualified with TapeStation D1000 ScreenTape. Sequencing was performed on the MiSeq^TM platform (Illumina). The Illumina sequencer generates raw images utilizing sequencing control software for system control and base calling through an integrated primary analysis software called RTA (Real Time Analysis). The BCL (base calls) binary files are converted into FASTQ files using the Illumina package bcl2fastq. Data analysis was conducted using the following tools: Cutadapt for adapter and primer trimming, and DADA2 for quality filtering, denoising, merging, and chimera removal. These tools ensured accuracy and reliability in the sequencing data processing [24].

The correlation analysis was conducted using the R program (version 4.3), employing Pearson correction coefficients, which are commonly used to measure the strength of relationships between variables. The analysis included data on physicochemical factors (e.g., pH, temperature, salinity, and conductivity) and the read abundance ratios of PAOs and GAOs. A coefficient of +1.0 indicates a perfect positive correlation, 0 indicates no correlation, and −1.0 indicates a perfect negative correlation. A correlation range of 0.9–1.0 indicates a very strong correlation, while a range of 0.7–0.89 suggests a strong correlation [25].

The facilities under study employ various modified forms of the activated sludge method for wastewater treatment, and their specific characteristics were monitored monthly over seven months, as summarized in Table 1. The geographical locations of the six WWTPs are as follows: YC, MG, US, and AD are located inland, at distances ranging from 65 km to 109 km from the coast. In contrast, YD and YH are situated within 2 km of the coast, directly adjacent to the sea. To explore the relationship between geographical location and the characteristics of influent and microbial communities in the WWTPs. Table 1 provides an overview of influent characteristics and process parameters for each facility. This includes geographic location, process type, daily treatment capacity, influent composition, solids retention time (SRT), temperature, salinity, and conductivity. The influent sources include not only domestic sewage but also livestock manure treatment water and treated food wastewater, which are mixed in the flow equalization tank before biological treatment stage. Notably, conductivity values were high in the inland facilities, such as YC and MG, despite their distance from the coast. These various operational parameters significantly influence microbial communities, making it important to consider these values for a more accurate interpretation of the treatment processes. With the exception of AD, all WWTPs in this study operate using modified versions of the conventional activated sludge process. The specific process types used in each plant are also listed in Table 1. The treatment efficiency in Table 1 was calculated based on the influent and effluent BOD values, as derived from sewage statistics provided by the Korea Environment Corporation. Except for YH (98.3%), all facilities demonstrated treatment efficiencies exceeding 99.0% for the final effluent. The DePinho process in AD is the only process that included media.

The operating temperature of the bioreactors averaged 21.4 ± 4.3°C, which is within the typical rage. This deviation was not significantly different from the results of a previous study that analyzed samples from 18 wastewater treatment plants, where temperatures ranged from 9°C in winter to 18°C in summer [10]. In terms of pH, it was generally within the neutral range, with an overall pH of 7.4 ± 0.28 (Table 1). Salinity, conductivity, and chloride ions in the influent were not correlated with geographical characteristics. Four WWTPs were located inland (YC, MG, US, and AD), and two were situated on the coast (YD, and YH). YD, located on the coast, includes various other effluents, such as those from seafood processing industries, in its influent. However, in the case of YH, despite its coastal location, all values of salinity, conductivity, and chloride ions were low. These parameters were more significantly influenced by the characteristics of the influent than by geographic characteristics.

To assess the overall phosphorus removal efficiency in the process, both biological and chemical phosphorus removal were investigated (Table 2). Except for YH, 30% of the influent TP was removed in the primary sedimentation tanks of all WWTPs. The biological phosphorus removal was calculated based on the difference in TP values measured after the primary sedimentation and before the secondary sedimentation (during the activated sludge treatment process). In contrast, the chemical phosphorus removal was calculated from the difference in TP values measured before and after the chemical treatment process for phosphorus removal. The removal efficiency of biological phosphorus did not correlate with the removal quantity. In the case of YD, although it had the highest removal quantity of TP (5.804 ± 2.588 mg/L), the removal efficiency was the lowest at 48.0 ± 14.4%. Chemical phosphorus treatment was calculated based on the difference in TP concentration between the effluent from the secondary sedimentation tank and the final effluent, with YD (0.802 ± 0.205 mg/L) showing the highest removal quantity. The highest dependency on chemical phosphorus removal was observed at YH (9.6 ± 5.5%). At the YH facility, on-site operators faced challenges in managing the biological processes. This was further supported by microbiological analyses, which showed a relatively high presence of non-essential microorganisms in the biological reactor compared to other regions.

The COD:N:P ratio of influent is a key parameter in the operation of WWTPs (Table 3). The recommended ratio for biological reactors, as reported in the literature, is 100:5:1 [26]. According to studies on WWTPs in Portugal and Denmark, when P was set to 1, the COD_Cr ranged from 40 to 206, and N ranged from 4 to 20, showing a wider range compared to our study. These six treatment plants were found to treat nitrogen and phosphorus with lower carbon sources compared to WWTPs in Portugal and Denmark [27]. This study utilized the COD_Mn method, which complicates direct comparisons with other research due to its tendency to underestimate organic carbon levels in municipal and industrial wastewater. While COD_Mn typically oxidizes only 30–60% of organic surrogates, both COD_Cr and TOC can oxidize 80–100% of organic compounds [28]. There has been a call in South Korea to replace COD_Mn with COD_Cr or TOC due to the low oxidation efficiency of the permanganate method. In 2022, South Korea addressed this issue by legally transitioning its industrial wastewater discharge regulations to TOC, leading to the inclusion of the TOC:N:P ratio in Table 2. The TOC/P ratio was highest at YC with a value of 26, followed by US at 18. Other locations had ratios ranging from 13 to 14 (Table 3). During the study period, the YC wastewater treatment plant also exhibited the highest biological phosphorus removal efficiency at 88.8% (Table 2). In contrast, YD showed the lowest TOC/P ratio of 9 and the lowest biological phosphorus removal efficiency at 48%. This indicates that the TOC/P ratio is a relevant factor in optimizing phosphorus removal processes. Carbon substrates have long been recognized as key indicators for optimizing phosphorus removal processes [16, 29]. While a higher carbon ratio improves phosphorus removal efficiency, studies suggest that excessively high COD/P ratios can increase the proportion of GAOs, which compete with PAOs, reducing phosphorus removal efficiency [30]. As considering the COD:N:P ratio, it can be observed that MG treats a larger amount of N and P with less carbon source compared to other WWTPs. This suggests that there may be specific operational characteristics in the treatment process, such as wastewater return rates or influent COD/N ratio, rather than simple influent and effluent indices [31, 32]. The impact of sludge recirculation on removal efficiencies is supported by studies demonstrating that treatment processes with optimized return rates achieve higher pollutant removal, including reductions in COD and nitrogen [31].

The microbial community in the aerobic stage of the targeted wastewater treatment facilities was analyzed through NGS of the 16S rRNA gene. The composition of bacterial clusters in each sample was examined at the class and genus levels. The Shannon index, an indicator of taxonomic diversity differences between samples, showed the highest value (9.0) in the MG sample and the lowest value (7.6) in the AD sample, with other samples ranging between 8.6 and 8.9. In the case of AD, the mixed-liquor suspended solids in the aerobic stage were relatively low (Table 3). However, the DeNiPho method employs soccer ball-sized carriers, which supports a thriving microbial population on the surface of the carriers.

Fig. 1 represents bacteria detected at the class level with an abundance of 1% or more. Genera with less than 1% abundance were aggregated into the category “known others”, while “unknown others” represents those not classified at the species level. The highest proportion of bacteria was found to be Alphaproteobacteria in AD media, representing 25.93%. Specifically, the genus Hyphomicrobium accounted for 3.60%. Hyphomicrobium spp. are primarily identified as key microorganisms in denitrification systems [33]. Given the role of the media in nitrogen removal, AD media showed a high distribution of nitrogen-related microorganisms, including Nitrospira spp. (2.88%). The sample with the highest number of unclassified species was YH. A notable feature of YH is that Actinomycetia (class) within Actinobacteria was present at only 0.92%, whereas in other regions, it occupied a high proportion of 6.6–16.9%. Among the lower genera of Actionmycetia, Mycolicibacterium (0.94–2.87%) was present at over 1% in five wastewater treatment plants, but in YH, it was found at a very low proportion of 0.16%.

Microorganisms related to phosphorus removal are classified into PAOs and GAOs. PAOs refer to bacteria capable of storing intracellular phosphorus, and they cannot be defined by a single species. Table 4 summarizes the bacteria detected in metagenomic analysis among the phosphorus-removal related microorganisms outlined by previous research [10]. Among various PAOs, two genera known for their excellent phosphorus accumulation capabilities are Tetrasphaera and Ca. Accumulibacter. However, only Tetrasphaera was detected in the metagenomic analysis results.

Microorganisms classified as PAOs and detected at levels of 1% or higher were identified in the following order of genera: Tetrasphaera, and Dechloromonas. The presence of PAOs related to phosphorus removal efficiency was most significant in the AD (1.83%) and the YC (1.50%) (Fig. 2). Tetraspheara was found to be dominant in YC and YD, while Dechloromonas was most dominant in AD The Tetrasphaera spp. detected in the NGS analysis included Jenkinsii, Vanveenii, and Veronensis. According to Table 4, Jenkinsii accounted for 0.89% only in YC, Vanveenii was identified in all processes except YH, and Veronensis was detected only in YD.

In other research findings, Ca. Accumulibacter has been reported to exist in domestic and international full-scale wastewater treatment plant activated sludge at levels of 1–3% [10, 14]. In a side-stream process, Ca. Accumulibacter was identified at levels of 0.007–0.114% and Dechloromonas at 1.642–5.921% as PAOs [21]. However, in our study, Ca. Accumulibacter was not detected in the final analysis (data not shown.). To ensure effective DNA extraction, the bead-beating duration during the genomic DNA extraction process was extended. The extraction time was increased from 10 minutes to 15 minutes to improve DNA recovery. Despite these efforts, DNA quality control tests repeatedly failed. Even in successful quality control tests, the Adapter & Primer Trimming values decreased from an average of 120,000 to 80,000, indicating fewer cells were present and confirming the absence of Ca. Accmulibacter.

In this paper, the analysis was performed using the DADA2 program, which processes data based on amplicon sequence variants (ASVs). However, other studies have used both BIOiPLUG and UCLUST, which employ OTUs for data analysis [21]. ASVs and OTUs differ in their sequencing data processing methods: ASVs offer higher resolution and accuracy, while OTUs represent clusters of taxa for greater representativeness. These differences can impact the assessment of species richness and diversity [34, 35].

Competing GAOs were found to exist at high rates in the YC (4.96%) and the US (1.30%). Amoung Alphaproteobacteria, Defluviicoccus was identified as the most prevalent GAOs in these two WWTPs, with only the species Defluviicoccus vanus being detected. Another GAOs, Propionivibrio limicola from Betaproteobacteria, was found in YC at a rate of 1.41%. Conditions favoring GAO dominance include high organic carbon sources, low influent phosphorus, a water temperature above 20°C, and a long hydraulic retention time [16]. Previous research has indicated that there are no significant changes in the composition of activated sludge communities between summer and winter in full-scale wastewater treatment plants.

Therefore, this study focused only on samples collected in July. In the YC, the operating temperature of the aerobic stage in July was recorded at 27°C. Additionally, the influent BOD of the primary treated water entering the biological reactor was significantly higher compared to other WWTPs. These factors suggest that conditions favorable for GAO dominance may exist in the YC WWTP.

Correlation coefficients were calculated to evaluate the factors affecting the distribution of phosphorus removal-related microbes. Fig. 3 shows the correlations between physicochemical factors during microbial sampling and the read abundance ratios of PAOs and GAOs. The variables used in the correlation analysis included the removal efficiency (%) of key pollutants, the removal mounts (mg/L) of these pollutants, and the operational conditions of aerobic reactors, such as pH, temperature, salinity, and conductivity. This comprehensive approach allowed for a thorough evaluation of the factors influencing the distribution of phosphorus removal-related microbes. The removal efficiency (%) reflects the effectiveness of the treatment system in reducing pollutant concentrations, while the removal amounts (mg/L) provide insight into the actual mass of pollutants eliminated from wastewater. This distinction is crucial for understanding both the performance and the broader impact of the treatment processes.

The analysis evaluated phosphorus removal efficiency, removal quantity, and the operational conditions of aerobic reactors, including pH, temperature, salinity, and conductivity. According to the results (Fig. 3), the read abundance of PAOs (%) exhibited a strong correlation with both COD removal efficiency (%) and COD removal amount (mg/L), with a correlation coefficient of 0.94. However, no significant correlation was observed between the read abundance of PAOs (%) and either TP removal efficiency (%) or TP removal amount (mg). The correlation between BOD₅ and the read abundance of PAOs (%) showed a moderate correlation of 0.70 for the removal amount (mg/L). This indicates that a higher amount of high-quality biodegradable organic matter provides favorable conditions for the dominance of PAOs, a finding consistent with previous studies [15, 36, 37]. In this study, special attention was given to Tetrasphaera, a key PAO. Previous research suggests that Tetrasphaera–PAOs can uptake acetate and other substrates during the anaerobic phase but lack the ability to synthesize polyhydroxyalkanoates (PHA) as storage compounds [12]. Due to this metabolic characteristic, it is likely that high COD or BOD₅ values, rather than TP values, play a more significant role in promoting the dominance of Tetrasphaera–PAOs. On the other hand, the read abundance of GAOs (%) exhibited a very strong correlation of 0.97 with TOC removal. Similarly, the removal efficiency (%) of MLSS demonstrated a very strong negative correlation of −0.96, supporting this correlation. Among the operational conditions of the aerobic reactors, only salinity showed a strong correlation of 0.81 with the read abundance of PAOs (%).

Correlations between the abundance of predominant bacterial groups and operational parameters were analyzed, revealing significant relationships [21]. The influent COD exhibited a very strong positive correlation (r = 0.95) with Alphaproteobacteria, while a strong negative correlation was observed between BOD₅ sludge loading and Actinobacteria abundance (|r| > 0.8). However, other correlations between the influent wastewater parameters and bacterial abundance were limited, suggesting that the relative abundance of Chloroflexi, Actinobacteria, and Betaproteobacteria may be more influenced by substrate concentrations than broad operational metrics like COD or BOD₅. This underscores the importance of considering organic matter in distinct fractions rather than as a single parameter.

In our study, we found that specific operational parameters were strongly correlated with the abundance of related microorganisms. For PAOs, the correlation was particularly strong with the removal amount of COD (mg/L). Conversely, for GAO, a strong relationship was observed with the removal amount of TOC (mg/L). These findings highlight the distinct characteristics of TOC and COD as indicators of organic matter in WWTPs. While COD is a measure of the total oxygen demand required to completely oxidize organic matter, TOC directly quantifies the organic carbon content in wastewater. This distinction is crucial, as PAOs typically utilize readily biodegradable substrates more efficiently, while GAOs may thrive under different organic conditions. These findings emphasize the need to consider the specific contributions of these compounds when analyzing microbial dynamics in biological treatment systems.