Samples were provided by Aker BP (Lysaker, Norway) and originated from four North Sea oil fields located 200–300 km offshore of Stavanger, Norway: Alvheim, Vilje, Volund, and Bøyla (Fig. 1; Table 1). The samples consisted of reservoir fluid (RF) samples and produced water (PW) samples, all of which were obtained directly at the wellhead.

For each production well, paired PW and RF samples were collected using standard 8-L miniseparators designed for representative sampling of multiphase production streams. The miniseparators were connected directly to the production wells only during sampling, with no intermediate processing equipment between the wellhead and the sampling point. Before sampling, the miniseparator was thoroughly flushed with reservoir fluids for approximately 2 minutes at an inlet pressure of 20 bar to ensure rapid fluid turnover and to minimize carryover of potential contaminants. The sampling cylinder was drained and refilled before the final sample was collected. Consequently, any influence of the sampling equipment on microbial composition is considered minimal.

Reservoir fluid samples contained varying amounts of an aqueous phase (0.1 to 10% of total RF volume), as well as minor amounts of suspended solids. This aqueous fraction was not removed prior to downstream processing and was analytically treated as part of the RF sample, whereas PW samples represented the separately collected aqueous phase of the same well stream. The water content and amounts of suspended solids varied among replicate samples.

Temperature was not actively controlled during sampling but occurred at well stream temperatures close to the reservoir conditions.

For all four oil fields, in-situ reservoir temperature varied between 63 °C and 66 °C at producing depths ranging from 2,100 to 2,200 meters below seafloor (mbsf). Water depths ranged from 120 to 130 meters. All four fields had the pressure maintained by underlying aquifers, produced water injections, or both. To the best of our knowledge, no seawater injection was carried out.

The Alvheim field comprises six discoveries, of which three were sampled in this study (East Kameleon [E. Kam.], Kneeler B [KNB] and Boa); Vilje, Volund, and Bøyla are separate fields tied back to the Alvheim FPSO. Their discovery years, production start dates, distances from Alvheim, and pressure-maintenance methods are given in Table 1.

Samples were collected from each of the Volund, Bøyla, Vilje, E. Kam., Boa, and KNB reservoirs and oil fields (Fig. 1). All samples were collected aboard the Alvheim FPSO vessel directly into sterile, 0.5-litre glass bottles capped with polytetrafluoroethylene (PTFE) lined caps without headspace to maintain anoxic conditions, kept at room temperature during the 3 days needed for transport, and stored for approximately one week at 4 °C until analysis. Immediate freezing was not feasible due to logistical and safety constraints associated with offshore sampling and helicopter transport.

Upon arrival, samples were opened and subdivided into reaction tubes in an anoxic glove box to avoid oxygen ingress. Samples were stored with a nitrogen headspace to ensure strict anoxia and avoid the growth of microaerophilic organisms. All samples were delivered in duplicate bottles, each of which was subdivided into duplicate aliquots, resulting in four technical replicates per sample.

RF samples consisted of a multiphase matrix containing crude oil, suspended solids, and a variable aqueous fraction. During DNA extraction, no physical separation or removal of the aqueous fraction from RF samples was performed prior to processing. In contrast, PW samples represented the isolated aqueous phase of the same reservoir fluids and were processed directly without oil-phase treatment.

In an anoxic glove box, RF samples were mixed 1:1 with sterile phosphate-buffered saline (PBS) containing 0.1 vol% Tween 80. Following vigorous shaking, samples were stored overnight at 55 °C. The next day, the polar fraction of the solution was filtered through Sterivex filters (pore size 0.22 µm; Merck KGaA, Darmstadt, Germany). To avoid contamination with airborne microbes, handling of the filter and subsequent DNA extraction were carried out in a laminar flow cabinet. Filter casings were cracked open, and filters were transferred to bead tubes. Before continuing with the manufacturer’s protocol for the PowerSoil Pro DNA extraction kit (DNEasy PowerSoil Pro kit, Qiagen, Hilden, Germany; Cat. No. 47014), 40 µl of 10% sodium dodecyl sulphate (SDS) was added to the bead tubes. The addition of SDS increased DNA yield from crude oil samples (Alibrandi et al. 2023). Unlike RF samples, PW samples required no cell/particulate separation from the oil matrix and, therefore, no PBS-Tween 80 treatment was required. For reproducibility, PW samples were incubated overnight at 55 °C (as for RF samples), then filtered and processed as described above. DNA extractions on blank Sterivex filters were conducted as negative controls. DNA was quantified using a Qubit 2.0 device and a dsDNA HS assay kit (Invitrogen, Carlsbad, USA).

All samples were processed in triplicate. Bacterial and archaeal 16S rRNA gene fragments (V4 hypervariable region) were PCR amplified with the universal barcoded primer pair 515F (5'-GTG TGY CAG CMG CCG CGG TAA-3') and 806R (5'- CCG GAC TAC NVG GGT WTC TAA T-3'). The final volume of each reaction mixture was 50 µl, consisting of 2 µl DNA template, 0.5 µl Taq DNA polymerase, 2 µl dNTP, 2 µl MgCl₂, 5 µl 10× polymerase buffer, 0.5 µl BSA, 2.5 µl of each primer, and 33 µl PCR-grade water (EurX, Gdansk, Poland). PCR amplification was run at 95 °C for 5 min of initial denaturation, followed by 32 cycles of 30 sec at 95 °C (melting), 30 sec at 56 °C (annealing), and 1 min at 72 °C (elongation), with a final elongation of 7 min at 72 °C. PCR products were cleaned using AMPure magnetic beads (Beckman Coulter, Brea, USA), and barcoded samples were normalized to 20 ng of DNA and pooled. Amplicon sequencing was performed on an Illumina MiSeq platform using 2 × 300 base pair (bp) reads at Eurofins Genomics (Ebersberg, Germany).

Read demultiplexing was performed using Cutadapt v. 3.5 (Martin 2011) with the following parameters: -e 0.2 -q 15,15 -m 150 --discard-untrimmed. Amplicon sequence variants (ASVs) were generated using trimmed reads and the DADA2 package v. 1.20 (Callahan et al. 2016) with R v. 4.1, applying the pooled approach with the following parameters: truncLen = c(220,180), maxN = 0, rm.phix = TRUE, minLen = 160. Taxonomic assignment was performed using DADA2 against the SILVA 16S rRNA SSU database release 138 (Quast et al. 2013). ASVs representing chloroplasts, mitochondria, and/or singletons were removed. Partial 16S rRNA gene sequences were aligned using SINA online v.1.2.11 (Pruesse et al. 2007).

Statistical analyses of alpha and beta diversity were conducted using the Vegan community ecology package in RStudio v. 4.3.1. Read counts were naturally low across these low-biomass, oil-rich samples, so the 16S rRNA gene dataset was rarefied to 1,500 reads, the count of the lowest sample, to retain all samples and compare them on an equal basis. This depth ensures comparability at the cost of reduced sensitivity to rare taxa. All downstream statistical analyses used this rarefied dataset. Two metrics were applied to evaluate microbial composition: relative read abundance of ASVs and α-diversity, including observed richness and Shannon index. To assess regional differences across reservoirs, we compared microbial beta-diversity using Bray-Curtis dissimilarities, visualized with Principal Coordinate Analyses (PCoA), and tested statistically with a non-parametric one-way PERMANOVA. Sample groups were defined by reservoir and sample type (e.g., Vilje RF, Vilje PW) to enable pairwise comparisons. PERMANOVA significance was evaluated using 9,999 permutations, and p-values were adjusted for multiple comparisons using the FDR method for the pairwise tests and the Bonferroni correction for the pooled RF-versus-PW comparison. PCoA revealed patterns of community dissimilarity, while the PERMANOVA assessed the statistical significance of these differences.

Beta-diversity between sample types (RF vs PW) was tested on Bray-Curtis dissimilarities (relative ASV abundance) and on Jaccard distances (presence/absence) by permutational multivariate analysis of variance (PERMANOVA; 9,999 permutations). To account for reservoir identity, a PERMANOVA was additionally run with permutations restricted within reservoir (reservoir as a permutation stratum) and as a two-way model (reservoir + sample type), partitioning variation between the two factors. Homogeneity of multivariate dispersion between RF and PW was assessed with PERMDISP (distance to the group centroid in principal-coordinate space; 9,999 permutations). Community structure was ordinated by principal coordinates analysis (PCoA) with 95% confidence ellipses per sample type. Shared and exclusive taxa were quantified at ASV and genus level, with a taxon being scored present in a sample type when detected in at least one replicate. One reservoir-fluid replicate (Vilje) that contained no reads after rarefaction was excluded, leaving 49 samples (26 RF, 23 PW).

Indicator species analysis (IndVal; Dufrêne and Legendre 1997) was used to identify ASVs indicative of PW or RF, using the multipatt function in the indicspecies package (Cáceres and Legendre 2009) with 9,999 permutations; ASVs with p < 0.05 were considered significant indicators. Analyses used the same 49 samples as the other 16S analyses (PW, n = 23; RF, n = 26). Before analysis, ASVs assigned to genera commonly reported as reagent or laboratory contaminants (Salter et al. 2014) were removed. Indicator ASVs left unassigned at the generic level by SILVA were identified by BLASTn against the NCBI nucleotide database.

Differential abundance between RF and PW was additionally assessed at the generic level with ALDEx2 (Fernandes et al. 2014). ASV counts were summed by genus; genera present in fewer than three samples or with fewer than ten total reads were removed (72 genera retained), and counts were centered-log-ratio transformed over 128 Monte Carlo Dirichlet instances. Genera were compared between RF and PW with a Wilcoxon rank-sum test and Benjamini-Hochberg correction; those with an expected BH-corrected p below 0.05 were considered differentially abundant.

DNA extracts were sent to CeGaT GmbH (Tübingen, Germany) for metagenomic sequencing. Libraries were prepared using the Nextera XT DNA Library Preparation kit (Illumina), and sequencing was performed on a NovaSeq 6000 Illumina platform at CeGaT, aiming for 50 million read pairs (2 × 150 bps). The DNA concentration from the KNB-reservoir fluid samples was insufficient for metagenomic sequencing. Therefore, only KNB’s produced water is included in the metagenomic analysis.

BBDuk (v. 37.09; Bushnell 2014) was used to remove Illumina artefacts and adapters from raw reads. Quality control and trimming were accomplished using Sickle (v. 1.33; Joshi and Fass 2011), and reads were de-duplicated with BBMap (v. 37.09; Bushnell 2014). Sourmash (v. 4.8.7; Pierce et al. 2019) was run on quality-checked forward reads with a k-mer (short sequence of length k) of size 31 to calculate Jaccard similarities between metagenomic samples. A Jaccard similarity of one indicates identical samples, while a zero means that no identical k-mers are found in the read sets. Quality-checked reads were assembled with metaSPAdes (v. 3.15.5; Nurk et al. 2017), and we then removed all scaffolds smaller than 1,000 bp from FASTA files with pullseq (v. 1.0.2; Thomas et al. 2015; github.com/bcthomas/pullseq). Prodigal (v. 2.6.3; Hyatt et al. 2010) was run in meta mode to predict and translate genes into amino acid sequences.

The 228,972 nucleotide sequences belonging to predicted genes from all metagenomes were clustered with MMSeqs2 (v. 15.6f452; Steinegger and Söding 2017) in cluster-mode 2, coverage-mode 1, with a minimum breadth of 95% and minimum sequence identity of 95%. Quality-checked reads were mapped against representative sequences from the resulting 77,195 clusters with bowtie2 (v. 2.3.5.1; Langmead and Salzberg 2012) in sensitive mode. Mappings were filtered to exclude all reads with more than five mismatches, and sequences were deemed present in samples if their minimum coverage breadth was greater than 95%.

Coverage values were normalized based on the bp counts of the forward read sequences. Only gene sequence clusters present in more than one sample were used for the PCoA calculation (to preclude the influence of under-sampling; 36,831 non-singleton clusters remained). Bray-Curtis distances were calculated from the normalized mean coverage depth using the R vegan package. MRPP and PERMANOVA (adonis2) tests were run in R to identify significant differences between sample types.

To confirm the diversity in community composition detected with 16S rRNA genes, we extracted the ribosomal protein S3 marker gene (rpS3) from the metagenomic data. RpS3 is found as a single copy within organisms, while 16S rRNA genes are often found in multiple copies (Sharon et al. 2015), making rpS3 a more reliable tool for resolving taxonomic assignments. Marker genes were identified via species-specific Hidden Markov Models (HMMs) and by comparing the amino acid sequences of predicted genes with diamond BLASTp (e-value cut-off of 1e-5; Buchfink et al. 2015) against the UniRef100 database (downloaded on 23.06.2021; The UniProt Consortium 2019). The rpS3 gene nucleotide sequences with 1,000 bp flanking regions were extracted from all samples and clustered with MMSeqs2 (v.15.6f452; Steinegger and Söding 2017) in cluster-mode 2, coverage-mode 1, with a minimum breadth of 95% and minimum sequence identity of 95%.

RpS3 sequences were taxonomically annotated by comparing them with rpS3 gene sequences extracted from the GTDB database (GTDB v. 214; Parks et al. 2022) with usearch -ublast (v. 10.0.240_i86linux64; Edgar 2010), detailing corresponding coverage, breadth, and taxonomy of extended rpS3 gene sequences (Suppl. material 1: table SS1). Quality-checked reads from all samples were mapped against representatives using bowtie2 (v. 2.3.5.1; (Langmead and Salzberg 2012) in sensitive mode. Reads with more than 10% mismatches were excluded. Mean coverage depths of extended rpS3 gene sequences were calculated for all sequences with a coverage breadth greater than 95%. ANOVA and Tukey HSD tests were run in R to identify significant differences in community composition between RF and PW samples.

As a targeted complement to the gene-cluster ordination above, we quantified the abundance of curated metabolic marker genes across the metagenomes.

All predicted protein-coding sequences were clustered using MMseqs2 (version 18.8cc5c; Steinegger and Söding 2017) with 90% nucleotide identity and a minimum 90% covered fraction in coverage mode 1 and cluster mode 2. A small mapping index was built from the resulting cluster representatives using Bowtie2 (version 2.5.5; Langmead and Salzberg 2012), and quality-controlled reads were mapped to the representative sequences in sensitive mode. The resulting SAM files were converted to BAM format and indexed using Samtools (version 1.22.1; Li et al. 2009). CoverM (version 0.7.0; Aroney et al. 2025) was used to calculate Reads Per Million (RPM; based on the TPM metric) as a quantification of cluster abundance and covered fraction, with a minimum read identity of 95% and a minimum aligned read fraction of 75%. RPM was set to zero if the covered fraction of a cluster representative was below 95% in the respective sample.

Cluster representatives were functionally annotated using three methods. First, a DIAMOND database (version 2.0.5; Buchfink et al. 2015) was generated from the Greening marker protein set, and representative sequences were searched against this database using DIAMOND blastx with an e-value threshold of 1 × 10^-5. Second, MetaCerberus (version 1.4.0; Figueroa III et al. 2024) was run in super mode using in-built Prodigal (Hyatt et al. 2010) and FragGeneScan (Rho et al. 2010) to predict amino acid sequences from the cluster representatives. Searches were performed with an e-value threshold of 1 × 10^-30 against the complete MetaCerberus HMM database (downloaded on 05 November 2025). Third, representative sequences were searched against the UniProt database (downloaded in 18 October 2024, The UniProt Consortium (2019) using DIAMOND blastx and an e-value threshold of 1 × 10^-5.

Annotations were assigned hierarchically. Hits to the Greening marker protein database were prioritized. If no Greening annotation was obtained, annotations based on MetaCerberus analyses using Prodigal predictions were used, followed by annotations based on FragGeneScan predictions. If no annotation could be assigned by these methods, the best-matching UniProt hit was reported. The table was filtered for marker genes displayed in Suppl. material 2: fig. S4. Finally, to test whether the functional composition of each pathway as a whole differed between sample types, the RPM profiles of all marker genes within a pathway were compared between RF and PW by permutational multivariate analysis of variance (PERMANOVA; adonis2 function, vegan package) on Bray-Curtis dissimilarities, using 9,999 permutations.

Length-filtered scaffolds were binned into MAGs using Abawaca (v. 1.0.0; Brown et al. 2015) and MaxBin2 (v. 2.2.7; Wu et al. 2016) with default parameters, bins were aggregated with DASTool (v. 1.1.6; Sieber et al. 2018), and the resulting selection was manually curated with uBin (v. 0.9.14; Bornemann et al. 2023). We then calculated completeness and contamination with CheckM2 (v. 1.0.1; Chklovski et al. 2023) and used GTDB-tk (v. 2.4.0; Chaumeil et al. 2022) with the Genome Taxonomy Database (v. 214; Parks et al. 2022) to assign taxonomy. The MAGs were placed into a de novo tree with GTDB-tk (v. 2.4.0; Chaumeil et al. 2022) de_novo_wf workflow, converted into an iTOLs usable format with GTDB-tk’s convert_to_itol workflow, and visualised with iTOLs (v. 6; Letunic and Bork 2024). We then applied Bowtie2 (v. 2.4.1; Langmead and Salzberg 2012) and samtools (v. 1.13; Li et al. 2009) to map all quality-checked reads against the combined MAG from all samples. Mappings were used to predict strain diversity and compare MAG presence/absence between RF and PW samples via inStrain (v. 1.8.1; Olm et al. 2021).

Coverage values of strains extracted from MAGs were normalized by the base-pair count of quality-controlled reads to correct for differences in sequencing depth between RF and PW samples. These normalized coverage values represent depth-corrected abundance estimates but were not further transformed into relative abundances (i.e., not rescaled to sum to unity within each sample). The low-abundance human contaminant Streptococcus oralis (Ren et al. 2024) was discarded from further analysis (inStrain coverage and breadth values in samples provided in Suppl. material 1: table S2). Figures displaying resulting data were constructed in R v. 4.2.2 and RStudio v. 2023.09.0+46.

We characterized the communities by 16S rRNA gene metabarcoding. Throughout, we refer to a taxon as dominant when it accounts for at least 10% of the relative read abundance (16S rRNA) or of the normalized coverage (metagenomic) within a sample. In the metabarcoding data, RF samples were dominated by taxa reported in the literature as typical of hot oil reservoirs (close to 80 °C), including Caminicella, Petrotoga, and the halophilic genus Flexistipes (Fig. 2A; Dahle et al. 2008; Vigneron et al. 2017; Christman et al. 2020). PW samples, by contrast, were dominated by marine-associated genera such as Marinobacter (Xu et al. 2008) and Malaciobacter (Kim et al. 2010).

Metagenomic profiling based on the rpS3 marker gene (11 libraries, five RF and six PW) reproduced this contrast (Suppl. material 1: table SS1). In RF, rpS3 identified taxa previously reported from reservoir fluids, including Haloanaerobium congolense (Ravot et al. 1997), Petromonas tenebris (Christman et al. 2020), Petrotoga spp. (Miranda-Tello et al. 2004), Thermoanaerobacter thermocopriae (Collins et al. 1994), and members of the Thermovirgaceae family (Dahle and Birkeland 2006). In PW, rpS3 recovered aerobic, typically seawater-associated genera, such as Marinobacter pelagius, which reached up to 98% of the 16S-inferred relative abundance in the Boa reservoir, and Malaciobacter marinus (Miller et al. 2018), which accounted for 61% of the relative abundance in KNBPW.

Taxa inferred to be potentially native to the reservoir were either absent or strongly underrepresented in PW samples. For example, species of Caminicella accounted for 86% relative abundance in RF but only 2% in PW from the Boa reservoir. Similarly, Syntrophotalea spp. accounted for 6% of RF in Vilje but were undetected in its PW. This pattern was consistent across all reservoirs (Fig. 3).

Overall, RF samples exhibited greater microbial diversity than PW, harboring a mean of 2.8 times more phyla (Suppl. material 1: table SS1). Bacillota_A (p = 0.009) and Desulfobacterota_F (p = 0.0042) were significantly more abundant in RF, whereas Pseudomonadota (p = 1.1 × 10^-11) were more abundant in PW. Shannon diversity indices derived from 16S rRNA gene sequencing data (Fig. 2B) further supported these observations, indicating higher diversity in RF. Among reservoirs, Vilje showed the highest microbial diversity, whereas Boa exhibited the lowest.

PERMANOVA analyses (Suppl. material 1: table S3) indicated significant differences across most sample groups before FDR correction, with uncorrected p-values ranging from 0.004 to 0.05, with a few exceptions. E. Kam.RF did not differ significantly from Bøyla, Vilje, or Volund RF. Similarly, Volund PW was comparable to its corresponding RF, and PW from Vilje and Boa showed no significant differences. After FDR correction, most comparisons were only near-significant (p = 0.054), likely reflecting the relatively small number of replicates (n = 4). When all PW and RF samples were grouped together, PERMANOVA detected highly significant differences (p = 0.0001), which remained robust after Bonferroni correction.

We reconstructed 26 MAGs from RF samples and 21 from PW samples (Suppl. material 2: fig. S1). PWMAGs were of higher quality than RFMAGs, with 16 of 21 exceeding 95% completeness and < 5% contamination, compared to 8 of 26 in RF (Suppl. material 2: fig. S2; Suppl. material 1: table S4); assembly statistics are given in Suppl. material 1: table S5.

Strains of anaerobic, halophilic, and thermophilic taxa, including Petromonas tenebris, members of the order Halanaerobiales, and Flexistipes spp. (Fig. 3; Suppl. material 2: fig. S1), were detected in both PW and RF samples, but were on average 10.6, 9.9, and 18.3 times more abundant in RF, respectively. In contrast, Marinobacter species were abundant in all PW samples but represented only a minor fraction of RF communities, being on average 6.4 times less abundant and undetected in the Bøyla RF sample.

While most strains occurred in at least one RF and one PW sample, five strains (Syntrophotalea, Malaciobacter, Marinobacter, and Halomonas) were detected exclusively in PW. In contrast, 16 strains, primarily belonging to Bacillota_A, Deferribacterota, and Synergistota, were found exclusively in RF.

Indicator species analysis (IndVal) reinforced this asymmetry (Suppl. material 2: fig. S3; Suppl. material 1: table S6). Of the 28 ASVs significantly associated with a sample type, 25 indicated RF, and only three indicated PW, and every PW indicator was a MarinobacterASV (IndVal up to 83%). The RF indicators were dominated by thermophilic, halophilic, and anaerobic taxa, the strongest being Flexistipes, Halomonas, Halanaerobium, Thermovirga, Dethiosulfatibacter, and Aeribacillus (IndVal > 75%). The reservoir thermophile Thermoanaerobacter was also recovered. This one-sided pattern, with many reservoir-associated indicators for RF against a single genus for PW, mirrors the differential-abundance result (Fig. 5B; Suppl. material 1: table S7).

A compositional differential-abundance analysis corroborated the indicator-species results and identified the genera driving the separation. Four genera differed significantly after FDR correction (BH < 0.05): Marinobacter, over-represented in PW, and Halanaerobium, Syntrophotalea and Flexistipes, over-represented in RF. A second tier of genera was nominally significant before correction, and all pointed toward RF (Dethiosulfatibacter, Halomonas, Thermovirga and Caminicella; Fig. 5B). Beyond Marinobacter, no genus was over-represented in PW even at the nominal level. This asymmetry, a broad set of reservoir-associated genera concentrated in RF against a single dominant genus in PW, is consistent with produced water representing a depleted and externally imprinted subset of the reservoir community.

Distinct clustering of microbial communities from RF and PW samples was confirmed by Bray-Curtis distances based on normalized 16S rRNA gene coverage (Fig. 4A), coverage of gene sequence clusters (Fig. 4B), and Jaccard similarities of quality-checked reads (Fig. 4C).

Although Jaccard’s similarity might be influenced by sequencing depth, it has been previously used to differentiate metagenomic samples from different sites (Bozzi et al. 2024). The highest intra-metagenomic Jaccard similarities were observed between Boa and Vilje PW samples, ranging from 0.14 to 0.197 (Fig. 4C). PCoA plots showed PW samples clustering together (Fig. 4A, B), indicating a high degree of similarity.

In contrast, RF samples were more widely dispersed, reflecting greater variability consistent with their higher diversity (Fig. 2B).

Sample type dominated community variation. RF and PW communities differed significantly on Bray-Curtis dissimilarities (PERMANOVA p = 0.0001), and this held when permutations were restricted within reservoir and in a two-way model, in which sample type explained more of the variation than reservoir identity (Suppl. material 1: table S8). The same separation appeared in presence/absence (Jaccard) data (p = 0.0001), so it reflects community membership, not only abundance. RF and PW differed within four of the six reservoirs (Boa, Bøyla, E. Kam. and Vilje, p < 0.05); the difference was not significant in Volund or KNB.

RF communities were also significantly more dispersed than PW (PERMDISP p = 0.002), confirming that PW is more homogeneous across reservoirs. This was not the case for presence/absence data (p = 0.20), indicating that the greater homogeneity of PW is driven by a few dominant taxa rather than by membership.

Consistent with this, RF was more even and diverse than PW (Shannon 1.9 versus 0.9, p < 0.001; Fig. 2B), while rarefied richness was comparable (66 versus 65), so the difference reflects evenness, with PW dominated by a few taxa such as Marinobacter. Only 41% of RF genera and 21% of RFASVs were recovered in the paired PW, whereas PW carried its own set of taxa absent from RF (67 genera and 928 ASVs). Produced water therefore both under-samples the reservoir community and adds a distinct set of aqueous-associated taxa.

In the gene-cluster ordination, RF samples spread along PCo1 in an arrangement broadly consistent with the geographical layout of the reservoirs (Fig. 4B). This gradient was weak: reservoir identity did not significantly structure the dataset as a whole (Suppl. material 1: table S10), and the pairwise separation reached significance only between the two most distant reservoirs, Vilje and Bøyla (MRPP delta = 0.005, A = 0.3753; PERMANOVA R² = 0.62, p = 0.03). We read this as a tentative geographic signal in RF rather than robust site resolution.

Functional marker-gene composition also differed between the two sample types. The abundance of marker genes for reservoir-relevant metabolic pathways differed significantly across all seven pathways tested (per-pathway PERMANOVA, all p < 0.05; Suppl. material 2: figs S4, S5, Suppl. material 1: table S9). These differences reflect encoded genomic potential rather than gene expression, and rest on a curated marker set. Sample type dominated both functional and taxonomic composition (R² = 0.67 and 0.59, both p < 0.001). Reservoir identity explained less and was not significant for either data type (R² = 0.22, p = 0.26 and R² = 0.29, p = 0.13; Suppl. material 1: table S10). Functional marker genes therefore separated RF from PW but did not distinguish reservoirs at this sample size. This does not mean the reservoirs are alike, only that any differences between them were too small to detect here.