Genome Mining and Comparative Pathogenomic Analysis of An Endemic Methicillin-Resistant Staphylococcus Aureus (MRSA) Clone, ST612-CC8-t1257-SCCmec_IVd(2B), Isolated in South Africa.

This study undertook genome mining and comparative genomics to gain genetic insights into the dominance of the methicillin-resistant Staphylococcus aureus (MRSA) endemic clone ST612-CC8-t1257-SCCmec_IVd(2B), obtained from the poultry food chain in South Africa. Functional annotation of the genome revealed a vast array of similar central metabolic, cellular and biochemical networks within the endemic clone crucial for its survival in the microbial community. In-silico analysis of the clone revealed the possession of uniform defense systems, restriction-modification system (type I and IV), accessory gene regulator (type I), arginine catabolic mobile element (type II), and type 1 clustered, regularly interspaced, short palindromic repeat (CRISPR)Cas array (N = 7 ± 1), which offer protection against exogenous attacks. The estimated pathogenic potential predicted a higher probability (average Pscore ≈ 0.927) of the clone being pathogenic to its host. The clone carried a battery of putative virulence determinants whose expression are critical for establishing infection. However, there was a slight difference in their possession of adherence factors (biofilm operon system) and toxins (hemolysins and enterotoxins). Further analysis revealed a conserved environmental tolerance and persistence mechanisms related to stress (oxidative and osmotic), heat shock, sporulation, bacteriocins, and detoxification, which enable it to withstand lethal threats and contribute to its success in diverse ecological niches. Phylogenomic analysis with close sister lineages revealed that the clone was closely related to the MRSA isolate SHV713 from Australia. The results of this bioinformatic analysis provide valuable insights into the biology of this endemic clone.

1. Introduction

The increased global detection of methicillin-resistant Staphylococcus aureus (MRSA) has led to its categorization as a high priorty pathogen in the World Health Organization’s (WHO) priority pathogen list for the research and development of new antibiotcs [1]. It is thus imperative to identify, profile and understand its diversity in various settings to provide insights into its evolution and spread [2]. MRSA has been categorised into three main types, namely healthcare-acquired (HA), community-associated (CA) and livestock-associated (LA) MRSA based on its involvement in nosocomial-, community- and livestock-associated infections, respectively. Although MRSA was assumed to have restricted human or animal hosts, recent reports indicate a blurring epidemiology of emerging MRSA clones between different hosts attributable to progressive clonal expansion, adaptability, and transmission between the various host types, thus making the traditional definitions indistinct [3,4,5].

MRSA is a leading foodborne bacterium implicated in the food chain with associated risk to occupationally exposed workers and food animals, mainly poultry, pigs and cattle [6,7,8,9,10]. A diversity of MRSA clones, including ST398, have been associated with infection and colonization in poultry and exposed personnel worldwide [11,12,13,14] with characteristics that indicate endemicity to a particular region and production system [14,15,16,17]. These survival features that favour their successful spread include the ability to express their pathogenicity through virulence determinants allowing them to adhere, lyse, invade, harm and escape from their host at many levels [18,19,20]. Furthermore, their ability to detect and react to various changes in environmental conditions related to stress (osmotic, oxidative and periplasmic), heat shock and toxification are also critical for their dominance and clonal expansion in diverse niches. Thus, a combination of virulence and environmentally acquired factors make S. aureus a highly efficient pathogen [21].

The ST612-MRSA clone has only been reported in Australia and South Africa in horses, companion animals, equine veterinarians, community-onset and nosocomial infections highlighting the blurring epidemiology of this endemic clone [15,22,23,24,25]. This parent clone has undergone clonal expansion into two main types, ST612-MRSA-IV-A and ST612-MRSA-IV-B, and it has been postulated that its global spread has occurred via travel to these popular tourist countries [14,26]. In South Africa, the clonal subtype ST612-CC8-t1257-SCCmec_IVd(2B) is dominant and has spread to different provinces across the country [16,27,28,29,30]. This multi-drug-resistant (MDR) local dominant clone was recently found in poultry and occupational farm workers with a diverse resistome to antibiotics commonly used in human and veterinary medicine with the concern of its spread across the poultry food-chain in the country [31]. Notwithstanding several reports on this endemic clone over the past decade, information on its pathogenic and specific features associated with dominance and persistence are lacking. Therefore this study undertook comparative pathogenomics analysis of the MRSA clone ST612-CC8-t1257-SCCmec_IVd(2B) isolated from the poultry farm system in South Africa using whole genome sequencing (WGS) to characterize its pathogenicity, virulome, environmental tolerance and persistence mechanisms.

2. Results

2.1. Demographics and Whole Genome Sequencing (WGS) In-Silico Analysis of Methicillin-Resistant Staphylococcus Aureus (MRSA)

A total of 145 S. aureus isolates were recovered across the ‘’farm to fork’’ continuum (farm; transport; slaughterhouse and retail points). Twelve (8.27%) out of the 145 S. aureus isolated were resistant to cefoxitin and were positive for mecA by PCR. WGS in-silico analysis revealed that 11 of the 12 MRSA isolates belonged to the same clone with sequence type (ST612), clonal complex (CC8), spa type (t1257) and harboured the SCCmec type (SCCmec_IVd (2B)). This endemic clonal subtype ST612-CC8-t1257-SCCmec_IVd(2B) was further investigated using bioinformatic analysis to gain genetic insights into their dominance in the food-chain. The critical points within the poultry continuum where the endemic clone was obtained included: animals on the farm (faecal samples; n = 4), occupational farm workers (nasal swabs, n = 3), slaughterhouse (carcass rinsate, n = 2) and retail outlets (whole carcass n = 2) (Table 1). The draft genomes of these isolates (deposited under Bioproject number PRJNA506780) were used for the bioinformatics analysis. The varied sources and sampling points across the food-chain allowed an understanding of the successful ecological dominance of the endemic clone in colonizing diverse niches via comparative analysis (Table 1).

Table 1

Demographics, genome features, classification, pathogenicity score, number and strain of pathogenic family linkage of the isolates belonging to the endemic clone and sister lineages.

Clone (ST612-CC8-t1257-SCCmec_IVd(2B) a					Genome Characteristics		Classification			e Pathogenicity Score(Number of Pathogenic Families)
Number	Strain ID	Point	Host	Source	Tandem Repeats	CRISPRs (Cas Cluster)	RM-System b	Agr type c	ACME d
1	SA1	Farm	Animal	Faecal	135	8 (1)	Type I and IV	Type I	II	0.926 (939) f
2	SA2	Farm	Animal	Faecal	127	8 (0)	Type I and IV	Type I	II	0.930 (895) f
3	SA3	Farm	Animal	Faecal	128	8 (1)	Type I and IV	Type I	II	0.925 (910) f
4	SA4	Farm	Animal	Faecal	133	9 (0)	Type I and IV	Type I	II	0.926 (878) f
5	SA5	Farm	Human	Nasal	136	6 (1)	Type I and IV	Type I	II	0.925 (931) f
6	SA6	Farm	Human	Nasal	129	6 (0)	Type I and IV	Type I	II	0.926 (907) f
7	SA7	Farm	Human	Nasal	145	7 (1)	Type I and IV	Type I	II	0.925 (912) f
8	SA8	Abattoir	Animal	Rinsate	143	8 (0)	Type I and IV	Type I	II	0.928 (870) g
9	SA12	Abattoir	Animal	Rinsate	129	7 (1)	Type I and IV	Type I	II	0.926 (903) f
10	SS1	Retail point	Animal	Carcass	133	8 (2)	Type I and IV	Type I	II	0.931 (841) f
11	SS2	Retail point	Animal	Carcass	258	7 (0)	Type I and IV	Type I	II	0.929 (828) f
ST8-CC8-t008-SCCmec_IVa (2B)
12	USA300	- h	Human	-	134	8 (0)	Type I and IV	Type I	I	0.924 (1094) f
ST8-CC8-t064-SCCmec_IVd (2B)
13	USA500	-	Human	-	133	7 (0)	Type I and IV	Type I	III	0.921 (1021) f
ST612-CC8-t064-SCCmec_IVd (2B)
14	SHV713	-	Animal	-	122	7 (0)	Type I and IV	Type I	II	0.924 (956) f

a Isolates belonged to the same clone with sequence type (ST612), Clonal complexe (CC8), spa type (t1257) and SCCmec type (SCCmec_IVd(2B)). b RM-System: Restriction-Modification System. c Agr: accessory gene regulator. d ACME-arginine catabolic mobile element: ACME type I (arcA+/opp3AB+), II (arcA+/opp3AB−) and III (arcA−/opp3AB+). e Pathogenicity score: Prediction of a bacteria’s pathogenicity towards the hosts using PathogenFinder. Strain of the closet pathogenic family linkage: f Staphylococcus aureus subsp. aureus USA300 (Accession number: CP000255) and g Staphylococcus aureus subsp. aureus JH9 (Accession number: CP000703). h Not applicable.

2.2. Phenotypic and Genotypic Characteristics of the 11 MRSA Belonging to the ST612-CC8-t1257-SCCmec_IVd(2B) and Sister Lineages

The antibiotic susceptibility testing (AST) results of the MRSA isolates are summarized in Supplementary Table S1. Putative acquired resistance genes and chromosomal mutations causing resistance to the different antibiotics were found (Table 2). Overall, the isolates belonging to the ST612-CC8-t1257-SCCmec_IVd(2B) clone exhibited a similar resistome profile which resonated with the AST (Table S1 and Table 2). However, resistome comparison with sister lineages showed marked differences. Interestingly, the endemic clone harboured more resistance (acquired and chromosomal) than all the other clones. The USA300 only harboured the mecA gene (Table 2).

Table 2

Genotypic characteristics of the MSRA isolates.

Isolate	Resistance Mechanisms
SA1	mecA, blaZ, aac(6’)-aph(2’’), erm(C), tet(M), gyrA(S84L), parC(S80Y), parE(D434N), rpoB (H481N)
SA2	mecA, blaZ, aac(6’)-aph(2’’), erm(C), mrs(A) tet(M), dfrC, gyrA(S84L), parC(S80Y), parE(D434N), rpoB (H481N)
SA3	mecA, blaZ, aac(6’)-aph(2’’), msr(A), tet(M), dfrC, gyrA(S84L), parC(S80Y), parE(D434N), rpoB (H481N)
SA4	mecA, blaZ, aac(6’)-aph(2’’), tet(M), dfrC, gyrA(S84L), parC(S80Y), parE(D434N), rpoB (H481N)
SA5	mecA, blaZ, msr(A), tet(M), dfrC, gyrA(S84L), parC(S80Y), parE(D434N), rpoB(H481N)
SA6	mecA, blaZ, msr(A), tet(M), dfrC, gyrA(S84L), parC(S80Y), parE(D434N), rpoB (H481N)
SA7	mecA, blaZ, aac(6’)-aph(2’’), tet(M), dfrC, rpoB (H481N)
SA8	mecA, blaZ, aac(6’)-aph(2’’), tet(M), dfrC, gyrA(S84L), parC(S80Y), parE(D434N), rpoB (H481N)
SA12	mecA, blaZ, aac(6’)-aph(2’’), erm(C), tet(M), dfrC, gyrA(S84L), parC(S80Y), parE(D434N), rpoB (H481N)
SS1	mecA, blaZ, aac(6’)-aph(2’’), erm(C), mph(C), msr(A), tet(M), tet(K), dfrC, gyrA(S84L), parC(S80Y), parE(D434N), rpoB (H481N)
SS2	mecA, blaZ, aac(6’)-aph(2’’), tetM, dfrC, gyrA(S84L), parC(S80Y), parE(D434N), rpoB (H481N)
USA300	mecA
USA500	mecA, blaZ, tetM, gyrA(S84L), parC(S80F)
SHV713	mecA, blaZ, tetM, rpoB (H481N)

Staphylococcus aureus subsp. aureus ATCC 25923 (Accession no. CP009361) was used as the wild type strain to elucidate fluoroquinolone and rifampicin resistance caused by chromosomal mutation.

2.3. General Features of Endemic Clonal Subtype ST612-CC8-t1257-SCCmec_IVd(2B)

The circular genome (CG) Viewer Server revealed a circular map for the isolates (Figure 1).

Figure 1

Graphic depiction of the circular map of the SA3 genome belonging to ST612-CC8-t1257-SCCmec_IVd(2B) clone. The two-outer circles show the open reading frame (ORF–in light blue). The inner circle shows the GC skew, with green and purple indicating positive and negative values, respectively. The GC content is indicated in black. This genome map was visualized using the CGView Server (http://stothard.afns.ualberta.ca/cgview_server/index.html).

Genome statistics of the clone are shown in Supplementary Table S2. The average size, GC content, number of contigs, N50, and L50, were 2.98 Mbp, 33.14%, 668, 11,425 bp, and 77 bp, respectively. Annotation with RAST (rapid annotation using subsystem technology) and PGAP (Prokaryotic Genome Annotation Pipeline) resulted in an average of 3235 protein-coding genes (CDSs) and 58 RNAs (Supplementary Table S2). The CDSs were assigned to 26 clusters of orthologous groups (COGs) and functional categories (Supplementary Table S3). The SEED subsystem distribution predicted amino acids and derivatives (300 genes), carbohydrates (203 genes), protein metabolism (181 genes), cofactors, vitamins, prosthetic groups, pigments (131 genes), membrane transport (60 genes) and DNA metabolism (79 genes) to be the most abundant functional categories in the clone (Supplementary Table S3).

2.4. In-Silico Prediction of Clustered, Regularly Interspaced, Short Palindromic Repeats (CRISPRs), Restriction-Modification System (R-M system), Accessory Gene Regulator (Agr) Type, and Arginine Catabolic Mobile Element (ACME)

The CRISPRCasFinder predicted clustered, regularly interspaced, short palindromic repeat (CRISPR) arrays (average = 7 ± 1) encoding the type I CRISPR-Cas systems on the nodes/contigs of the strains. Of note, the in-silico analysis revealed that all the isolates in the clone carried the same Restriction-Modification (R-M) system (type I and IV), Agr (type 1) and ACME (type II) (Table 1). Tandem repeats finder predicted an average of 145 repeats in the endemic clone (Table 1). Comparison with close lineages showed similar CRISPRCas system, R-M system (type I and IV) and Agr (type 1) except for their ACME types that differed (Table 1).

2.5. Pathogenicity and Virulome Insights of the Endemic Clone

The endemic clonal subtype ST612-CC8-t1257-SCCmec_IVd(2B) had a 0.927 mean probability (P score) of being pathogenic to hosts and was found to match with 892 pathogenic families. All the pathogenic families were linked to the family Staphylococcaceae of which the Staphylococcus aureus subsp. aureus USA300 strain (Identity: 100%) was the organism with the highest linkage. Furthermore, the P score of the ST612, its sister lineages (ST8), as well as other well-known epidemic MRSA clones (ST5, ST59, ST239) predicted were similar ranging from 0.921–0.933 without any significant differences (Table 1 and Table 3). The whole virulome analysis predicted a total of 78 putative virulence-encoded proteins belonging to eight major virulence determinants of Staphylococcus, namely: adherence factors, antiphagoctyosis, exo-enzymes, immune evasion, iron uptake, plasminogen activator, secretion system and toxins (Figure 2 and Supplementary Tables S4–S6). Of note, 68 conserved virulence factors were observed across the endemic clone. The minor differences in the virulence factors of the clone occurred in the possession of a biofilm operon system (icaABCDR), an extracellular adherence protein/major histocompatibility complex (MHC) analogous protein (eap/map), the acquisition of hemolysin (hld) and enterotoxins (sec, seb, selk, selp) (Figure 2 and Supplementary Tables S4–S6). Interestingly, the epidemic reference genome S. aureus subsp. aureus (USA300_TCH1516, 2,872,915 bp, NC010079) used for the comparison contained all the predicted putative virulence factors illustrated in Figure 2.

Table 3

Comparative in-silico prediction of the of the pathogenicity of known epidemic MRSA clones.

Isolate	Sequence Type(Clonal Complexe)	Pathogenicity Score(Number of Pathogenic Families)	Accession Number
ISU935	ST5 (CC5)	0.933 (899)	CP017090.1
JKD6008	ST239 (CC8)	0.927 (1036)	CP002120.1
M013	ST59 (CC59)	0.930 (337)	CP003166.2

Figure 2

Heatmap generated with phylogeny and distribution of virulence factors across the endemic clone; ST612-CC8-t1257-SCCmec_IVd(2B). The green colour represents the presence of the gene, and the yellow colour represents the absence of the gene. The virulence factors are represented by Roman numerals I: adherence factors; II: immune evasion; III enzymes; IV secretion systems; V anti-phagocytosis and VI: Toxins.

2.6. Genomic Encoding Mechanisms and Pathways of Bacterial Persistence and Tolerance

All the isolates across the food-chain contained conversed factors involved in stress tolerance (osmotic, oxidative and periplasmic) regulation (Table 4). The endemic strains possessed heat shock proteins (GroES, GroEL, S4 paralog and GrpE), bacteriocins (BceR, BceA, BceB and BceS) and putative systems for detoxification (Table 4). This was also observed in all the sister strains (USA300, USA500 and SHV713).

Table 4

In-silico identification and characterization of tolerance and persistence mechanisms in the endemic clone and sister clones (USA300, USA500 and SHV713).

Type	Associated Proteins/Enzymes/Genes
Stress response
Protection from Reactive Oxygen Species	Superoxide dismutase (Fe) (EC 1.15.1.1)
Oxidative stress	Superoxide dismutase (Mn) (EC 1.15.1.1)
	Ferric uptake regulation protein FUR
	Peroxide stress regulator (PerR), FUR family
	Organic hydroperoxide resistance protein
CoA disulfide thiol-disulfide redox system	CoA-disulfide reductase (EC 1.8.1.14)
Glutathione: Redox cycle	Glutathione peroxidase (EC 1.11.1.9)
Osmotic stress
Osmoregulation	Glycerol uptake facilitator protein (GlpF)
Choline and Betaine Uptake and Betaine Biosynthesis	Choline ABC transport system, permease protein (OpuBB)
	Choline ABC transport system, permease protein (OpuAC)
	Choline ABC transport system, permease protein (OpuAA)
	Betaine aldehyde dehydrogenase (EC 1.2.1.8)
	Choline ABC transport system, permease protein (OpuBC)
	Glycine betaine transporter (OpuD)
	Choline ABC transport system, permease protein (OpuBD)
	Glycine betaine ABC transport system, permease protein (OpuAB)
	Choline ABC transport system, permease protein (OpuBA)
	Choline dehydrogenase (EC 1.1.99.1)
Periplasmic stress
Periplasmic Stress Response	Intramembrane protease (RasP/YluC)
Others
SigmaB stress response regulation	Serine-protein kinase (RsbW) (EC 2.7.11.1)
	Serine phosphatase (RsbU), regulator of sigma subunit
	RNA polymerase sigma factor (SigB)
	Anti-sigma B factor antagonist (RsbV)
Bacterial hemoglobins	Hemoglobin-like protein (HbO)
High frequency lysogen (Hfl) operon	RNA-binding protein (Hfq)
Heat shock	GroES, GroEL, S4 paralog and GrpE
Bacteriocins	Two-component response regulator (BceR)Bacitracin export ATP-binding protein (BceA)Bacitracin export permease protein (BceB)Two-component sensor histidine kinase (BceS)
Detoxification	Nudix proteins (nucleoside triphosphate hydrolases):8-oxo-dGTPase Bsu (YtkD) and ADP-ribose pyrophosphatase (EC 3.6.1.13)
	Housecleaning nucleoside triphosphate pyrophosphatases:Nucleotidase (YfbR), HD superfamily and Dimeric dUTPase (EC 3.6.1.23)
	Nucleoside triphosphate pyrophospho-hydrolase (MazG)

2.7. Comparative Phylogenomic Insights of the Endemic Clone

The phylogenomic analysis depicted a clear distinction between different clones, although they belonged to the same clonal complex, CC8 (Figure 3).

Figure 3

A graphical view of the phylogenetic tree of the 11 MRSA isolates belonging to the endemic clonal subtype and its closet lineages belonging to the same clonal complex (CC) 8. The patterns on the cladogram are: purple colour depicting the endemic clonal subtype ST612-CC8-t1257-SCCmec_IVd(2B), from the critical points in the poultry farm system in South Africa, yellow indicating the ST612-CC8-t064-SCCmec_IVd(2B)-SVH7513 strain (Accession number: CP029166) (the first complete genome of the ST612; a livestock-associated MRSA isolated in Australia), black indicating the ST8-CC8-t064-SCCmec_IVd(2B)-S. aureus strain 2395-USA500 (Accession number: CP007499.1) (an epidemic community-associated MRSA strain) and red depicting ST8-CC8-t008-SCCmec_IVa(2B)-S. aureus USA300 strain TCH1516 (Accession number: CP000730) (a healthcare-associated MRSA strain).

The isolates of the endemic clone were more closely related to the Australian strain isolated from a horse followed by the USA500 strain which shared ancestry with these clones while the USA300 isolate was the least related strain (Figure 3). This collaborated with their molecular typing schemes (sequence, spa and SCCmec types). The average nucleotide (orthoANIu) value confirmed the phylogenetic analysis, predicting an average ≥99.90% between endemic clones which differed with an average nucleotide identity (ANI) value of 99.85%, 99.81% and 99.70% for SHV713, USA500 and USA300 respectively, in decreasing order of similarity. Metadata coupled with the SNPs tree of the endemic clone offered a clear visualisation of the conversed and variable virulence determinants associated with the clonal subtype (Figure 3).

3. Discussion

The ability of MRSA clones to colonise and cause infection in a niche has been attributed to specific features that favour their clonal spread and survival. The dominant and endemic clonal subtype ST612-CC8-t1257-SCCmec_IVd(2B) has been reported for its ecological success in South Africa across various sectors including the food-chain. However, studies on its pathogenicity as well as specific characteristics associated with dominance and persistence are lacking. Therefore, this study undertook comparative pathogenomics analysis of the MRSA clone ST612-CC8-t1257-SCCmec_IVd(2B) in South Africa using WGS and comprehensive bioinformatic tools.

There was a general similarity with regards to the phenotypic and resistance patterns in the MRSA isolates of human and animal origin in the clone. This is because all the MRSA isolates in the clone harboured the same SCCmec mobile genetic element (SCCmec_IVd(2B)) carrying mecA mediating resistance to methicillin and other semisynthetic penicillinase-resistant β-lactams that are frequently co-carried with genes conferring resistance to aminoglycosides, macrolide-lincosamide-streptogramin B (MLSB) and spectinomycin of the MRSA isolates. This co-selection phenomenon renders most of the antibiotics useless as treatment options [32].

Functional screening of the genomes revealed a wide array of similar central metabolic and biochemical networks (carbohydrate, protein, amino acid breakdown networks), regulatory models (regulation and cell signalling systems), cellular processes (cell wall and capsule formation, cell division and DNA uptake) within the endemic clone which are crucial for the survival of the clones in the microbial community (Supplementary Table S3) [33,34,35].

In-silico analysis of the clone predicted the presence of the CRISPR array (Table 1) responsible for adaptive immunity mechanisms in prokaryotes that protect bacterial clones against invading viral attacks [36,37,38]. The CRISPR immune system has been postulated as a critical force for the sustainability of bacteria in the microbial landscape, and this may have contributed to the dominance of this endemic clone in diverse settings [39]. Further analysis revealed other uniform defence systems in the endemic clone such as the type I and IV restriction-modification (R-M) system which provides protection against exogenous DNA and thus increases its survival [40]. The two R-M systems found in the strains were not peculiar as about 12.5% of bacteria normally harbour more than one R-M system [41]. More so, the typing of the accessory gene regulator (agr) system modulating the expression of temporal genes encoding virulence determinants in S. aureus predicted the agr-type I system in all the isolates [42]. This corroborates other reports on agr-typing in ST612 isolates identified using polymerase chain reaction (PCR) and highlights the ability of WGS to accurately predict different typing techniques [43,44]. Comparative analysis with close lineages showed similar CRISPRCas systems, R-M systems (type I and IV) and Agr (type 1) except for their ACME types that differed (Table 1). This indicates that the MRSA USA300 strain, though it contains less resistance genes, harboured the more hypervirulent ACME type I which enhances the colonization, transmissibility, and persistence of the on host contributing to the global success of this lineage compared to the ST612-CC8-t1257-SCCmec_IVd(2B) [45].

Interestingly, the agr type 1 system has also been reported in the ST8 (USA300 and USA500) strains, all of which belong to the clonal complex (CC)-8 and are the closest lineages of the endemic clonal subtype in South Africa [43,46]. This quorum sensing agr-system enables isolates to monitor their local population density through the detection and secretion of small auto-inducer substances to regulate gene expression of temporal factors aiding their survival in the different niches [47]. Furthermore, the possession of a type II_ACME system encoding putative virulence factors with low pathogenicity enhances their ability to grow and survive on their host [48]. Of note, this is the first report of an ACME system in ST612 although its presence has been reported in its sister lineage, ST8 [23,48,49]. Interestingly, all these defence systems (CRISPR, R-M, Agr and ACME) (Table 1) which are also employed in typing the strain can serve as targets for drug development and need to be explored further [47].

The pathogenic potential (Pscore) is a probabilistic value indicating the possibility of the strain been pathogenic to the host with the probability ranging from 0 to 1. Therefore, the P score is a machine learning algorithm used for the in-silico prediction of the pathogenic potential an isolate. Estimation of the pathogenic potential of the clone using these trained algorithms to distinguish between friend or foe strains, predicted a higher probability (average Pscore ≈ 0.927) of the clone being pathogenic to its host (animal and human) (Table 1). The fact that this clone has been associated with community and hospital infections, as well as colonising humans and animals in South Africa and Australia, require further studies to substantiate this claim [15,22,23,24,25]. More so the Pscore of the ST612, its sister lineages (ST8) as well as other well-known epidemic MRSA clones (ST5, ST59, ST239) were similar, ranging from 0.921–0.933 without any significant differences indicating their possibility of been pathogenic to their host (Table 1 and Table 3). However, given the complex interplay between host-pathogen interactions, further experiments are required to ascertain the practicality of this distinctive measure (pathogenic potential), and it should hence be interpreted with caution [50,51].

Comparative analysis of the clone revealed the possession of a battery of virulence factors which were mostly conserved across the isolates in the endemic clone. The clone contained adherence and immune evasion factors, predominantly belonging to the microbial surface components recognizing adhesive matrix molecules (MSCRAMMs) family (atl, ebh, clfA, clfB, ebp, eap/map, efb, fnbA, sdrC, sdrD, sdrE and spa) (Figure 3 and Supplementary Table S4), that facilitate attachment and enable them to invade their host [52,53]. MSCRAMMs act as adhesin receptors to facilitate the attachment and specific binding to the extracellular matrix. The isolates also harboured intercellular adhesion proteins (icaA, icaB, icaC, icaD, icaR) associated with attachment, proliferation and differentiation of micro-colonies into unique biofilm structures protecting them from the host defences such as antibodies or phagocytosis, making them extremely difficult to eradicate [54,55]. The association of dominant clonal lineages and biofilm formation has been reported in various studies on MRSA [56,57,58,59]. Naicker et al. reported a strong biofilm formation in ST612 compared to the other clones which may have contributed to its dominance in South Africa [60]. The clone contained only one capsular serotype, cap5, which enhances microbial virulence and offers protection against phagocytic uptake (Figure 2 and Supplementary Table S5) [61]. They also harboured the highly versatile type VII secretion system encoding the eight cluster genes associated with the release of secretions that cause pathogenesis in S. aureus [62,63].

The ST612-CC8-t1257-SCCmec_IVd(2B) contained an array of six putative toxins which may play a significant role in its pathogenesis and survival (Figure 2 and Supplementary Table S6). The possession of hemolysins (Hly/hla and hld) offers the strain the ability to induce cell membrane damage, triggering cytokine formation, and reducing or killing neutrophils [64]. The clone harboured the eta gene, a member of exfoliative toxins that cause bullous impetigo and staphylococcal scalded skin syndrome when expressed [65,66]. Moreover, they contained a set of exotoxins (set30, set31, set34, set35, set36, set37, set38, set39, set40) which are involved in host tissue damage and aid in diverting the immune response to the bacteria [67]. More so, they encoded enterotoxins (seb, sec, selk, selp, selq) that have been associated with symptoms of food poisoning raising an associated food safety threat [67,68]. Of note, the endemic clone possessed only the panton-valentine leucocidin (PVL) F component (lukF-PV) which forms a bi-component layer with gamma hemolysins (hlgA, hlgB, hlgC) to induce cell activation leading to a Ca2+ influx and apoptosis of the cell [69,70,71]. There have been two main reports on the association of PVL in ST612-MRSA-IV isolates indicating their presence and absence in the clone [44,72]. Thus, the linkage between SCCmec type IV and PVL production is still debatable. The differential expression of these putative virulence determinants probably confers a competitive advantage, contributing to its remarkable success as a pathogen. The detection of these virulence genes can delineate the most prevalent exposing proteins which can culminate in the development of novel vaccines for this local clone [73].

Bacteria deal with continuous stress and thus develop complex mechanisms to cope with various kind of environmental pressures and stressors to survive the most hostile environments [74]. We identified several mechanisms contributing to the survival of this clone in extreme environments (Table 4). The clone exhibited genomic signatures to regulate osmolarity and pH level in order to establish themselves in the diverse biological niches [75,76]. For example, to overcome dramatic changes in osmolarity (the concentration of a solution) or osmotic stress, they harboured the glycerol uptake protein facilitator (glpF) reported to play a pivotal role in glycerol efflux and influx during hyperosmotic conditions [77]. They also possessed the osmoprotectant, Choline-glycine betaine uptake and biosynthesis system, which operates the feedback mechanism via the up- and down- regulation of both molecules (choline and glycine betaine) to control high salt stress in the micro-environment [78,79]. Additionally, both choline and glycine betaine can serve as sources of carbon and nitrogen for the strains to grow and survive in a variety of niches [78,79]. For oxidative stress response, the clone contained super-oxidase dismutase, reductase and peroxidase enzymes that break harmful products by neutralizing them before they cause damage to essential cellular components, including DNA, membrane lipids and proteins [80]. The periplasmic sensing and signal transmission mechanisms, and the intramembrane metallo-protease (RasP/YluC) implicated in the cleavage of anti-sigma factor regulons to salvage protein misfolding caused in cell stress conditions were also present in the clone [81,82]. The clone was armed for heat shock response by inducing the expression of molecular chaperones (heat shock proteins) to combat the adverse effects on proteins caused by stressors such as heavy metals, oxidative stress heavy and increased temperatures in the microenvironment [83,84,85]. Comparative analysis of the other ancestry lineages (USA300, USA500 and SHV713) of the ST612-CC8-t1257-SCCmec_IVd(2B) clonal subtype revealed similar stress response parameters in all the genomes indicating the possibility of a conserved environmental tolerance and persistence mechanisms. However, stress biology is still a matter of active research. Hence, further studies will be needed to understand the actual molecular mechanisms regulating the tolerance and persistence phenotypes in the clone to aid in the identification of new targets for developing innovative anti-infective treatments.

Phylogenomic analysis depicted a clear distinction between the clone and sister lineages (USA300 and USA500) belonging to the same clonal complex (CC8) but shared the closet ancestry with the Australian strain isolated from a horse (SHV7513) (Figure 3). This is because ST612 is a double-locus variant of ST8 explaining why the USA300 and USA500 are distinct from this endemic clone [2,30,44]. Furthermore, the endemic clone was found on a different phylogenetic branch from the SHV7513, indicating its uniqueness and the clonal expansion of the ST612 parent clone (Figure 2). This corroborated the typing of the isolates and ANI value highlighting the ability of WGS to tentatively interpret genetic differences among strains [86,87]. Metadata coupled with the phylogenomic tree of the endemic clone offered a clear visualisation of the conserved and variable virulence determinants associated with the clonal subtype (Figure 3). Further studies are required to ascertain the implication of the minor differences in the virulence determinants (biofilm operon system, extracellular adherence protein/MHC analogous protein and acquisition of toxins) of the clone (Figure 2 and Supplementary Tables S4–S6). Further studies exploring the difference in the genetic compositions of ST612 and other major epidemic MRSA clones is recommended to understand the dynamics of this pathogen.

The results of this bioinformatic analyses would provide valuable and deeper insights into the ecological success and biology of this endemic clone. This would be a good step towards the development of strategies by the appropriate agencies to help curb this opportunistic pathogen which is on the rise

4. Materials and Methods

4.1. Study Design and Identification of MRSA

A longitudinal study was conducted in a poultry farm system in the uMgungundlovu District in KwaZulu-Natal (KZN), South Africa over a six-week period between 8 August and 14 September 2017. The poultry farm system (Farm A), its environments, with traced slaughterhouse and retail products served as the principal study sites for data and sample collection. A sampling scheme was set-up based on the World Health Organization Advisory Group on Integrated Surveillance of Antimicrobial Resistance (WHO-AGISAR) guidelines [88] to sample isolates across the farm continuum (animals on the farm, transport/holding, post-slaughter and retail products). On the fifth week, when the flock were ready to be slaughtered, swab samples from holding/transport (crates and trucks swabs) were collected during transportation of the target flock to the slaughterhouse of Farm A. At the slaughterhouse, upon the sacrifice of the flock, the carcass rinsate and caecal contents were collected. A portion of the poultry meat from the same flock, in the form of whole chicken, thigh and neck which are supplied to consumers in frozen packets, were purchased from Farm A at the retail point. All collected samples were immediately stored at 4 °C to maintain moisture, and cell viability. These samples were transported to the laboratory and processed within 4 hours of sampling.

4.1.1. Isolation, Identification and Molecular confirmation of S. aureus

All the samples were inoculated into tryptone soya broth (TSB) (Basingstoke, Hampshire, England) and incubated at 37 °C for 2 h while shaking (100 rpm). These samples were then streaked on HiCrome Aureus Agar Base (Himedia Laboratories, Mumbai, India) and incubated overnight at 37 °C in aerobic atmosphere. After incubation, colonies showing a unique brown black colour with a clear zone were streaked on mannitol salt agar (Himedia Laboratories, Mumbai, India) for further screening. Presumptive S. aureus colonies were examined for coagulase activity by the tube plasma agglutination test as well as DNAse tests [89]. The identified colonies were then confirmed using the API Staph kit (BioMérieux, Marcy-l’Etoile, France) and PCR using S. aureus species-specific primers for the nucA gene, which codes for thermostable nuclease [90].

4.1.2. Detection of MRSA and Antibiotic resistance testing

The phenotypic detection of MRSA isolates was performed as previously reported [91]. All isolates resistant to cefoxitin 30 μg (Oxoid, England) inhibition zone ≤ 21 mm were considered as MRSA [92]. MRSA was confirmed using real-time polymerase chain reaction (PCR) targeting the mecA gene [32,93]. S. aureus ATCC 25923 (susceptible to methicillin) and S. aureus ATCC 43300 (resistant to methicillin) were used as negative and positive controls, respectively. Antibiotic resistance profile of the MRSA isolates was determined by broth microdilution according to the European Committee on Antimicrobial Susceptibility testing breakpoints [94]. Methicillin-sensitive strain, S. aureus ATCC 29213 was used as a control. All MRSA isolates were the subjected to WGS sequencing and analysis in order to ascertain the dominant clonal lineage, its pathogenicity, virulome, environmental tolerance, persistence mechanisms and determine their possible putative threat for human health using a comprehensive bioinformatic analysis.

4.2. Whole Genome Sequencing (WGS) Analysis and Characterisation

4.2.1. Purification, Sequencing and Pre-Processing of Genomic Data

One colony-forming unit from a visibly pure culture of each MRSA isolate was selected for WGS. Genomic DNA (gDNA) of the isolates was extracted and purified using the GenElute Bacterial Genomic DNA kit (Sigma Aldrich, St. Louis, MO, USA) per the manufacturer’s instructions. Following extraction, quantification was performed on a Nanodrop 8000 (Thermo Scientific, Waltham, MA, USA) with verification by agarose gel electrophoresis [95]. A paired-end library was prepared using Nextera XT DNA Sample Preparation Kit and the whole-genome sequencing was carried out on an Illumina MiSeq machine (Illumina, San Diego, CA, USA). The sequenced reads were quality trimmed using Sickle version 1.33 (https://github.com/najoshi/sickle) and de novo assembled using SPAdes version 3.11 [96] and the CLC Genomics Workbench version 10.1 (CLC, Bio-QIAGEN, Aarhus, Denmark). All resultant contiguous sequences were then submitted to GenBank for gene annotation using the NCBI Prokaryotic Genome Annotation Pipeline [97].

4.2.2. WGS-Based Molecular Typing of the Obtained MRSA

The SCCmec type and its structural composition in the MRSA isolates were determined in-silico using the web-based tool SCCmecFinder, freely available at https://cge.cbs.dtu.dk/services/SCCmecFinder [98]. Spa typing of MRSA isolates was performed in-silico using the assembled genomic sequences on the online platform tool Spa Typer 1.0 [99]. Multilocus sequence typing (MLST) typing was performed in-silico using the WGS data online platform tool MLST 1.8 [100]. An eBURST [101] analysis was performed on all sequence types (STs) of S. aureus in the MLST database (http://saureus.mlst.net). STs were assigned to clonal complexes (CC) where they had six identical alleles with at least one other ST within the clonal complex.

4.2.3. In-Silico Resistome Profiling

The Comprehensive Antibiotic Resistance Database platform (https://card.mcmaster.ca/analyze/rgi) [102] and the ResFinder through the GoSeqIt tools web-platform (https://www.goseqit.com/web-services/) [103] were used for the prediction of resistance genes. To detect the molecular basis of resistance (chromosomal SNPs) against quinolones (gyrA, parC, parE) and rifampicin (rpoB), the nucleotide allele sequences were translated with tBLASTn to call SNPs in these genes the using the S. aureus ATCC 29212 (Accession no. CP009361) as the wild-type strain.

4.3. Genome Visualization and Gene Annotation

The raw reads were de-novo assembled using the SPAdes assembler [96]. The genomes of the strains were visualized using the CG Viewer Server (http://stothard.afns.ualberta.ca/cgview_server/index.html) [104] (Figure 1). NCBI Prokaryotic Genome Annotation Pipeline (PGAP) available at http://www.ncbi.nlm.nih.gov/ [97] and SEED subsystems in the RAST server (rapid annotation using subsystem technology) available at http://rast.nmpdr.org/ [105] were used to annotate these genomes. The size, GC content, number of contigs, N50, L50, average coverage and the number of RNAs and protein coding sequences obtained for each isolate. The annotated functional category of in-silico-predicted proteins was also visualised via the SEED subsystems in the RAST server. A progressive alignment algorithm implemented in MAUVE was used to determine the rearrangements in each genome [106]. The Tandem Repeat Finder available at https://tandem.bu.edu/trf/trf.html [107] was used to analyze the DNA sequences of the 11 isolates to predict repeats in the genome.

4.4. Detection of CRISPR Array, Restriction-Modification System (R-M system), Accessory Gene Regulator (agr) Type, Arginine Catabolic Mobile Element (ACME)

The CRISPRCasFinder available at https://crisprcas.i2bc.paris-saclay.fr/CrisprCasFinder/Index [108] was used to identify putative CRISPR loci and Cas cluster in the draft genomes. Annotations from the Pathosystems Resource Integration Center (PATRIC) online platform and Restriction Modification Finder at https://cge.cbs.dtu.dk/services/Restriction-ModificationFinder/ predicted the R-M system in the isolates. The alignment of fully annotated reference agr types (I, II, III and IV) was used to identify the specific agr type in all the isolates using an identity match and query length of ≥90%. A similar alignment of the ACME components (arcA and opp3AB) was used to predict the specific ACME type and the classification done as follows: ACME type I (arcA+/opp3AB+), II (arcA+/opp3AB−) and III (arcA−/opp3AB+), where positive and negative indicated presence and absence, respectively.

4.5. Genomic Insights of the Isolates in the Endemic Clone

4.5.1. Pathogenicity and Virulome Predictions

PathogenFinder available at https://cge.cbs.dtu.dk/services/PathogenFinder/ [109] was used to predict pathogenicity towards human hosts. Furthermore, the P score sister lineages (ST8) and well-known epidemic MRSA clones (ST5, ST8, ST59, ST239) were computed to offer a better comparison of the pathogenicity of this clone. Virulence determinants (sequences and functions) corresponding to different major bacterial virulence factors (adherence, antiphagoctyosis, exoenzyme, immune evasion, iron uptake, plasminogen activator, secretion system and toxin) associated with S. aureus were collected from GenBank. The predicted factors were then validated using virulence factors of pathogenic bacteria database, VFanalyzer (available at http://www.mgc.ac.cn/VFs/). The known epidemic S. aureus subsp. aureus USA300_TCH1516, 2,872,915 bp, NC_010079) which is also close relative of the clone (ST612-CC8-t1257-SCCmec_IVd(2B)) was used as the reference genome for this inference. Furthermore, the Victors virulence factors search program (available at http://www.phidias.us/victors/) [110] and PATRIC_VF tool [111], were used to overrule the inherent shortfalls of all the tools.

4.5.2. Genomic Prediction of Mechanisms of Bacterial Persistence and Tolerance

The SEED subsystems in the RAST server and PATRIC database platform were used to identify and profile genes putatively associated with tolerance and persistence such as the stress response (osmotic, oxidative and periplasmic), heat shock, dominance and sporulation, bacteriocins and detoxification in the clone from diverse sources.

4.6. Comparative Phylogenomic Analysis and Metadata Insights

An all-by-all BLAST phylogenomic comparison was performed with the endemic clonal subtype ST612-CC8-t1257-SCCmec_IVd(2B) genome downloaded from the NCBI database (https://www.ncbi/) via their accession numbers to investigate phylogeny. The genomic sequences of three close relative strains of the endemic clone, which are members of clonal complex 8 (CC8), were rooted as references for the comparative phylogenomic analysis. The three strains were: S. aureus USA300 strain TCH1516 (Accession number: CP000730), S. aureus strain 2395, USA500 (Accession number: CP007499.1) and S. aureus strain SVH7513 (Accession number: CP029166) (the first complete genome of the ST612; a livestock-associated MRSA strain). Nucleotide sequences of all strains were respectively aligned using the default parameters/settings of classical sequence analysis of the CLC Genomics Workbench (version 10.1.1) to generate an aligned file (CLC Genomics Workbench 11.0.0; https://www.qiagenbioinformatics.com/)). The created aligned file was used to draw the maximum likelihood phylogenetic tree to infer the evolutionary relationship using optimized parameters (Construction method: unweighted pair group method with arithmetic mean (UPGMA), nucleotide substitution model: Jukes cantor, protein substitution model: WAG, transition/transversion ratio: 2, estimate substitution rate: yes, number of substitution rate: yes, number of substitution rate: 4, perform Bootstraps analysis: Yes, Replicates: 1000) of the CLC Genomics Workbench (version 11.0.0) [112]. The average nucleotide identity (ANI) was used to provide a robust measurement of genetic distance among the genomes of the three reference strains and the endemic clone, for the conserved genes of the genomes using the ANI calculator (https://www.ezbiocloud.net/tools/ani) [113,114]. Additionally, an SNP-based phylogenomic tree was generated with the endemic clone via the CSI Phylogeny-1.4 (https://cge.cbs.dtu.dk/services/CSIPhylogeny-1.2). The obtained phylogenomic tree was downloaded in Newick format and annotated, visualized or edited using an interactive tree of life (ITOL) (https://itol.embl.de/). The edited trees were coupled with their metadata (genomic profiles of putative virulence determinants) via ITOLs to generate heatmaps [115].

5. Conclusions

The genetic insights into the dominance of the MRSA endemic clone ST612-CC8-t1257-SCCmec_IVd(2B) revealed a battery of highly conserved defence systems, putative virulence determinants, environmental tolerance and persistent mechanisms which enable it to withstand endogenous and exogenous lethal threats and contribute to its ecological success in diverse biological niches. The myriads of genomic signatures and mechanisms possessed by this dominant clone are potential targets for drug development necessitating further studies on transcriptomics and functional genomics to inform the design of novel therapeutic strategies to curb this pathogen.

6. Ethical Considerations

Ethical approval was received from the Animal Research Ethics Committee (Reference: AREC 073/016PD) and the Biomedical Research Ethics Committee (Reference: BCA444/16) of the University of KwaZulu-Natal. The study was further registered with the South African National Department of Agriculture, Forestry and Fisheries (Reference: 12/11/1/5 (879)). Human samples were obtained from participants 18 years or older upon explicit, voluntary, written informed consent. All additional information obtained from the farm (herein, noted as Farm A) were kept confidential as part of the memorandum of understanding (MOU) between the Antimicrobial Research Unit (ARU) and the farm.