Phenology and ecological role of aerobic anoxygenic phototrophs in freshwaters

Background Aerobic anoxygenic phototrophic (AAP) bacteria are heterotrophic bacteria that supply their metabolism with light energy harvested by bacteriochlorophyll-a-containing reaction centers. Despite their substantial contribution to bacterial biomass, microbial food webs, and carbon cycle, their phenology in freshwater lakes remains unknown. Hence, we investigated seasonal variations of AAP abundance and community composition biweekly across 3 years in a temperate, meso-oligotrophic freshwater lake. Results AAP bacteria displayed a clear seasonal trend with a spring maximum following the bloom of phytoplankton and a secondary maximum in autumn. As the AAP bacteria represent a highly diverse assemblage of species, we followed their seasonal succession using the amplicon sequencing of the pufM marker gene. To enhance the accuracy of the taxonomic assignment, we developed new pufM primers that generate longer amplicons and compiled the currently largest database of pufM genes, comprising 3633 reference sequences spanning all phyla known to contain AAP species. With this novel resource, we demonstrated that the majority of the species appeared during specific phases of the seasonal cycle, with less than 2% of AAP species detected during the whole year. AAP community presented an indigenous freshwater nature characterized by high resilience and heterogenic adaptations to varying conditions of the freshwater environment. Conclusions Our findings highlight the substantial contribution of AAP bacteria to the carbon flow and ecological dynamics of lakes and unveil a recurrent and dynamic seasonal succession of the AAP community. By integrating this information with the indicator of primary production (Chlorophyll-a) and existing ecological models, we show that AAP bacteria play a pivotal role in the recycling of dissolved organic matter released during spring phytoplankton bloom. We suggest a potential role of AAP bacteria within the context of the PEG model and their consideration in further ecological models. Supplementary Information The online version contains supplementary material available at 10.1186/s40168-024-01786-0.


Introduction
Recurrent seasonal changes of aquatic microbial communities are among the best-studied phenomena in freshwater lakes and reservoirs.The Plankton Ecology Group (PEG) model initially described the dynamic interactions between phytoplankton and zooplankton [1] and was later amended with the eutrophic and oligotrophic scenarios and role description of heterotrophic protists [2].Subsequently, the importance of bacterioplankton was revealed, especially during the spring phytoplankton bloom [3][4][5], increasing our understanding of the contribution of microorganisms to the functioning of limnic ecosystems [6].Bacteria represent an important part of aquatic microbial communities.They generate fresh particulate organic matter by utilizing dissolved organic carbon (DOC) and render it accessible to organisms at upper trophic levels [7].However, the role of photoheterotrophic bacteria, which present a significant part of bacterial biomass and activity [8], remains overlooked.
Aerobic anoxygenic phototrophic (AAP) bacteria are a functional group of photoheterotrophs that rely upon external sources of organic carbon and supplement their metabolism with energy obtained from light through bacteriochlorophyll-a (BChl-a) type II reaction centers.The capacity to harvest light energy enables AAP bacteria to reduce their respiration and increase biomass yield [9,10].Moreover, the AAP community shows higher growth rates, larger cell sizes, and greater activity than heterotrophic bacteria [11][12][13][14][15][16].Photoheterotrophy by AAP bacteria increases carbon transfer efficiency, enlarging the availability of biomass for upper trophic levels and reducing CO 2 emitted to the atmosphere [17].However, little is known about the phenology of the AAP community and the absence of exhaustive seasonal sampling hampers the understanding of their role in lakes.AAP bacteria peak during spring in lakes, when they may account for up to 22% of bacteria [18,19].Their abundances and diversity dynamics correlate with irradiance, temperature, chlorophyll-a (Chl-a), oxygen, and DOC [18,[20][21][22][23][24].
One of the obstacles in the study of AAP bacteria is the fact that they do not represent a monophyletic group.On the contrary, phototrophic genes have been gained and lost multiple times in closely related species [25,26].Therefore, AAP species cannot be identified based on the most common marker used in community studies, the 16S rRNA gene.Instead, the pufM gene, which encodes the subunit M of the anoxygenic type-II reaction center, has been widely employed to study AAP communities [15,[27][28][29][30][31][32].However, these studies were unsuccessful in providing a taxonomic assignment for abundant pufM phylotypes.This is caused by the low taxonomic resolution of the short amplicon sequences and the lack of a curated reference database.The increased availability of metagenome-assembled and single-cell amplified genomes (MAGs and SAGs) has expanded our knowledge of metabolic potential within multiple bacterial lineages and should allow us to establish a comprehensive pufM database for amplicon assignment.
To improve the taxonomic assignment, we designed a novel primer set targeting a larger 450 bp region of the pufM gene and compiled an extensive database of 3633 non-redundant pufM gene sequences from existing genome and metagenome sequence datasets.We applied this novel metabarcoding assay to 215 samples from 3 years, collected from meso-oligotrophic freshwater Cep lake (Czechia) at biweekly intervals from multiple depths.We hypothesized that the AAP community would show a recurrent seasonal succession.Due to their specific metabolism, this succession would exhibit different abundance patterns than overall heterotrophic bacteria.Specifically, we expected that AAP abundance would peak during the spring phytoplankton bloom.Additionally, since they are a functional, taxonomically diverse group, we surmise that the spring AAP peak is orchestrated by specific phylotypes, rather than the involvement of the entire AAP community.

Sampling and measuring environmental variables
Samples were collected biweekly from April 2017 to December 2019 from the freshwater Cep lake (48°92′49.24″N,14°88′68.11″E).This meso-oligotrophic lake is located in the Třeboň Basin Protected Landscape Area, Czechia, and has an area of 130 ha and a maximum depth of 12 m.Five liters of water were collected from 0.5, 2, 5, and 8 m using a 3-L Ruttner water sampler (KC Denmark A/S, Denmark) and transported to the laboratory in closed plastic containers in a cooler box, which were pre-rinsed three times with the sampled water.Temperature and oxygen profiles were taken with an EXO1 multi-parameter probe (YSI Inc., Yellow Springs, USA).Total and AAP bacterial abundances were counted in Zeiss Axio Imager.D2 epifluorescence microscope equipped with Collibri LED module illumination system (Carl Zeiss, Jena, Germany).Microscopic samples were excited at 325-370 nm, 450-490 nm, and 545-565 nm and the AAP cells were detected by their BChl-a autofluorescence at wavelengths > 850 nm in epifluorescence.Since Chl-a autofluorescence spectra emission also encompasses the 890 nm wavelength, cells that show fluorescence emission at 690 and 573 nm, corresponding with Chl-a and phycoerythrine, respectively, were not counted as AAP bacteria [17,33,34].Concentrations of Chl-a and BChl-a were determined in organic solvent extracts by reversed-phase high-performance liquid chromatography [31].The quantification of environmental nutrients was performed as described in Procházková, 1959 (nitrate); Murphy and Riley, 1962 (phosphate); Kopáček and Hejzlar, 1993 (total phosphorous); Kopáček and Procházková, 1993 (ammonia) and Shabarova et al. 2021 (DOC) [35][36][37][38][39].
All pufM sequences were pooled and duplicated sequences were removed.Protein sequences were aligned with MAFFT v7.453 (-maxiterate 1000 -localpair) [53] and a maximum likelihood tree was calculated using iqtree2 [54] with automatic model selection performed by ModelFinder [55], and 1000 iterations of ultrafast bootstrapping with 1000 rounds of SH-aLRT testing (-alrt 1000 -B 1000) [56].pufL sequences were identified as they formed a long branch.Bona-fide pufM sequences were retained (Supplementary file S1), and alignment and phylogenetic trees were redone and visualized using iTOL [57].The environmental origin of each sequence was obtained manually from source databases (Supplementary file S2).

AAP community analysis by pufM gene amplicon sequencing
Between 300 and 1460 ml of water was filtered through sterile 0.2 µm Nucleopore Track-Etch Membrane filters (Whatman ® , Maidstone, United Kingdom) that were immediately placed inside sterile cryogenic vials (Biologix Group Limited, Jinan, Shandong China) containing 0.55 g of sterile zirconium beads, flash-frozen in liquid nitrogen and stored at − 80 °C until DNA extraction (max.6 months).Total nucleic acids were chemically extracted according to Griffiths et al. 2000 [58] with modifications [59], re-suspended in 35 µl of DNase and RNase-free water (MP Biomedicals, Solon, OH, USA), and stored at -20 °C.Concentration and quality of the extracts were checked using NanoDrop (Thermo Fisher Scientific).
To improve the accuracy of the taxonomic assignation and reduce the number of unclassified amplicon sequence variants (ASVs), a new primer pair for pufM gene was designed.pufM_uniF primer [27] was used as a reverse (pufM_UniFRC in the current study, 5´-RAANGGR TTR TARWANARRTTNCC-3') and pufM_longF was designed ~ 450 bp upstream (5'-YGGSCCG WTC TAYST-SGG-3') using a pre-existing database of 1500 sequences [31].The specificity and coverage of the new primer pair were tested in comparison to the commonly used pufM primers [27,32] against the new pufM database.The analysis was done in Geneious Prime (v2023.0.1) with up to three mismatches in the binding region and in both forward and reverse directions.The primers' specificity was also tested separately for Pseudomonadota (formerly known as Proteobacteria), Alphaproteobacteria, Gammaproteobacteria, Gemmatimonadota, Chloroflexota, Myxococcota, and Eremiobacterota based on alignments done in Geneious Prime by MUSCLE alignment (v5.1.).
The PCR conditions were optimized using genomic DNA from Gemmatimonas phototrophica (Gemmatimonadota), Sphingomonas glacialis (Alphaproteobacteria) and Congregibacter litoralis (Gammaproteobacteria), and environmental DNA from the current sampling.The final conditions were as follows: initial denaturation for 3 min at 98 °C, 35 cycles of 98 °C for 15 s, 52 °C for 30 s, 72 °C for 18 s, and final elongation at 72 °C for 5 min.Triplicate PCR reactions (20 μL) using Phusion ™ High-Fidelity PCR MasterMix (Thermo Fisher Scientific, USA) were pooled and the amplicons of ~ 450 bp were purified from 1.5% agarose (MP Roche, Germany) gel using the Wizzard SV Gel and PCR clean system (Promega, USA) and quantified with Qubit dsDNA HS assay (Thermo Fisher Scientific, USA).Samples were randomly distributed within two runs to account for the batch effect and sequenced on Illumina Miseq 2 × 300 bp PE (Macrogen, South Korea).
The taxonomic assignment was done through phylogenetic placement using The Evolutionary Placement Algorithm v0.3.5 [65] that placed the ASVs into the phylogenetic tree calculated from the new reference database sequences that were back-translated from protein alignments using trimAl [66].The taxonomic assignation was handled according to the ASV phylogenetic position using Gappa [67] (Supplementary file S3, Taxonomy sheet).

Phytoplankton community analysis based on 16S rRNA gene amplicons
The V3-V4 region of the bacterial 16S rRNA gene was amplified using 341F and 785R primer pair [68] as described in Piwosz et al. 2022 [17] The subset of sequences assigned to Chloroplast was extracted and their taxonomy was further affiliated using a curated reference database of the plastidial 16S rRNA gene: Phy-toRef [69].Bar plots were visualized using ggplot v3.4.3 [70].
Relationships between environmental data and AAP community were analyzed using distance-based linear models (DistML) [78,79] in the PERMANOVA + add-on package of the PRIMER7 software [80] (Primer Ltd., Lutton, UK).From strongly correlated environmental variables (correlation coefficient > 0.6) only one was selected for further analysis.The model was calculated on the CLR-transformed relative abundance data of AAP bacteria [74], using a stepwise selection procedure.The best model was selected based on statistical significance (9999 permutations) and the value of Akaike's Information Criterion (AICc).

AAP core and network community analyses
ASVs present in more than 80% of the samples from the 3 years and four depths, were considered the AAP lake core microbiome.The percent contribution of each core ASV to their respective maximum percent contribution was calculated and plotted with bubble plots using Phyloseq v1.30.0 [75] and ggplot2 v3.3.6 [70].
SparCC analysis was applied to calculate lake community co-occurrence correlations from the compositional data [81].Only correlations with pseudo p value < 0.02 and stronger correlation than ± 0.7 were selected.The network was plotted using Cytoscape v3.9.1 [82].

Time series and trend lines
Interannual trend analysis was done in TTR package v0.24.3 (https:// github.com/ joshu aulri ch/ TTR) in R v 4.3.0.Raw data on total and AAP bacterial abundances, temperature, and Chl-a concentrations were averaged for months and depths and transformed into time series assuming annual frequency.They were decomposed into trend, seasonal, and random components using decompose with default settings.Spearman correlation between interannual trends of AAPs abundance and temperature and Chl-a concentrations was done for the extracted trend component.

New database and longer amplicons enhance taxonomic assignments of pufM gene ASVs
In order to improve the taxonomic assignment of the pufM gene amplicons, we constructed a new reference database containing 3633 pufM sequences (> 646 bp) from Pseudomonadota (synonym for Proteobacteria), Gemmatimonadota, Chloroflexota, Eremiobacteriota, and Myxococcota (Supplementary file S1).The database includes 529 genera, 114 families, 53 orders, and 9 classes (Supplementary Figure S1) from cultured species and MAGs originating from a wide variety of habitats (Supplementary file S2), mostly from freshwater (2140 sequences) and marine environments (381 sequences).
The newly designed primer set covers 80.5% of the new database at maximum of three mismatches (Supplementary file S4).The older, highly degenerated primer pairs, UniF + UniR, pufMF + pufMR and UniF + pufM_WAW, cover 98.9%, 85% and 96.2%, respectively.However, the amplicon length of the novel primers is about two times longer (~ 450 bp) allowing for the proper taxonomic assignation of more than 95% of the alphaproteobacterial and above 75% of the gammaproteobacterial reads at the order level.Additionally, 38 alphaproteobacterial, 36 gammaproteobacterial, and 6 Gemmatimonadota genera were detected (Supplementary Figures S2, S3, and S4).

Seasonal changes in Cep lake
Environmental conditions in Cep Lake showed seasonal dynamics typical for a temperate freshwater lake (Supplementary file S5).In January and February, the lake was partially frozen and stratified from April/May until September with maximum temperatures around 24 °C in July and August.The metalimnion was located between 5 and 8 m depth in 2017, and between 2 and 5 m in 2018 and 2019.In all 3 years, the autumnal mixing, characterized by higher values of dissolved oxygen and lower temperatures, was initiated in October (Supplementary Figure S5A-B).
Chlorophyll-a measurements varied throughout the year with seasonal maxima representing spring and autumn phytoplankton blooms (Supplementary Figure S5C).The spring phytoplankton bloom terminated at the onset of stratification, and was composed, according to 16S amplicons affiliated to plastids, mostly by Bacillariophyta and Chrysophyceae (Supplementary Figure S6).

Seasonal dynamics of AAP community composition
The maximal AAP abundances (3.42-5.50× 10 5 cells mL −1 ) corresponded to 15-20% of the total bacteria (Supplementary Figure S5D-E) closely following the spring phytoplankton blooms.A second AAP bacterial peak occurred towards the end of summer, before the autumn phytoplankton peak.Alpha diversity of the AAP community was lower during the spring abundance peaks and rose during the second part of the year (Supplementary Figure S5F).
The AAP community followed an annual recurrent pattern during the three consecutive years, with a distinction between the epilimnion and hypolimnion communities during the stratified period (Fig. 1A).Samples from the same season in different years were more similar to each other than samples from different seasons in the same year, indicating the persistent temporal succession of the community and similar interannual community ).For 2018 and 2019, 2 years for which nutrient data is available (Supplementary file S5), phosphorus and ammonia increased this explanation by 6.75%, up to 29.86%.
Interestingly, in addition to a seasonal cycle, we observed an interannual variation in AAP community composition (Fig. 1B).The decomposition of time series on monthly averaged values for the whole water column showed increasing interannual trends in temperature and AAP abundance during 3 years of sampling, and a decreasing trend in Chl-a concentration (Supplementary Figure S7).Trends of temperature and AAP abundance were significantly correlated (Spearman correlation coefficient rho = 0.8, p value < 0.0001).
The AAP community was dominated by Gammaproteobacteria over the whole water column, with an average relative contribution exceeding 50% and reaching up to 90% during stratification (Supplementary Figure S8).Alphaproteobacteria was the second most abundant class, showing the maxima contributions in spring and autumn, reaching over 50% in the spring of 2018.Classes Gemmatimonadetes (Gemmatimonadota) and Myxococcia (Myxococcota) made up 5% and 2% of the AAP community, respectively.Unclassified Pseudomonadota and unclassified Myxococcota showed transient contributions of < 9% and < 1% of the AAP community, respectively.Chloroflexota and Eremiobacteriota were not detected.
The 100 most abundant ASVs (based on their average relative abundances) comprised 75% of the reads and exhibited seasonal recurrence, peaking every year at specific times of the year (Fig. 2).The majority of ASVs demonstrated a transient contribution and were generally absent outside their maxima (e.g., ASV28 (Rhizobiales)),   while a few showed pronounced relative abundance throughout the year (e.g., ASV4 (Rhodoferax)).Interestingly, ASV67 (Limnohabitans), whose relative abundance was on average < 0.36%, was the sole ASV detected in every sample.Out of 1588 pufM ASVs (Supplementary file S3, Reference ASV sheet), the stable part of the AAP community (defined here as ASVs present in > 80% of the samples) consisted of only 22 ASVs (Supplementary Figure S9): 8 Rhodoferax, 4 Limnohabitans, 2 Aestuariivirga, 1 Methylobacterium, 1 Rubrivivax, and 6 other Burkholderiales.These core ASVs varied largely in their contribution from the most abundant ASV2 (Aestuariivirga, with an average relative abundance of 4%) to the least abundant ASV237 (unclassified Burkholderiales, with an average relative abundance of 0.06%).Their seasonal dynamics differed substantially throughout the year and distinct relative abundance patterns were observed even for ASVs from the same genus.For instance, ASV5 and ASV49, both Rhodoferax, peaked in autumn and spring, respectively.Similar differences were observed for two Aestuariivirga: ASV2 peaked during the spring mixing period (from March to May), while ASV62 showed its highest contribution during the summer stratification.It is noticeable that the core AAP community also included ASVs outside the 100 most abundant, such as ASV115 (Rhodoferax) and ASV237 (unclassified Burkholderiales), which had a low but steady contribution during the whole sampling season.Furthermore, the dynamics of some core ASVs, such as ASV31, were different each year.
To identify phylotypes most contributing to the difference in AAP community composition during spring and autumn peaks, we carried out an analysis of compositions of microbiomes with bias correction (Supplementary file S7).Composition of genera and orders contributing to the AAP bacterial peaks was different (Fig. 3).Spring peak consisted of a higher prevalence of Alphaproteobacteria versus Gammaproteobacteria genera (9 vs 5), while during autumn, the community was more diverse and included also Gemmatimonadota and Myxococcota.The highest genera contributors to the spring peak (log fold change > 2) were RFPW01, UBA1936, Cypionkella, and Rhodoferax, while UBA964 and UBA5518 (Gammaproteobacteria: Steroidobacterales and Pseudomonadales, respectively) contributed most to the late summer peak.The only genera substantially contributing to both peaks was Aestuariivirga.

Network analysis
To study possible interactions between AAP bacteria, we performed a network co-occurrence analysis.The network calculated for the entire lake concluded that 99 ASVs presented 139 significant interactions (92% positive; Fig. 4).Thirteen highly connected nodes with more than six co-occurrence correlations were identified as hubs.They belonged to the genus UBA964, Steroidobacterales, Rubrivivax, and unclassified Burkholderiaceae from Gammaproteobacteria.The most connected nodes from Alphaproteobacteria belonged to Aestuariivirga and Rhodobacteraceae (4 edges each).The majority of correlations (~ 60%) occurred between ASVs from the same genus, family, or order, creating groups of densely connected nodes (e.g.genus UBA964).Furthermore, some ASVs peaking in spring correlated positively with each other (Aestuariviirga and unclassified Burkholderiaceae) in contrast to Rhodoferax and Cypionkella.A similar pattern was observed for ASVs peaking in summer-autumn.

Discussion
Seasonal succession of planktonic communities in temperate lakes has been intensively studied [2].The PEG model defines the key phases in the annual development of ecological succession as well as the interactions between different organisms in aquatic ecosystems.
The annually recurrent phenomena in freshwater lakes include a spring phytoplankton bloom followed by a zooplankton-induced clear-water phase in early summer, a late-summer phytoplankton bloom, and a period of low productivity in winter.Recently, seasonal dynamics of heterotrophic bacteria have been incorporated into the model [5], as they rapidly respond to the transitions in the lakes' pelagic functioning, especially during the phytoplankton bloom [3,5,83].Whilst AAP bacteria have been shown to substantially contribute to the bacterial abundance, biomass, and activity in freshwater lakes [13,17], they are not considered in the PEG model [1,2].Dynamics of microbial communities are typically investigated by sequencing of the 16S rRNA gene and such analyses do not allow for disentangling the metabolic functionality of bacteria, especially when traits of interest follow a heterogenic pattern of presence within the same taxonomic ranks.This is the case for AAP bacteria, where members from the same genus might or might not contain the ability to carry out anoxygenic photosynthesis [26].Amplicon analysis of a functional gene may overcome this hindrance, but it requires a comprehensive database with taxonomically assigned reference sequences.Currently, some phylogenetic clades (A-L) cannot be assigned to the genus or even order level [84,85].Moreover, different researchers assemble databases for their environment of interest, with different quality thresholds and criteria [17,86], which hampers direct comparison between different studies.Finally, short amplicons obtained with the  4 AAP community co-occurrence correlation network.The radius of the nodes (ASVs) is directly proportional to the number of significant correlations.Colour of the nodes shows the taxonomy at the class level with Alphaproteobacteria in red-orange colour scale, Gammaproteobacteria blue colour scale, unclassified Burkholderiaceae (grey), and Gemmatimonadota (brown).Blue lines indicate positive and red lines have negative correlations.Rectangles with solid lines indicate nodes that interact with other nodes from the same taxonomic rank (genus or family) and positive log fold change for Spring AAP peak are highlighted in yellow rectangles fill and yellow node stroke while green rectangles fill and green node stroke for autumn AAP peak most commonly used primer combinations, puf_UniF-puf_ UniR, puf_UniF-pufM_WAW, and pufMF-pufM_WAW [27,32] often do not allow for taxonomic assignment below the class or order level, resulting in a substantial number of unclassified reads [18,29,86,87].

Comprehensive database and longer amplicons allow for improved taxonomic assignments of pufM ASVs
We constructed the largest curated pufM database to date that includes sequences from essentially all environments, with a special high representation of freshwater lakes (Supplementary Figure S1, Supplementary file S1 and S2).Phylogenetic trees based on pufM gene do not match 16S rRNA phylogeny due to horizontal gene transfer events [88][89][90].In contrast, the taxonomic assignment of genomes and MAGs based on the whole genome or 120 selected marker genes is more accurate and consistent [91].Thus, we included only sequences originating from taxonomically assigned genomes and MAGs, excluding all environmental pufM sequences classified into phylogroups based on phylogenetic analysis [84].During our quality control, sequences originating from Bdellovibrionota, Verrucomicrobiota, Omnitrophota, Planctomycetota, or Bacteroidota were excluded from the final database.Some members of these phyla were reported to encode pufM [87].However, our manual inspection revealed that their pufM gene was present in short contigs and often found as the only phototrophic gene, which does not warrant phototrophic functionality.Furthermore, pufM genes from these phyla did not form a monophyletic clade, suggesting dubious multiple and independent events of horizontal gene transfer.Thus, only members of Pseudomonadota, Chloroflexota, Gemmatimonadota, Eremiobacteriota, and Myxococcota were included.Our choice of sequences ensures high quality and allows for future extensions the database as more metagenomic data is produced and more AAP bacteria are cultured, enhancing its fidelity and functionality.Moreover, as our database contains entire or almost entire pufM gene sequences, it can be used for amplicon taxonomy assignment independently of the actual primers used.
The pufM gene is one of the most conserved genes from the puf operon which codifies the genes for the synthesis of the anoxygenic photosynthesis apparatus [92,93].Commonly used puf_UniF-puf_UniR, puf_UniF-pufM_WAW, and pufMF-pufM_WAW primer pairs [27,32] hybridize on the most conserved regions at the end of the gene, separated by ~ 110-160 bp.They show high coverages (Supplementary file S4) but produce short amplicons that hamper taxonomic assignation below the class or order level resulting in a high fraction of unclassified reads [18,29,86,87].Thus, we designed a new primer set producing longer amplicons to increase the taxonomic resolution, as has been shown for other genes [94].The new primer pair has lower in silico coverage against our new database than puf_UniF-puf_UniR [27] (80.5 vs 98.9%).Some groups, such as Chloroflexota, Aquidulcibacter, and Polynucleobacter, were poorly covered, which may explain their absence in our amplicons.Nevertheless, the number of sequences identified at every taxonomic level increased compared to a previous study in the same lake [18]: 95% of the alphaproteobacterial and above 75% of gammaproteobacterial reads were classified at the order level.Additionally, the number of newly detected genera was substantially higher (80 vs 12) and the Shannon index showed a wider range of diversity (Supplementary Figure S5F) since longer amplicons enable us to detect more nucleotide variations, and thus revealing higher diversity, advancing our knowledge on AAP community composition.

Phenology of AAP bacteria and consideration of the PEG model
The seasonal succession of the AAP community revealed differential strategies of adaptation to environmental conditions, unveiling generalist AAP bacteria appearing most of the time, whereas specialists or opportunists showed a transient contribution to the AAP community (Fig. 2).Within the generalists, we identified the core AAP community that consistently contributed throughout the seasons for three consecutive years and across all depths (Supplementary Figure S9).The coexistence of the core AAP community, the dominance of positive correlations in the networks, and the numerous correlations between ASVs of similar taxonomic ranks (Fig. 4) suggest partial metabolic redundancy within some closely related AAP bacteria that maintain functional kinship.In contrast, the complex seasonal succession pattern indicates that the lake's AAP community is extremely diverse, with over 100 reported genera of AAP bacteria (Supplementary Figure S1, S2, and S3).AAP communities represent a large functional repertoire (even within the same genus) allowing for niche speciation via temporal succession, facilitating their geographical coexistence.Finally, the 3-year recurrence of the AAP community (Fig. 1A) documents its indigenous character in this freshwater lake and the higher importance of selection over the environmental drift and dispersal processes at a short temporal scale [95].Changes in the AAP abundance coincided with shifts in their community composition indicating that abundance peaks were caused by specific phylotypes.These phylotypes differed between both abundance peaks (Fig. 3).Generally, the AAP community was dominated by Gammaproteobacteria except for spring abundance peaks, when Alphaproteobacteria, which already have shown higher phototrophic activities in spring [31], increased their contribution.Additionally, the directional interannual variation of the AAP community (Fig. 1B) signifies the evolution of AAP populations, potentially influenced by changes in environmental and biological variables such as temperature or Chl-a (Supplementary Figure S7).
While the PEG model has enhanced our understanding of seasonal patterns, it still does not encompass all aquatic components, such as viruses or specific functional bacterial groups.This includes AAP bacteria, which are characterized by a heterogenic behavior but still represent an important functional group, fulfilling valuable ecological and biochemical processes in the aquatic environment.For that reason, amending them into the present PEG model will certainly improve our understanding of aquatic community functioning.
Cep Lake is representative of a meso-oligotrophic temperate lake in the northern hemisphere [96], thus our conclusions might be applied to other similar lakes.AAP bacteria played an important role during and shortly after the spring bloom when their abundance and contribution to the total bacterial community were recurrently the highest.This spring AAP abundance peak preceded that of the overall heterotrophic bacteria (Fig. 5).These results are consistent with the previous study in the same lake [18].The faster response of AAP bacteria highlights that photoheterotrophy confers a distinct metabolic importance as food for bacterivores in microbial food webs is well documented [97][98][99] and due to their, on average, larger cell size and higher activity than other heterotrophic bacteria [13,14], they might contribute disproportionally to the carbon cycling despite their relatively low abundances [19,99].Additionally, Chl-a concentration has been identified as a variable explaining the dynamics of the AAP bacterial community and it is plausible to assume that the spring abundance peaks of AAP bacteria are triggered by the excess of carbon released by the phytoplankton bloom (mostly diatoms; Supplementary Figure S6) and the lack of grazing pressure after winter.Moreover, AAPs have a highly efficient photoheterotrophic metabolism [17] increasing secondary bacterial production and disposing of more carbon to higher trophic levels via the microbial loop.This emphasizes the urgent need for more quantitative studies to further decipher carbon transfers along microbial and classical food webs.The AAP bacterial peak is terminated by selective and extensive grazing of bacterivorous protists and mesozooplankton that are also present in summer [83,98,100].The absence of a pronounced AAPs peak following the phytoplankton bloom in autumn (Fig. 5) could be attributed to the higher grazing pressure, distinct phytoplankton composition (Supplementary Figure S6), decreasing temperatures (Supplementary file S6), and/or the decreasing light availability at shorter day length [2,18,30,101].Altogether, we propose for the first time the inclusion of the AAP bacteria into the PEG model, encouraging other microbial ecologist to account for their role in other lakes of different trophic status.Finally, our study provides novel insight into the ecology of phototrophic Myxococcota.While their average contribution was low, they were detected during stratification over three consecutive years (Supplementary Figure S8) and constituted a member of the summerautumn AAPs peak, emphasizing their potential significance in microbial communities during summer as they showed a potentially predatory and photoheterotrophic metabolism [102].

Conclusions
Our study revealed annual recurrent seasonal patterns of AAP bacteria in a freshwater lake, supporting the potential inclusion of this important functional group into the PEG model.The high abundance of AAP bacteria during the spring phytoplankton bloom highlights their crucial role in recycling phytoplankton-derived dissolved organic matter and their role in aquatic food webs, which needs to be further quantified and better understood.Differential contribution patterns of the core community and temporal succession of the AAP community indicate strong competition within AAP bacteria communities, which forces them to conduct temporal niche partitioning in order to geographically coexist.In contrast, positive co-occurrence correlations between closely related AAP bacteria indicated their functional redundancy.Our findings provide unprecedented insights into the phenology of AAP bacteria in a temperate freshwater lake, blazing a trail for future studies to verify the proposed role in other types of lakes.and technical support for the LIMNOS sequencing data set and all related environmental data.

Fig. 1
Fig. 1 Development of AAP community structure.Principal component analysis of centered log-ratio transformed AAP community composition.Each point represents a sample with 0.5 m (red circle), 2 m (orange triangle), 5 m (green square) and 8 m (blue cross) in PC1 and PC2 axis (A) and 0.5 m (circle), 2 m (triangle), 5 m (square) and 8 m (cross) coloured according to the date of sampling in PC1 and PC3 axis (B).Dashed line and arrows in panel A indicate AAP community succession following an annual chronological direction

Fig. 2
Fig. 2 Recurrence of 100 most abundant ASVs.Relative abundance is individually normalized in 0.5 m, 2 m, 5 m, and 8 m during the different seasons of the 3-year sampling.Coloured box on the left indicate the taxonomic assignation at class level of each ASV, red for Alphaproteobacteria and blue for Gammaproteobacteria. Colour bar on the top and bottom indicate the seasons of the year

Fig. 3
Fig.3 Community composition of AAP abundance peaks.Analysis of compositions of microbiomes showing log fold change values at order and genus level in Spring (A) and in Summer-Autumn peaks (B) with green colour scale for Burkholderiales, blue colour scale for other Gammaproteobacteria, red colour scale for Alphaproteobacteria, yellow for Gemmatimonadota and grey for Myxococcota orders Fig.4 AAP community co-occurrence correlation network.The radius of the nodes (ASVs) is directly proportional to the number of significant correlations.Colour of the nodes shows the taxonomy at the class level with Alphaproteobacteria in red-orange colour scale, Gammaproteobacteria blue colour scale, unclassified Burkholderiaceae (grey), and Gemmatimonadota (brown).Blue lines indicate positive and red lines have negative correlations.Rectangles with solid lines indicate nodes that interact with other nodes from the same taxonomic rank (genus or family) and positive log fold change for Spring AAP peak are highlighted in yellow rectangles fill and yellow node stroke while green rectangles fill and green node stroke for autumn AAP peak

Fig. 5
Fig.5 Proposed inclusion of AAP bacteria in the PEG model.Annual succession patterns of microbial communities for (A) Phytoplankton and zooplankton according to original PEG model in oligotrophic scenario(Sommer et al., 1986), and (B) monthly averaged annual succession pattern of AAP abundance in CEP lake, phytoplankton dynamics through Chl-a, and of bacterial abundance.Trend lines are normalized to maxima and minima values for each variable, light transparent areas indicate 95% confidence interval