Transcriptome Analyses of Prophage in Mediating Persistent Methicillin-Resistant Staphylococcus aureus Endovascular Infection

Persistent methicillin-resistant Staphylococcus aureus (MRSA) endovascular infections represent a significant subset of S. aureus infections and correlate with exceptionally high mortality. We have recently demonstrated that the lysogenization of prophage ϕSA169 from a clinical persistent MRSA bacteremia isolate (300-169) into a clinical resolving bacteremia MRSA isolate (301-188) resulted in the acquisition of well-defined in vitro and in vivo phenotypic and genotypic profiles related to persistent outcome. However, the underlying mechanism(s) of this impact is unknown. In the current study, we explored the genetic mechanism that may contribute to the ϕSA169-correlated persistence using RNA sequencing. Transcriptomic analyses revealed that the most significant impacts of ϕSA169 were: (i) the enhancement of fatty acid biosynthesis and purine and pyrimidine metabolic pathways; (ii) the repression of galactose metabolism and phosphotransferase system (PTS); and (iii) the down-regulation of the mutual prophage genes in both 300-169 and 301-188 strains. In addition, the influence of different genetic backgrounds between 300-169 and 301-188 might also be involved in the persistent outcome. These findings may provide targets for future studies on the persistence of MRSA.


Introduction
Methicillin-resistant S. aureus (MRSA) is a major cause of life-threatening endovascular infections, including bacteremia and infective endocarditis (IE) [1,2]. Persistent MRSA bacteremia (PB; defined as ≥5 days of positive blood cultures in the presence of antibiotic therapy) represents~15 to 30% of such infections [3,4]. In addition, it is very worrisome that most PB isolates appear to be susceptible in vitro to gold-standard anti-MRSA antibiotics (e.g., vancomycin (VAN) and daptomycin (DAP)) by the Clinical and Laboratory Standards Institute (CLSI) breakpoints [4][5][6], yet persistent in vivo. Thus, PB represents a uniquely vital variant of traditional antibiotic resistance mechanisms. This problem underscores an urgent need to understand the mechanism(s) of specific factors driving this syndrome.
The current study aimed to define the impact of φSA169 on genetic factors which may contribute to the PB phenotypes by RNA sequencing (RNA-seq) using PB 300-169 wild type (WT), RB 301-188 WT, and φSA169 lysogenized RB 301-188 (301-188::φSA169) strains. The transcriptomic analyses emphasized genetic factors that might contribute to the PB outcomes and provided clues for future studies on molecular mechanisms of PB outcomes.

Bacterial Strains, Plasmids, and Growth Medium
Three MRSA strains, including PB 300-169 WT (300-169), RB 301-188 WT (301-188), and 301-188 WT φSA169 lysogenization (301-188::φSA169), were used in our previous [11] and current studies. The PB 300-169 strain was isolated from a patient with 16 days of persistent MRSA bacteremia, while the RB 301-188 strain was obtained from a patient with 2 days of MRSA bacteremia [4]. In addition, all the three study strains have a minimum inhibitory concentration (MIC) to VAN of 0.5 µg/mL and are susceptible to VAN in vitro based upon the CLSI breakpoints [11]. The strains were routinely grown at 37 • C in tryptic soy broth (TSB; Becton Dickinson and Company, NJ, USA) or on tryptic soy agar (TSA) plates if not otherwise specified.

RNA Isolation
RNA isolation was performed following the method described in previous studies [13,14]. In brief, overnight cultured cells of the study strains were pelleted by centrifugation and resuspended in Buffer RLT from RNeasy kit (Qiagen, Germantown, MD, USA), and then transferred into lysing matrix B (MP Biomedicals, Irvine, CA, USA) containing 0.1 mm silica spheres for mechanical lysis using Fastprep (Thermo Fisher, Waltham, MA, USA). Total RNA was isolated according to the manufacturer's instructions of the RNeasy kit. DNA in the samples was removed using a TURBO TM DNase kit (Thermo Fisher, Waltham, MA, USA) [11]. Biological duplicates from two different experiments were prepared for each study strain. RNA samples with concentrations ≥ 100 ng/µL and 260/280 ratio between 1.9 and 2.0 were submitted to the Novogene Corporation Inc (Sacramento, CA, USA) for RNA-seq.

RNA-Seq and Data Analyses
RNA degradation, purity, integrity, and quantitation were checked prior to the RNAseq. RNA-seq libraries were constructed using NEBNext ® Ultra TM RNA Library Prep Kit for Illumina ® (NEB, Ipswich, WA, USA). The index-coded samples were clustered using the PE Cluster Kit cBot-HS (Illumina, San Diego, CA, USA) on a cBot Cluster Generation System. Then, the samples were sequenced, and paired-end reads were obtained. For data analyses, RNA-seq reads were mapped to the genome of the PB 300-169 strain (Accession: JASL00000000) [12] using Bowtie2 [15]. Analyses of differential expressions between any two study strains (two biological replicates per study strain) were performed using DESeq2 R package based on a negative binomial distribution. The resulting p values were adjusted using Benjamini and Hochberg's approach for controlling the false discovery rate. The genes with an adjusted p value (p adj) ≤ 0.05 and |log 2 (fold change)| > 0 were defined as differentially expressed genes (DEGs), indicating the genes had significantly different expression levels in the two strains comparison. The DEGs list generated from the comparison of transcriptomic profiles between the isogenic strain set (301-188 and 301-188::φSA169) indicated the impact of φSA169. In addition, comparisons of 300-169 vs. 301-188 and 300-169 vs. 301-188::φSA169 were also performed to further investigate the role of the distinct genetic backgrounds on the transcriptional changes. The DEGs were Genes 2022, 13, 1527 3 of 13 classified using the Kyoto Encyclopedia of Genes and Genomes (KEGG) mapper tool with the ST45 mode strain of MRSA CA-347 [16].

Verification of RNA-Seq Results by qRT-PCR
The expression levels of selected genes from the DEGs listed above were confirmed by qRT-PCR as described previously [11,17,18]. The expression of gyrB was used as a wellstudied host gene to normalize transcripts levels, and relative expression was calculated by the ∆∆C T method [5]. The relative expression level was then used to calculate the fold changes in the selected genes in strain comparisons.

Global Analyses of Gene Expression
Each sample yielded a high percentage of exon-mapped reads (85.3-90.1%) that covered over 2000 genes, indicating the abundance of mRNA and low interference from non-coding RNAs. More than 86% of the mapped genes had at least one fragment per kilobase of transcript sequence per million (FPKM), suggesting that the transcriptional profiles covered most of the genes in the study strains. Principal component analysis (PCA) was performed to assess the overall differences in the gene expression of the study strains ( Figure 1). Among the study strains, 300-169 had a different genetic background vs. 301-188, while 301-188 and 301-188::φSA169 were isogenic strain-set with the only difference in the absence/presence of φSA169. The strains 301-169 and 301-188 had the most distant locations on the PCA biplot, indicating the most significant genetic variation, while 301-188 and 301-188::φSA169 had the closest locations suggesting minor variation, which might be due to the same genetic background (Figure 1). defined as differentially expressed genes (DEGs), indicating the genes had significantly different expression levels in the two strains comparison. The DEGs list generated from the comparison of transcriptomic profiles between the isogenic strain set (301-188 and 301-188::φSA169) indicated the impact of φSA169. In addition, comparisons of 300-169 vs. 301-188 and 300-169 vs. 301-188::φSA169 were also performed to further investigate the role of the distinct genetic backgrounds on the transcriptional changes. The DEGs were classified using the Kyoto Encyclopedia of Genes and Genomes (KEGG) mapper tool with the ST45 mode strain of MRSA CA-347 [16].

Verification of RNA-Seq Results by qRT-PCR
The expression levels of selected genes from the DEGs listed above were confirmed by qRT-PCR as described previously [11,17,18]. The expression of gyrB was used as a wellstudied host gene to normalize transcripts levels, and relative expression was calculated by the ΔΔCT method [5]. The relative expression level was then used to calculate the fold changes in the selected genes in strain comparisons.

Global Analyses of Gene Expression
Each sample yielded a high percentage of exon-mapped reads (85.3-90.1%) that covered over 2000 genes, indicating the abundance of mRNA and low interference from noncoding RNAs. More than 86% of the mapped genes had at least one fragment per kilobase of transcript sequence per million (FPKM), suggesting that the transcriptional profiles covered most of the genes in the study strains. Principal component analysis (PCA) was performed to assess the overall differences in the gene expression of the study strains ( Figure 1). Among the study strains, 300-169 had a different genetic background vs. 301-188, while 301-188 and 301-188::φSA169 were isogenic strain-set with the only difference in the absence/presence of φSA169. The strains 301-169 and 301-188 had the most distant locations on the PCA biplot, indicating the most significant genetic variation, while 301-188 and 301-188::φSA169 had the closest locations suggesting minor variation, which might be due to the same genetic background ( Figure 1). :φSA169 strains. The X-axis represents the first principal component (PC1) that displays the maximum variation through the data, while the Y-axis represents the second principal component (PC2) that displays the next highest variation. Each dot represents a biological duplicate of a study strain. The replicates of each study strain were well clustered, indicating the good reproducibility of the samples in each strain. RNA-seq results of the study strains were scattered in the graph, which indicated the significant genetic variations of the study strains. The X-axis represents the first principal component (PC1) that displays the maximum variation through the data, while the Y-axis represents the second principal component (PC2) that displays the next highest variation. Each dot represents a biological duplicate of a study strain. The replicates of each study strain were well clustered, indicating the good reproducibility of the samples in each strain. RNA-seq results of the study strains were scattered in the graph, which indicated the significant genetic variations of the study strains.
The transcriptome profiles of each study strain were compared to identify the DEGs ( Figure 2, Table 1). There were 153 DEGs in 301-188::φSA169 vs. 301-188 (Figure 2a), while over 1200 DEGs were found in 300-169 vs. 301-188 ( Figure 2b) and 300-169 vs. 301-188::φSA169 ( Figure 2c). In the strain 301-188::φSA169, 77 and 76 DEGs were significantly up-and down-regulated, respectively, compared to the parental 301-188 (Table 1). Over half of the up-regulated DEGs (49 out of 77) were the genes of φSA169 (Table S1), while more than one-third of the down-regulated DEGs (24 out of 76) belonged to the mutual prophage in both 300-169 and 301-188 (Table S2). The high log2(fold change) values of the 49 φSA169 genes (Table S1) indicated the absence in 301-188. In the 300-169 strain, 666 and 633 DEGs were significantly up-and down-regulated, respectively, compared to 301-188 ( Table 1). The detailed up-and down-regulated DEGs in the comparison of 300-169 vs. 301-188 are presented in Tables S3 and S4, respectively. In the comparison of 300-169 vs. 301-188::φSA169, a total of 637 and 613 DEGs were significantly up-and down-regulated, respectively ( Table 1). The detailed up-and down-regulated DEGs are presented in Tables S5 and S6, respectively.  Prophage φSA169 was initially identified in PB 300-169 and transduced into RB 301-188 to construct the 301-188::φSA169 strain. Therefore, φSA169 was an exogenous genomic element for the 301-188 chromosome despite the similar genetic background between 300-169 and 301-188 (e.g., CC45, agr I, and SCCmec IV); thus, the gene expression of φSA169 may differ in the 300-169 vs. 301-188::φSA169. There were 58 out of a total of 67 annotated genes in φSA169 detected in the current RNA-seq results ( Figure 3). The plotted expression levels of φSA169 genes in both 300-169 and 301-188::φSA169 are presented in Figure  3. Bacteriophage (phage) genes are highly mosaic and grouped into different modules based on the functions of the gene products (18). In general, φSA169 genes in the modules of lysogeny, packing and morphogenesis, and lysis were highly expressed, while genes in the replication module had low expression ( Figure 3). In addition, the transcriptional profiles of φSA169 were similar in both strains. However, some φSA169 genes, especially in the packing and morphogenesis module, had different expression levels in the two strains, which might imply the impact of the distinct genetic backgrounds.  Prophage φSA169 was initially identified in PB 300-169 and transduced into RB 301-188 to construct the 301-188::φSA169 strain. Therefore, φSA169 was an exogenous genomic element for the 301-188 chromosome despite the similar genetic background between 300-169 and 301-188 (e.g., CC45, agr I, and SCCmec IV); thus, the gene expression of φSA169 may differ in the 300-169 vs. 301-188::φSA169. There were 58 out of a total of 67 annotated genes in φSA169 detected in the current RNA-seq results ( Figure 3). The plotted expression levels of φSA169 genes in both 300-169 and 301-188::φSA169 are presented in Figure 3. Bacteriophage (phage) genes are highly mosaic and grouped into different modules based on the functions of the gene products (18). In general, φSA169 genes in the modules of lysogeny, packing and morphogenesis, and lysis were highly expressed, while genes in the replication module had low expression ( Figure 3). In addition, the transcriptional profiles of φSA169 were similar in both strains. However, some φSA169 genes, especially in the packing and morphogenesis module, had different expression levels in the two strains, which might imply the impact of the distinct genetic backgrounds.

The Impact of φSA169 on Transcriptional Profiles
The 301-188::φSA169 and 301-188 formed an isogenic strain set; thus, the DEGs from the comparison of the two strains were likely caused by φSA169. On the other hand, 300-169 and 301-188 strains had distinct genetic backgrounds; thus, the DEGs profile of these two strains might be affected by both φSA169 and their genetic backgrounds. Therefore, the overlapping DEGs between the two comparisons (301-188::φSA169 vs. 301-188 and 300-169 vs. 301-188) might indicate the specific impact of φSA169. There  (Table 2). Over half of the down-regulated DEGs (24 out of 45) belonged to the mutual prophage in both 300-169 and 301-188 strains, and the remaining 21 DEGs were the MRSA host genes, including lacABCDEF, treP, and pfkB (Table 3). The transcriptional profiles of φSA169 in 300-169 and 301-188::φSA169 were similar, and genes from the packing and morphogenesis module had higher expression levels than genes from other modules. Expression levels of some φSA169 genes, especially the genes from packing and morphogenesis, were significantly higher in 300-169 compared to 301-188::φSA169. * p adj < 0.05, ** p adj < 0.01, *** p adj < 0.001

The Impact of φSA169 on Transcriptional Profiles
The 301-188::φSA169 and 301-188 formed an isogenic strain set; thus, the DEGs from the comparison of the two strains were likely caused by φSA169. On the other hand, 300-169 and 301-188 strains had distinct genetic backgrounds; thus, the DEGs profile of these two strains might be affected by both φSA169 and their genetic backgrounds. Therefore, the overlapping DEGs between the two comparisons (301-188::φSA169 vs. 301-188 and 300-169 vs. 301-188) might indicate the specific impact of φSA169. There were total of 65

DEGs Impacted by Both φSA169 and MRSA Genetic Backgrounds
There were 36 (Figure 4a, Table S9) and 23 (Figure 4b, Table S10) DEGs up-and downregulated in all three comparisons, respectively. It indicated that these DEGs were affected by both φSA169 and the genetic backgrounds of 300-169 and 301-188. The up-regulated DEGs included 26 genes in φSA169 and 10 other staphylococcal genes (Table S9). The down-regulated DEGs consisted of 10 genes in the mutual prophage and 13 MRSA genes (Table S10).

Global KEGG Analyses of DEG Profiles
To understand the gene functions and pathways associated with the persistent outcomes, we classified the DEGs using the KEGG pathways mapper tool ( Figure 5). In 301-188::φSA169, a significant number of genes were down-regulated compared to 301-188 (e.g., carbohydrate metabolism and membrane transport; Figure 5a). In 300-169, genes involved in carbohydrate and amino acids metabolisms, metabolism of cofactors and vitamins, and membrane transport were mainly differentially expressed vs. 301-188 (Figure 5b). Some pathways were found up-regulated in 300-169 vs. 301-188 (e.g., glycan biosynthesis and metabolism, transcription, and drug resistance; Figure 5b). The KEGG analysis profile of 300-169 vs. 301-188::φSA169 (Figure 5c) was similar to 300-169 vs. 301-188 (Figure 5b), suggesting the significant differences may be due to the different genetic backgrounds.

φSA169-Specific KEGG Analyses
The overlapping DEGs of 301-188::φSA169 vs. 301-188 and 300-169 vs. 301-188 might represent the genes regulated explicitly by φSA169 (Figure 4). The KEGG profile of the overlapping DEGs indicated that most of these genes were involved in metabolic pathways ( Figure 6). For instance, the DEGs of fatty acid biosynthesis (fabFH), purine metabolism (purA), and RNA degradation (AS94_08925) were up-regulated by φSA169. Among the down-regulated DEGs by φSA169, many of them belonged to galactose metabolism (lacABCDEF) and phosphotransferase system (PTS) (treP, pfkB) ( Figure 6). Figure 6. KEGG analysis of the DEGs impacted by φSA169. Fatty acid biosynthesis had the most genes up-regulated, compared to the other pathways, while galactose metabolism and phosphotransferase system (PTS) were the pathways that had most genes down-regulated.

Verification of the Selected DEGs
DEGs that were up-/down-regulated in both comparisons 301-188::φSA169 vs. 301-188 and 300-169 vs. 301-188 were thought to be significantly impacted by φSA169. The expression of four DEGs (fabH, purA, lacF, and treP) involved in different KEGG pathways was selected to verify the RNA-seq results using qRT-PCR. Primers for the selected genes are listed in Table  S11. Genes fabH/purA and lacF/treP represented significantly up-and down-regulated DEGs by φSA169, respectively. The fold changes of the four genes determined by the qRT-PCR were similar to the results obtained in the RNA-seq assays (Figure 7).

Verification of the Selected DEGs
DEGs that were up-/down-regulated in both comparisons 301-188::φSA169 vs. 301-188 and 300-169 vs. 301-188 were thought to be significantly impacted by φSA169. The expression of four DEGs (fabH, purA, lacF, and treP) involved in different KEGG pathways was selected to verify the RNA-seq results using qRT-PCR. Primers for the selected genes are listed in Table S11. Genes fabH/purA and lacF/treP represented significantly up-and down-regulated DEGs by φSA169, respectively. The fold changes of the four genes determined by the qRT-PCR were similar to the results obtained in the RNA-seq assays (Figure 7).

φSA169-Specific KEGG Analyses
The overlapping DEGs of 301-188::φSA169 vs. 301-188 and 300-169 vs. 301-188 might represent the genes regulated explicitly by φSA169 (Figure 4). The KEGG profile of the overlapping DEGs indicated that most of these genes were involved in metabolic pathways ( Figure 6). For instance, the DEGs of fatty acid biosynthesis (fabFH), purine metabolism (purA), and RNA degradation (AS94_08925) were up-regulated by φSA169. Among the down-regulated DEGs by φSA169, many of them belonged to galactose metabolism (lacABCDEF) and phosphotransferase system (PTS) (treP, pfkB) ( Figure 6). Figure 6. KEGG analysis of the DEGs impacted by φSA169. Fatty acid biosynthesis had the most genes up-regulated, compared to the other pathways, while galactose metabolism and phosphotransferase system (PTS) were the pathways that had most genes down-regulated.

Verification of the Selected DEGs
DEGs that were up-/down-regulated in both comparisons 301-188::φSA169 vs. 301-188 and 300-169 vs. 301-188 were thought to be significantly impacted by φSA169. The expression of four DEGs (fabH, purA, lacF, and treP) involved in different KEGG pathways was selected to verify the RNA-seq results using qRT-PCR. Primers for the selected genes are listed in Table  S11. Genes fabH/purA and lacF/treP represented significantly up-and down-regulated DEGs by φSA169, respectively. The fold changes of the four genes determined by the qRT-PCR were similar to the results obtained in the RNA-seq assays (Figure 7).

Discussion
Many phages carry virulence factors that significantly contribute to genome variation, pathogenesis, and antibiotic resistance in S. aureus [7,19,20]. Despite the obvious importance of phages, studies on the interactions between phage and MRSA persistent outcome are limited. Recently, we demonstrated that the lysogenization of clinical RB 301-188 strain with phage φSA169 resulted in persistent phenotypes in vitro and in an experimental endocarditis model [11]. Thus, the current study was designed to determine the impact of φSA169 on genetic factors that may contribute to persistent MRSA endovascular infections.
The RNA-seq results revealed that MRSA host genes up-regulated by φSA169 were mainly involved in fatty acid biosynthesis (fabF and fabH), purine (purA), pyrimidine (AS94_12220), and RNA degradation (AS94_08925). Both fabF and fabH encode essential enzymes for fatty acid biosynthesis in many pathogens, including S. aureus [21]. Fatty acids are crucial hydrophobic components of membrane lipids and are important metabolic energy sources in bacteria [22]. It has been reported that defected unsaturated fatty acid biosynthesis in Streptococcus mutans results in attenuated virulence (e.g., less transmissible, less carious lesions) in a rodent model of dental caries [23]. In addition, fatty acid biosynthesis contributes to virulence in Group B Streptococcus (GBS) [24]. Importantly, fatty acid biosynthesis pathway inhibition has been investigated as a possible antimicrobial agent in bacteria [25]. In the current study, significantly higher expressions of fabF and fabH were observed in the φSA169-carrying strains, which may result in survival advantage and consequent persistence.
As a member of pur regulon, purA encodes the enzyme that catalyzes the conversion of inosine-5-phosphate (IMP) to adenylosuccinate [26]. We and others have previously shown that purine biosynthesis promotes virulence and persistence in S. aureus [14,[26][27][28]. For instance, the inactivation of purA causes the lower expression of a broad spectrum of genes (e.g., energy production and conversion) and attenuates the ability of S. aureus to cause kidney infection in mice [27]. Li et al. reported that higher purine biosynthesis production correlates with persistent outcomes in an experimental MRSA endocarditis model [14].
In addition, several studies demonstrated that the inactivation of purine biosynthesis repressor, purR, leads to a greater amount of secreted virulence factors and hypervirulence in the murine model of S. aureus bacteremia model [26,28]. In the current study, the purine biosynthesis gene, purA, was found to be significantly up-regulated by φSA169. Therefore, φSA169-related higher purA expression might contribute to the persistent outcomes we observed in our recent study [11].
It is also interesting that φSA169 significantly down-regulated several genes related to the galactose metabolism. Galactose is a common monosaccharide used by organisms [29]. S. aureus employs lac operon to import and metabolize galactose [30]. In a previous study, the down-regulation of lac operon was observed in a rpoB (A621E) mutant S. aureus strain that had decreased susceptibility to vancomycin compared to the parental strain [31]. Therefore, down-regulated lac operon in the φSA169-carrying strains might contribute to the persistent outcomes with VAN treatment in vivo [11]. However, more research into galactose metabolism and its role in pathogenesis and persistence in S. aureus is needed. The RNA-seq displayed down-regulation of the phosphotransferase system (PTS) by φSA169. It has been demonstrated that the PTS plays an important role in carbohydrate transport, and the regulation of sugar utilization genes, which further contributes to overall metabolic efficiency in Gram-positive bacteria [32,33]. Gera et al. reported that deleting ptsI that encodes cytosolic enzyme I (EI) (∆ptsI) in group A Streptococcus (GAS) strains resulted in a hypervirulent phenotype compared to their respective wild-type strains (e.g., significantly increased skin lesion severity and size) in a murine model of disseminating skin and soft tissue infection [33]. Thus, PTS appears to reduce the virulence of GAS skin infection. However, a conflict phenotype of interrupted ptsI in S. aureus was reported with an attenuated virulence compared to its wild-type strain in a systemic infection model [34]. We suspect this discrepancy is possibly due to (i) the PTS regulation of virulence factors in GAS vs. S. aureus and (ii) the animal models used (skin and soft tissue infection vs. systemic infection). Importantly, galactose is one of the carbohydrates that utilizes PTS [35]. Thus, down-regulated PTS in φSA169-carrying strains might be correlated with the lower expression of galactose metabolism genes. Detailed studies are needed to define the specific role of PTS, and the interaction between PTS and galactose, in persistent MRSA endovascular infection.
In this study, we also observed that some genes within the mutual prophage in both 300-169 and 301-188 strains were negatively impacted by φSA169, which suggested that the mutual prophage genes might be another φSA169-derived genetic factor that participated in the PB outcomes. It has been reported that the pathogenesis of S. aureus Newman requires the participation of its all four prophages [7]. Thus, multiple prophages might have combined effects on virulence and pathogenesis in S. aureus. Therefore, φSA169 might contribute to the PB outcomes by mediating the gene expression of the mutual prophage.
Besides the impact of φSA169 on genetic factors in the MRSA host genes and the mutual prophage, the different genetic backgrounds between 300-169 and 301-188 strains might also play a role in the persistent outcomes ( Figure S1). We have previously demonstrated that key global regulators were differently expressed in 300-169 and 301-188 [11,14]. These differences may impact downstream virulence factors, subsequently contributing to the persistent outcome.
We recognize that there were some significant limitations in the current study. For instance, we only studied a PB 300-169 WT (300-169) containing φSA169, RB 301-188 WT (301-188) in the absence of φSA169, and 301-188 WT with φSA169 lysogenization (301-188::φSA169) in the current and previous research [11]. It would be important to verify the genetic impact of φSA169 using φSA169 deletion in the PB 300-169 strain background. In addition, it would be interesting to study the combinational effect of VAN with φSA169 on the MRSA host and φSA169 genes, which may demonstrate how φSA169 mediates the response to VAN treatment in the IE model [11]. Therefore, further investigations are needed to address these limitations.

Conclusions
In this study, we explored the impact of prophage φSA169 on genetic factors, which might play an essential role in MRSA-persistent endovascular infection. The results highlighted that φSA169 contributed to PB outcomes mainly through mediating metabolisms, especially the up-regulation of fatty acid biosynthesis and down-regulation of galactose metabolism and PTS. In addition, the mutual prophage in both 300-169 and 301-188 strains and different genetic backgrounds of these two strains might also be the genetic factors that contribute to the PB outcomes.