All the methodological workflow (from the original KMC consortium until the genome sequencing of isolated strains) is graphically depicted as Figure 1. The original kombucha microbial culture was obtained from the collection of microorganisms of the Institute of Molecular Biology and Genetics (Kyiv, Ukraine). It was maintained in filter-sterilized black tea (Lipton, 1.2%, w/v) with white sugar (3%, w/v) (BTS) at 28°C, prior to preparation as dry organo-mineral KMC pellicle fragments samples previously described in Podolich et al. (2019). The KMC samples were placed in a three-level tray hosting four samples (independent sampling units or biological replicates) on either level (top, middle, and bottom) and mounted on an EXPOSE-2 platform outside the ISS as described in de Vera et al. (2019). These samples were maintained for 18 months (1.5 years) within the tray under a simulated Mars atmosphere (95.55% CO₂, 2.70% N₂, 1.60% Ar, 0.15% O₂, ∼370 ppm H₂O), under pressure of 980 Pa. No UV filters were used in the kombucha experiments, and space varied between 4.5 × 10⁶ and 8.4 × 10⁶ kJ/m² (de Vera et al., 2019). A set of analogous samples were kept in the laboratory at room temperature in darkness. Returned and ground reference KMC samples were cultivated for reactivation in filter-sterilized sugared (7%) black tea infusions at 28°C until a cellulose-based pellicle appeared, meaning that bacteria from the genus Komagataeibacter were active. In reference samples, pellicles appeared in 7 days, and returned samples formed pellicles after 60 days of a culturing. After KMC sample reactivation, K. oboediens IMBG185 was isolated from the exposed top-level KMC samples and opened to UV irradiation, and K. oboediens IMBG180 from initial reference KMC (iKMC) was kept in the lab during the spaceflight experiment, desiccated and stored at room temperature. Both samples were sequenced together to avoid variation due to different sequencing runs. Both strains were cultured in the HS medium (Hestrin and Schramm, 1954) at 28°C for 3 days before extraction of genomic DNA.

Total DNA from the planktonic cells of KMC was isolated by using the innuSPEED Bacteria/Fungi DNA isolation kit (Analytik Jena AG, Germany) during the stage of biofilm formation (7th day). The nucleic acids were quantified and qualified by a NanoDrop ND-1000 spectrophotometer (NanoDrop Technologies, Wilmington, DE, United States). A 450 bp library was prepared from genomic DNA with the NEBNext Fast DNA Fragmentation and Library Preparation Kit (New England Biolabs, Ipswich, MA, United States) following the manufacturer’s instructions. A library of 450 bp was prepared with the NEBNext Fast DNA Fragmentation and Library Preparation Kit (New England Biolabs). Library quality was checked with Agilent 2100 Bioanalyzer. The genomes were sequenced using Illumina HiSeq 2500 (Illumina, San Diego, CA, United States) to obtain paired-end reads. The obtained final genome coverage of IMBG185 (715.0x) and IMBG180 (1964.0x).

De novo genome assemblies of K. oboediens strains were carried out with a customized pipeline (Figure 2). Using this pipeline, the read quality was checked using FastQC v0.11.8¹, followed by seven different assembly strategies, executed with distinct software and parameters. The first assembly was carried out by Edena v3.131028 (Hernandez et al., 2008), in which the overlap-layout-consensus approach was employed to assemble the raw paired short reads solely exploiting the overlaps between them. The second assembly was executed using SPAdes v3.13 (Bankevich et al., 2012), which breaks the raw paired short reads into k-mers (substrings of size k) and applies De Bruijn graphs, which were generated to perform genome assemblies. The k-mers used for this second assembly were the best three predicted by KmerStream version 1.1 (Melsted and Halldórsson, 2014). The third assembly is a hybrid between Edena and SPAdes, where the contigs produced in the first assembly were employed as trusted contigs in a new assembly performed by SPAdes, again using the raw paired short reads and the three best k-mers detected by KmerStream. In the fourth assembly, SPAdes was executed again after quality control of the short reads performed by AdapterRemoval v2 (Schubert et al., 2016) and k-mer prediction using KmerStream. AdapterRemoval was employed to collapse paired short reads that overlap by at least 11 base pairs, remove adapters, and trim base pairs from the 3′ end with quality lower than Phred 20 and, ultimately, exclude short reads that were left with less than 31 nucleotides after trimming. Trimmed paired reads were passed as parameters to SPAdes in conjunction with the singleton and collapsed short reads. The singleton and collapsed short reads were also used in the fifth assembly, but this time, the raw paired reads were used as parameters to SPAdes besides the contigs produced by Edena, which were used as trusted contigs. The k-mers used were the top 3 predicted by KmerStream. In the sixth assembly, which was also performed by SPAdes, we employed the singleton, collapsed, and paired short reads enhanced by AdapterRemoval along with the contigs assembled by Edena, which were used as trusted contigs. The top 3 k-mers predicted by KmerStream were also used for the last assembly. The seventh assembly was performed by Unicycler v.0.4.5 and Spades v.3.9.1. Afterward, we ran QUAST (Gurevich et al., 2013) to evaluate these seven distinct assemblies for each genome, considering the N50 score, lesser numbers of contigs, and the expected genome size. Finally, the genome obtained from Unicycler v. 0.4.5 was considered as the best assembled for each of the two genomes (IMBG180 and IMBG185).

The genome sequences were annotated using Prokka v. 1.14.5 prokaryotic genome annotation (Seemann, 2014), as well as using BUSCO v. 4.0.1 to predict the completeness of single-copy orthologs (Seppey et al., 2019), and MAUVE v. 2.4.0 (Darling et al., 2004) to verify the synteny of the genes between the two genomes. The predicted contigs were also analyzed with GO FEAT (Araujo et al., 2018), an on-line platform for functional annotation and enrichment of genomic data integrated with UNIPROT², INTERPRO³, PFAM⁴, NCBI⁵, and KEGG⁶ databases.

The draft genomes of K. oboediens strains, both the space/Mars-like conditions (IMBG185) and the reference ground (IMBG180), were deposited in NCBI with the following accession numbers: K. oboediens IMBG185 (BioProject: PRJNA632016, BioSample: SAMN14942470) and K. oboediens IMBG180 (BioProject: PRJNA633374, BioSample: SAMN14942824, Accession: JABLUV000000000), respectively.

For species identification, a multilocus sequence analysis (MLSA), a small-scale phylogenetic analysis, was performed based on concatenated partial sequences of the housekeeping genes: 16S rRNA, dnaK, rpoB, and groEL, as previously described (Liu et al., 2018). These gene sequences were extracted from the strains IMBG185 and IMBG180 of deposited genomes and other 16 Komagataeibacter species (Table 1) from NCBI and aligned with MAFFT v7.310 (Katoh and Standley, 2013). Acetobacter aceti was used as the outgroup. The phylogeny was reconstructed using the maximum likelihood method implemented in IQ-TREE v1.6.12 (Nguyen et al., 2015), and the robustness of the branches was evaluated by bootstrap analysis based on 1000 iterations. For large-scale phylogenomic analysis, MUMmer (NUCmer) was used to align the genome sequences, and Pyani tool version 0.2.x was used to calculate the percentage nucleotide identity for the analysis of average nucleotide identity (ANI) and to perform the biplot dendrogram associated to a color matrix.

The prediction of genomic islands and transposases in the genomes were performed using the Genomic Islands Prediction Software GIPSy (Soares et al., 2016), using K. saccharivorans strain CV1 as a reference genome. This tool provides information about the genes present in these islands (symbiosis, pathogenicity and virulence, antibiotic resistance, and metabolic) and their respective location in the genomes. Only genomic islands with normal and strong classifications were considered. In addition, GIPSy also detects transposase genes based on HMMER3 (Eddy, 2011) with E-value ≤ 1 × 10^–4. The online tool PAPA (Prion aggregation prediction algorithm) (Toombs et al., 2012) was used to predict candidate prion aggregating regions in the sequences annotated as “Transcription termination factor Rho” in both genomes.

Prophage Hunter tool (Song et al., 2019) was used for prediction and classification of prophages. Circular contigs with high depth coverage, predicted by Unicycler (Wick et al., 2017), were considered plasmids. Non-circular contigs were assessed using a combination of three different tools: (i) Mob-recon (reference-based approach, version 2.1.0) (Robertson and Nash, 2018), (ii) HyAsP (hybrid approach, version 1.0.0) (Müller and Chauve, 2019), and (iii) Gplas (version 0.5.0) (Arredondo-Alonso et al., 2020). Non-circular contigs were considered plasmids only when predicted by the majority of these three tools.

The prediction of Cas proteins and classification of CRISPR-Cas system subtype in both genomes were performed via CRISPRCas Identifier (Padilha et al., 2020).

The KAAS web server (Moriya et al., 2007) was used to predict KEGG Orthology (KO) identifiers for both the sequenced strains. The representative set of organisms was selected as the default for prokaryotes, in addition to the closely related species Gluconacetobacter diazotrophicus (three-letter code: gdi), Acetobacter pasteurianus (apt), and Granulibacter bethesdensis (gbe). The bidirectional best-hit approach was used as an assignment method. In order to complement KEGG analysis, the MetaCyc (Caspi et al., 2020) ontology was also employed to perform metabolic prediction. For this analysis, the replicons and associated annotation for both genomes were input to Pathway Tools (Karp et al., 2002), which creates a Pathway/Genome Database (PGDB) containing the predicted metabolic pathways for the organisms. Metabolic reconstruction included determining gene-protein-reaction associations, which are primarily based on the corresponding enzyme commission (EC) numbers. The pathway score yielded by Pathway Tools was used as a proxy for the presence or absence of a given pathway.

Pathway Tools was also employed to thoroughly study the key cellulose synthesis bcs operon, which was mined in K. oboediens using the bcs operon identified in K. xylinus as template (PMID: 31983038; PMID: 30761107). BLAST analyses against the UniProt database (PMID: 15215342) were also performed to search for orthologs related to this metabolic function.

Stress-related genes and pathways were identified by manual literature searches as well as by the retrieval of curated proteins deposited in UniProtKB/Swiss-Prot, which were associated to the “Stress response” keyword (KW-0346) and belonged to either Archaea or Bacteria domains. In order to reduce the redundancy of this protein set (originally composed of 7,724 sequences), the proteins were clustered at 90% sequence identity using UniRef90 clusters, yielding a final set of 2,687 stress-related sequences. Subsequently, for retrieving orthologs in K. oboediens, BLASTp was performed using each predicted K. oboediens proteome set as a query against the reference stress-related proteins, followed by manual curation of the alignments.

A principal component analysis (PCA) was also conducted to assess the metabolic profiles of K. oboediens and related bacteria. In order to perform this task, the enzymatic content (as predicted by KEGG) was used as input to the prcomp function in R⁷, which employs the singular value decomposition method to extract principal components. The compared bacteria (and corresponding labels, when applicable) were Acetobacter aceti (Aac), A. cerevisae (Ace), A. malorum (Ama), A. pasteurianus (Apa), A. pomorum (Apo), Clostridium kluyveri, Escherichia coli, E. albertii, E. fergusonii, Komagataeibacter cocoi (Kco), K. diospyri (Kdi), K. europaeus (Keu), K. hansenii (Kha), K. intermedius (Kin), K. kakiaceti (Kka), K. maltaceti (Kma), K. medellinensis (Kme), K. nataicola (Kna), K. oboediens LMG 18849 (Kob 18849), K. rhaeticus (Krh), K. saccharivorans (Ksa), K. sucrofermentais (Ksu), K. swingsii (Ksw), K. kylinus (Kxy), Salmonella bongori, S. enterica, Shigella boydii, Sh. flexneri, and Sh. sonnei.

A comparison between coding sequences (CDS), as well as corresponding plasmid CDS, of K. oboediens in space and terrestrial conditions was performed to search for SNV between their genomes and their plasmids. An in-house Python script (see Availability of data and materials) was used to find the corresponding genes in both genomes based on the following criteria: (i) both genes had to present exactly the same 60 first nucleotides and (ii) have the same gene length. These criteria were applied assuming that the space stress would not have caused severe mutations that would make these criteria to fail the analysis, as such stress would possibly be lethal to the bacteria.

Furthermore, a list of corresponding pairs of genes present in both genomes was generated, as well as a list of genes without a match in the other genome. The list of genome-specific, putative exclusive genes was then further manually curated. The lists generated are available in Supplementary Material 1.

In order to evaluate the K. oboediens pyruvate phosphate dikinase (ppdk) single mutation amino acid, we performed ppdk homology modeling and molecular dynamics for discovering the influence of the amino acid SER287 on the ATP/ADP site. Firstly, the K. oboediens PPDK amino acid FASTA sequence was submitted to The Swiss Model Workspace (Arnold et al., 2006) for automated modeling using the PDB template 5LU4⁸, because it is complexed with the ADP molecule, as well as Mg⁺⁺ cofactor. Additionally, we have considered the template resolution below 3.0 Å, identity, and sequence coverage above 50%. After QMEAN validation, we submitted both the ppdk model and its template to a structural alignment using the UCSF Chimera (Pettersen et al., 2004) program to identify the catalytic and other important binding site positions, as well as to verify if the ppdk mutation could affect the enzyme function. Considering the ATP/ADP binding site’s high identity between both the ppdk model and its 5LU4 template, and, since ppdk SER287 was aligned with 5LU4 GLY285, we decided to generate a single mutation in 5LU4 (GLY285 to SER285) in order to use this crystallographic structure complexed with ADP and Mg⁺⁺ ligands and performed a molecular dynamics simulation for describing how this mutation would affect the ATP/ADP binding site, which would be important for enzyme activation and its transition state (Minges et al., 2017b). In this case, using a homologous crystallographic structure avoids a possible bias, which could happen only if a modeled protein docked with its respective ATP/ADP ligand was used.

For the Molecular Dynamic experiments, we used the GROMACS 2020 program. The mutated complex 5LU4, bound to ADP and Mg⁺⁺, was subjected to a stability MD in water solvent for 50 ns. Subsequently, we followed the GROMACS tutorial for protein–ligand complexes (Lemkul, 2019) with adaptations, as well as with the program manual (Berendsen et al., 1995) for all the simulation steps and analyses. The MD graphs were generated using the Grace 5.1.25 program⁹, and the output data were used for creating a video of the simulation using the VMD program (Humphrey et al., 1996).

A total of two distinct and highly supported clades (bootstrap values: 99–100%) were retrieved in multigene phylogenetic analysis of the genus Komagataeibacter: one of them (clade 1), including the well-studied and representative species, K. hansenii, and the other one (clade 2) comprising another much studied and representative species, K. xylinus. Clade 2, in turn, encompassed two differently and highly supported clades: clade 3, with 76% bootstrap, including, besides K. xylinus and other related species, K. medellinensis and K. saccharivorans, and clade 4, encompassing K. oboediens, K. intermedius and related species (Figure 3A). Moreover, the analysis clearly revealed that both strains (IMBG180 and IMBG185) are very closely related to K. oboediens (Figure 3A), thus confirming that they belong to the same species K. oboediens. Furthermore, this finding of the small-scale phylogenetic analysis (MLSA) was fully corroborated with large-scale phylogenomic analysis based on global genomic similarity (ANI), showing that both strains (IMBG180 and IMBG185) and all the K. oboediens available genomes are the most identical (the large red cluster) amongst all the other genomes of Komagataeibacter species (Figure 3B). Moreover, ANI (average nucleotide identity) index between K. oboediens IMBG180 and IMBG185 were 99.9%, suggesting that they are almost identical.

A remarkable similarity between IMBG185 (space/Mars-like conditions) and IMBG180 (ground reference) strains of sequenced genomes (Figure 4 and Supplementary Material 2) was noticed. Very small differences in size (3,560,605 and 3,562,447 bp), total gene prediction (3,282 and 3,275), total number of protein-coding genes (3,230 and 3,222), and other features between samples IMBG185 and IMBG180, respectively, may be because they are not of the same strain or because of minor changes in IMBG185 due to the stressful conditions it was exposed to. Besides, both strains displayed exactly the same G + C content (61.85%) and the same number of rRNA genes (3) (Table 2).

In the genome completeness analysis, 432 orthologous genes were found in the two analyzed genomes (IMBG185 and IMBG180), according to the following annotation databases: tBLASTn 2.2+, Augustus 3.2, Prodigal, and HMMER3.1. The genomes also showed 427 complete orthologous genes divided by categories: singles-copy (425), duplicated (2), fragmented (1), and missing (1). These results revealed the very high similarity between IMBG185 and IMBG180 sequenced genomes (Supplementary Material 2).

The number of symbiosis (4), metabolic (4), and virulence + pathogenicity islands (7), and antibiotic resistance (6) islands found in both IMBGs were identical. Furthermore, considering only the identifications with identity ≥ 30% and E-values ≤ 1 × 10^–7, the number of protein-coding genes inside all the genomic islands were also similar (Figure 4 and Supplementary Material 3). Moreover, the prediction of transposases retrieved 39 potential protein sequences for both IMBGs. For Prion-like proteins, the transcription terminator Rho gene had been detected in both IMBGs, but the Cb-Rho cPrD (prion-forming) domain was not found in the IMBG180 and IMBG185 genomes. This may mean that Cb-Rho does not behave as a prion in IMBGs.

A total of 96 candidates for prophage regions in the IMBGs were predicted, and these regions have 92 inactive, three active, and one ambiguous prophage (Figure 3 and Supplementary Material 4). Regarding the regions containing active prophages, 62 genes found in IMBGs displayed high similarity to Rhodovulum phage RS1 (score = 0.87), and 52 genes from postflight (IMBG185) and 60 genes from reference (IMBG180) genomes showed high similarity to Ruegeria phage DSS3-P1 (score = 0.97). Moreover, 18 genes from the genome IMBG185 contained similar sequences to Streptomyces phage (score = 0.92), whereas in the IMBG180 genome, 19 prophage genes in unrelated regions with other bacterial species were also characterized.

Plasmid prediction of the assembled contigs of K. oboediens IMBGs genomes resulted in two plasmids (1 and 2) on each strain, with a variation of only 20 nucleotides in one of the contigs of plasmid 2. Although only plasmid 1 was confirmed as circular, plasmid 2, which contains two contigs, was characterized as a mob conjugative plasmid (Table 3).

Komagataeibacter oboediens IMBGs strains genome exhibited two CRISPR-Cas cassettes, which were classified as CAS-VI-B and one as CAS-III-D (or Type I, according to Makarova et al., 2020) subtype. The only difference between the isolates was the position of some cas genes along the genome sequence (Supplementary Material 4).

The comparison of IMBG genomes resulted in a list of 3,234 corresponding genes present in both samples. This represents 98.5% (3,234/3,282) of IMBG185 (space/Mars-like) and 98.7% (3,234/3,275) of IMBG180 (ground reference) coding sequences. The IMBG185 genome exhibited 48 putative genome-specific genes and the IMBG180 genome displayed 41 putative specific genes. Furthermore, among the corresponding genes, only two of these genes presented SNVs: pyruvate phosphate dikinase (ppdk) (1 non-synonymous SNV), and sugar phosphatase (yidA) (3 synonymous SNVs) genes (Table 4). The plasmid CDS SNV analysis did not find any variation between the strains, indicating that the CDS plasmid sequences of both strains were identical. The non-synonymous substitution SER285 (transition) that was identified in the ppdk gene of IMBG185 was further analyzed with protein modeling and molecular dynamics in order to evaluate if this SNV could affect the predicted binding capacity of the ADP molecule and pyruvate metabolism.

The K. oboediens PPDK 3D model was constructed and structurally aligned with its respective template (5LU4). The region where the mutation occurred is well-aligned with its template corresponding amino acids (5LU4 ALA 255 to HIS 390), which are involved with binding to ADP and Mg⁺⁺ ligands (Minges et al., 2017a,b; Figure 5). RMSD between K. oboediens PPDK chain B and 5LU4 chain A was 0.93 Å with overall identity of 57%. Figure 6 clearly shows the structural alignment between K. oboediens PPDK (magenta) and 5LU4 (dark gray) with the amino acids SER287 and GLY285 (respectively on nearby positions) as they appear surrounding the ADP/MG⁺⁺ binding site. Additionally, the amino acids responsible for binding to these ligands in the ADP/MG⁺⁺ binding site are in the same positions between the two proteins.

After 50 ns of solvent-accessible molecular dynamics of the complex 5LU4 mutated G285S, the behavior of the ADP ligand and the influence of serine on the nucleotide-binding positioning were evaluated. The root mean square fluctuations (RMSF) of all the chain B amino acids were analyzed (Figure 7A), displaying that the amino acid S285 had an RMSF of 0.3 nm, which indicated that its sidechain possibly had freely rotational movements. This observation was further confirmed by the MD video (see Supplementary Material 5¹⁰), showing the rotational movements of this amino acid. Additionally, the binding stability of the ADP molecule to the protein backbone (Figure 7B) was also evaluated, revealing that the root mean square deviation of this ligand was below 0.4 nm, which indicates that the ligand is stable inside its binding site. Furthermore, the MD video also confirmed the ADP positioning inside its 5LU4 binding site for all the 50 ns time of the simulation. The Coulombic interaction energy between ADP and 5LU4 (Figure 7C) showed that they had a high affinity with each other, with an average of −191.586 KJ/mol. Therefore, it is unlikely that they have been affected by the mutated SER285 amino acid. Furthermore, ADP formed a high number of hydrogen bonds (Figure 7D) with the enzyme active amino acids, varying from 1 to 8, during all the simulation. Altogether, it is unlikely that the mutated SER285 can affect the binding capacity of the ADP molecule or the enzyme metabolism.

The genome-scale metabolic reconstructions of K. oboediens species comprised 1,894 enzymatic reactions associated with 1,111 genes and 1,418 metabolites (small-molecules in the MetaCyc ontology), covering 33.9% of the annotated open reading frames (Supplementary Material 6). Furthermore, we evaluated differences of the metabolic repertoire between the two isolates of K. oboediens IMBG180 and IMBG185 using KEGG and MetaCyc ontologies for metabolic reconstructions.

The genome-based analysis revealed the absence of gene content differences that could impact the metabolic capacities of post-flight K. oboediens IMBG185, since the predicted content of metabolic genes of the ground reference (IMBG180) and strains exposed to space/Mars-like conditions (IMBG185) was quite similar, independent of the metabolic ontology applied (KEGG or MetaCyc) (Supplementary Material 6). Subsequently, we compared the metabolic functions within the Komagataeibacter genus, as well as how it compares to other bacteria. As expected, there is extensive conservation, in terms of metabolic capacities, within the genus (Figure 8). Furthermore, bacteria with similar metabolic profiles are also closely related to Komagataeibacter, particularly the acetic acid bacteria within the Acetobacter genus (Figure 8). We also included in these analyses representatives from Enterobacteriaceae, which clustered together in a single group but have strong metabolic dissimilarities when compared to the acetic acid bacteria (Figure 7), probably reflecting different lifestyles of these organisms. It is worth noting that K. kakiaceti was the only Komagataeibacter representative that did not group with the remaining bacteria of this genus (Figure 8). Next, we evaluated pathways of interest coded in the genomes of K. oboediens IMBGs that could be important to the overall functioning of the KMC.

The main pathways found in IMBGs as glycolysis, pyruvate oxidation, Krebs cycle, pentose phosphate pathway (oxidative and non-oxidative phases), phosphoribosyl diphosphate (PRPP – an important intermediate in cellular metabolism), and Entner-Doudoroff pathway had their KOs completely matched (100%) with the reference pathway of KEGG. Six of seven of the gluconeogenesis potential proteins were detected (86%), and only one protein was not complete (86%) considering the pentose phosphate pathway (7 proteins). Other carbohydrate pathways found were the De Ley-Doudoroff pathway (D-galactonate degradation), nucleotide sugar biosynthesis, and UDP-N-acetyl-D-glucosamine biosynthesis (Supplementary Material 7).

For both IMBGs, several fermentative complete pathways were identified through the genome potential enzymes. Key enzymes of the fermenting processes were found, such as phosphate acetyltransferase (pta) and acetate kinase (ackA) for acetate production, pyruvate decarboxylase (pdc) and alcohol dehydrogenase (adh) for ethanol production, acetolactate synthase (ilv) for 2-acetolactate production, α-acetolactate decarboxylase (budA) and diacetyl reductase (budC) for acetoin production, and finally the enzyme diacetyl reductase (budC) to complete the acid homeostasis in cells (Supplementary Material 7).

The putative presence of genes related to the fatty acid biosynthesis (FAB) route for MCCA synthesis such as acetyl-CoA carboxylase (accABCD) (all the subunits) and ortholog coding for malonyltransferase (fabD) (at least one), ketoacyl-ACP synthase (fabF; fabH), ketoacyl-ACP reductase (fabG), hydroxyacyl-ACP dehydratase (fadZ), and enoyl-ACP reductase (fabI; fabK) enzymes was detected (Supplementary Material 8).

In order to disclose general stress-related genes coded in IMBGs genomes, we firstly conducted a large-scale analysis by leveraging information on all the stress-associated, manually reviewed proteins available in UniProtKB/Swiss-Prot and compared both sequenced strains. This resulted in the identification of highly conserved genes in IMBGs genomes (Supplementary Material 9), and no differences in gene content between both IMBGs were observed in this analysis. Genes with chaperone activities, such as groL, groS, dnaJ, dnaK, and hfq were detected, as well as the products of recA and uvrABC systems involved against DNA damage and other resistance properties. Protease activities coded by hslVU, lon, and clp were detected in the genomes, and also hrcA and lexA, both of which regulate the expression of part of these genes. Trehalose biosynthesis pathway (otsA/otsB) genes associated with ROS scavenging (katB, oxyR) were also identified in both genomes (Supplementary Material 9).

The cellulose synthesis key bcs operons following the K. xylinus (a close member of K. oboediens) as described by Liu et al. (2018), Gullo et al. (2019) were constructed and can be accessed on Gullo et al. (2019), La China et al. (2020). There was no topological change in any of the bcs operon genes between postflight and reference genomes, and all the genes of the four operons had 100% identity at both the DNA and amino acid levels, as recently characterized by Lee et al. (2021).

A total of five protein-coding genes associated with biological nitrogen fixation was retrieved in both genomes of IMBG185 (space/Mars-like conditions) and IMBG180 (ground reference) strains (Table 5). One gene (nifU) product is involved in the assembly and incorporation of iron into the metalloclusters of a nitrogenase, whereas the other four genes either encode transcription factors that regulate (positively and negatively) nitrogenase (such as fixK and vnfA) or code for proteins that regulate the transcription factors of nitrogenase (such as the PII proteins glnB and glnK).

There are five enzyme-coding genes directly associated with the biosynthesis of hopanoid lipids in both genomes (Table 6). One of the genes (sqc) codes for the enzyme responsible for the synthesis of hopanoid lipids, whereas another gene (hpnR) encodes an enzyme that modifies these hopanoids. Moreover, the other three genes (hpnC, hpnD, and hpnE) code for enzymes that jointly participate in the conversion of all trans-farnesyl diphosphate to squalene.