Inside environmental Clostridium perfringens genomes: antibiotic resistance genes, virulence factors and genomic features

Clostridium perfringens (previously known as Bacillus welchii) is a strictly anaerobic, Gram-positive, rod-shaped bacterial pathogen (Hassan et al. 2015). This rapidly-growing bacterium has the ability to form endospores that are extremely resilient to toxic chemical and environmental stressors, such as heat or radiation (Paredes-Sabja et al. 2008). As a result, C. perfringens inhabits diverse environments including soil, water, sewage, food and even the gastrointestinal tract of mammals (Sathish & Swaminathan 2009). In the gastrointestinal tract, Clostridium perfringens does not infect healthy cells directly but rather acts by producing toxins and enzymes that cause a number of enteric and systemic diseases in the infected host (Gross et al. 1989; Kiu & Hall 2018).

Pathogenic C. perfringens strains have been found to produce and secrete more than 20 different extracellular toxins and/or hydrolytic enzymes (Revitt-Mills et al. 2015). They are classified into seven groups (toxinotypes A–G) according to the combination toxin types they produce (α-toxin, β-toxin, ɛ-toxin and ι-toxin, enterotoxin (CPE) and NetB) (Kiu & Hall 2018). Clinically, the toxinotype A strain is the most common disease-causing agent responsible for gas gangrene (clostridial myonecrosis), necrotic enteritis, as well as mild diarrhoea (Brynestad & Granum 2002). Although these toxins and enzymes each play a specific role in the disease process, it is the synergistic actions of both on the host that could potentially be the key virulence factors involved in its pathophysiology (Rood 1998; Revitt-Mills et al. 2015).

Beside virulence factors in Clostridium perfringens, another trait such as antibiotic resistance poses a danger to clinically treated patients, mainly due to their ability to render the antibiotic treatment of severe and life-threatening infections ineffective (Kiu & Hall 2018). Therefore, extensive research has been conducted on antibiotic resistance profiles in C. perfringens using phenotypic methods (e.g. testing of minimal inhibitory concentration (MIC); Koch et al. 1998; Akhi et al. 2015). However, there are only a few studies that have used genomic methods (i.e. WGS, whole-genome sequencing) to investigate antibiotic resistance genes (ARGs; Kim et al. 2017; Li et al. 2017).

It is considering the clinical relevance of Clostridium perfringens holds and its wide distribution in nature; the current state of knowledge on the pathogenicity of C. perfringens or its non-clinical environment is still poor when compared with other more well-known pathogens (Shimizu et al. 2002). New technologies, such as next-generation sequencing (NGS), have become the new gold standard in in silico analysis. This is because NGS helps to understand bacterial pathogenesis, identify and characterise genes coding for virulence factors, toxins or antibiotic resistance in pathogens, especially C. perfringens strains which are hard to grow under laboratory conditions (Bakour et al. 2016). Thus, the aim of this study was to expand on the existing knowledge of Clostridium perfringens through the analysis of assembled genomes from surface water and identify the presence of virulence and ARGs.

Clostridium perfringens was isolated from water of the Schoonspruit River (26°40′46.0″S, 26°34′58.7″E) in the North West Province (South Africa) in March 2016. A modified Fung's double tube method, along with tryptose–sulfite–cycloserine agar (Oxoid, UK), was used for the bacterial isolation from collected samples (Fourie et al. 2019). The isolates were purified through streaking and grown under anaerobic conditions at 44 °C and subjected to antibiotic susceptibility testing as recommended by the Clinical and Laboratory Standards Institute (CLSI 2016). The MICs of five clinically relevant antibiotics (ampicillin, tetracycline, clindamycin, chloramphenicol and metronidazole) were determined by the agar dilution method. Thereafter, multidrug-resistant (MDR) isolates were selected for genomic DNA extraction. gDNA extraction from overnight cultures was done using the NucleoSpin^® Tissue kit (Macherey-Nagel, Germany) according to the instructions provided by the manufacturer. Species confirmation was done by amplifying the 16S rRNA gene, followed by Sanger sequencing of the gene fragment. Subsequently, three MDR C. perfringens isolates with high DNA purity were selected for whole-genome sequencing.

DNA libraries were performed using the Nextera XT kit (Illumina Inc., San Diego, CA) according to the instructions provided by the manufacturer and sequenced using Illumina's MiSeq300 in paired-end reads. All short (<50 bp) or low quality value (QV <15) fragments were removed. The quality assessment and trimming were done in Trimmomic (v.0.36; Bolger et al. 2014), followed by de novo assembly in SPAdes (v.3.9.0) to generate scaffolds and contigs (Bankevich et al. 2012). Gene prediction and annotation were conducted using the NCBI Prokaryotic Genome Annotation Pipeline (v.4.3) along with the Rapid Annotation using Subsystem Technology (RAST) server (Overbeek et al. 2014; Tatusova et al. 2016). Subsequently, the obtained draft genome sequences used in this study were deposited into the NCBI database (Fourie et al. 2019).

ARGs and virulence factors were identified within the genomes of all three C. perfringens strains by using BLASTx comparison to the virulence factors database (VFDB) and the ARGs database in deepARG (v.2.0; Chen et al. 2016; Arango-Argoty et al. 2018). The presence of genomic islands (GIs) within the assembled genomes was predicted using the VRprofile (v.2.0; Li et al. 2018). Prophage regions were identified and characterised and located using PHASTER (PHAge Search Tool – Enhanced Release; Arndt et al. 2016). Default parameters were used for all software unless otherwise specified.

Average nucleotide identity (ANI) was determined by using the OrthoANI algorithm (v.1.4) to analyse genomic relatedness between the three analysed and referenced C. perfringens strains. OrthoANI percentages were calculated, and a UPGMA dendrogram was constructed (Lee et al. 2015). Annotated protein sequences were submitted to OrthoVenn2 to identify unique and/or shared orthologous clusters among three C. perfringens genomes, as well as three reference genomes (Xu et al. 2019). Reference strains used for comparison were Clostridium perfringens: ATCC 13124 (Genbank: NC_008261), Str 13 (Genbank: BA000016) and FORC_025 (Genbank: NZ_CP013101).

The three assembled Clostridium perfringens strains (SC4-C13, SC4-C17 and SC4-C24) that were originally obtained from surface water produced genome sizes of 3.6, 3.5 and 3.4 Mbp, respectively (Table S1) and were organised in 205 (both SC4-C13 and SC4-C17) and 110 (SC4-C24) contigs (Fourie et al. 2019). The GC content between the three strains averaged at around 28.1%. Annotation of the genomes indicated the presence of open reading frames (ORFs) ranging between 3,079 and 3,245 for all three strains (Table S1), which is consistent with published data on this bacteria. The average genome size known for C. perfringens varies from 3.0 to 4.1 Mbp, with a GC content ranging between 27 and 28% (Kiu et al. 2017; Kiu & Hall 2018). The number of coding sequences present in their genome averaged between 2,500 and 3,600 (Kiu & Hall 2018).

RAST was used to functionally annotate the open reading frames of each genome into 25 subsystem categories (Figure 1). Based on the categories generated, all three strains showed similar amounts of assigned ORFs for each category. The majority was allocated to ‘carbohydrates’, thus showing the importance of carbohydrates for these pathogens, since they play a vital role in their growth and sporulation (Sacks 1983). This was then followed by ‘protein metabolism’, ‘amino acids and derivatives’ and ‘cofactors, vitamins, prosthetic groups and pigments’. The less abundant ORFs belonged to the ‘metabolism of aromatic compounds’ category, along with ‘motility and chemotaxis’, ‘transposable elements’ and ‘nitrogen metabolism’.

The deepARG tool applies a deep-learning approach to identifying and annotating ARGs from genomic data. It does this by creating a dissimilarity matrix based on ARG categories and comparing it against the three major datasets currently available, namely Comprehensive Antibiotic Resistance Database (CARD), Antibiotic Resistance Genes Database (ARDB) and Universal Protein Resource (UNIPROT; Arango-Argoty et al. 2018). This approach identified a total of 11 ARGs in the genomes of the three Clostridium perfringens strains (Table 1). Strains SC4-C13 and SC4-C17 both possess genes that encode for resistance against seven different classes of antibiotics (glycopeptide, β-lactamase, macrolide–lincosamide–streptogramin B (MLSB), tetracycline, trimethoprim, kasugamycin and bacitracin), whereas strain SC4-C24 only harboured resistance genes against six classes, with the absence of the β-lactam resistance gene (bla2). Furthermore, all three strains hold two or more genes that encode for multidrug resistance.

β-lactam antibiotics are the most commonly used class of antibiotics in human health medicine and agriculture (Price et al. 2019). The prolonged misuse of these antibiotics have resulted in the more frequent occurrence of β-lactam resistance in bacteria. This, in turn, has diminished its effectiveness in fighting against bacterial infections that were once easily treated. The bacterial mechanism of β-lactamase resistance is primarily due to the production of β-lactamase enzymes that hydrolyse the β-lactam ring, thereby inactivating the antibiotic (Bush & Bradford 2016). The presence of the bla2 gene which encodes for β-lactamase in both SC4-C13 and SC4-C17 is worrisome, since most β-lactam resistance is not commonly associated with Clostridium species, with only a few reported incidences of resistance in clinical settings (Finegold et al. 2005; Mishra et al. 2016). Based on the UNIPROT database, the bla2 gene was originally identified in Firmicutes bacterium CAG:212, which could suggest that the two C. perfringens strains acquired this gene from an external source and incorporated it into their genomes.

Glycopeptides encompass a group of antibiotics that are vital in combating infections caused by antibiotic-resistant bacteria, with vancomycin still the first choice of treatment against Gram-positive bacteria since its introduction 50 years ago (Blaskovich et al. 2018). Vancomycin works by binding to the d-Ala–d-Ala C-terminal of the growing pentapeptide chain during cell wall synthesis and inhibits further elongation and cross-linking of the peptidoglycan chain blocking through transglycosylation and transpeptidation (Sanakal & Kaliwal 2011). However, resistance to this antibiotic has become a common occurrence, especially in vancomycin-resistant enterococci and even in methicillin-resistant Staphylococcus aureus (Willems et al. 2005). The resistance in these bacteria has been confirmed by the presence of several vancomycin resistance genes (van). Our strains possessed many of these genes (vanHD, vanTrL, vanRB, vanRG and/or vanRI), which encode for enzymes that produce synthetic precursors that replace the d-Ala–d-Ala C-terminal, altering the vancomycin-binding site (Sanakal & Kaliwal 2011). Quite a few studies have reported vancomycin resistance genes in Clostridium species (Peltier et al. 2013; Chia et al. 2017). However, Clostridium perfringens has not been one of the few to show such resistance since they are normally susceptible to vancomycin (Citron et al. 2005; Camacho et al. 2008; Chia et al. 2017). The presence of these genes could suggest a new resistance trait in these species if they are expressed.

A further study of the results presented in Table 1 shows that all three C. perfringens genomes harboured the same gene that encodes resistance for MLSB (ermQ), tetracycline (tetB(P)), trimethoprim (dfrK), kasugamycin (KsgA) and bacitracin (bacA). MLSB antibiotics work by inhibiting protein synthesis in both Gram-positive and Gram-negative bacteria. Resistance to this antibiotic is due to 23S rRNA methylation (encoded by erm genes) that prevents MLSB antibiotics from binding to the ribosome. Several erm genes have been detected in Clostridium species, such as C. perfringens and C. difficile (Spigaglia & Mastrantonio 2002). However, the ermQ gene was first identified in C. perfringens and found to be a distinct class of MLSB resistance determinant (Berryman et al. 1994; Shoemaker et al. 2001).

Tetracyclines are broad-spectrum antibiotics and the second most used antibiotics group after β-lactams. Bacterial resistance to tetracycline is achieved by harbouring one or more of the 36 known tet genes (Sheykhsaran et al. 2019). These tet genes follow one of three resistance mechanisms: energy-dependent efflux pumps, ribosomal protection proteins, or drug target alteration and enzymatic inactivation (van Hoek et al. 2011). The tetB(P) genes, found in our strains along with tetA(P) (not present), are originally two overlapping genes within the Tet P determinant of C. perfringens (Sloan et al. 1994). However, because the Tet P determinant has only been found in Clostridium species and is often associated with conjugative and non-conjugative plasmids, it has shown the capacity to spread rapidly across the whole genus (Roberts 2011; Vidor et al. 2019).

Clostridium species are known to show endogenous resistance to trimethoprim (Huovinen et al. 1995). However, the presence of the drfK gene in the three C. perfringens genomes suggests that it could have been acquired from an external source since it has mainly been associated with enterococci and staphylococci species (López et al. 2012). While the drfK gene may not be vital for trimethoprim resistance in Clostridium, its presence may have implications for further dissemination in the environment, or to other bacteria that are primarily susceptible to trimethoprim.

Antibiotics such as kasugamycin and bacitracin are widely used in agriculture and animal husbandry as growth promoters, with bacitracin also being used as a topical ointment for skin infections (Duffin & Seifert 2009; Charlebois et al. 2012). Resistance to the aminoglycoside antibiotic kasugamycin is due to the presence of the KsgA gene, which produces a predicted dimethyltransferase, whereas an ATP-dependent ABC-type efflux system consisting of three BcrA, -B and -C proteins is responsible for bacitracin resistance (El Ghachi et al. 2004; Duffin & Seifert 2009). Studies have shown that the use of antibiotics in agriculture and growth promoters in animal breeding have a direct correlation to the development and dissemination of antibiotic-resistant bacteria and genes into the surrounding water sources. This happens by means of surface runoff, where antibiotic resistance can be further intensified via horizontal gene transfer (HGT; Kümmerer 2004).

Genes for efflux pumps (mepA), ATP-binding cassette transporters (vgaB) and system regulators (arlR), which are associated with multidrug resistance, were found in our genomes as well. The arlR is part of a two-component regulatory system, a transmembrane sensor and its associated response regulator (arlS-arlR). This regulator was found to increase the expression of efflux pumps associated with Staphylococcus aureus (Fournier et al. 2000). The vgaB gene is one of the ATP transporter genes and can increase the level of resistance to pristinamycin, a mixture of streptogramin A and streptogramin B compounds via efflux pumps (Roberts 2002; Chesneau et al. 2005). The mepA encodes for a multidrug resistance efflux pump in S. aureus, which has reported low-level resistance to quaternary ammonium compounds, such as chlorhexidine, pentamidine and tigecycline (Costa et al. 2013).

When comparing the genotypic characteristics with the phenotypic profiles (Table S2 in Supplementary Materials) of these three strains, some conclusions can be drawn. All the strains were shown to be susceptible to chloramphenicol, along with the absence of any associated resistance genes to support this. Two of the strains, SC4-C13 and SC4-C17, showed phenotypical resistance to ampicillin (β-lactam) and possessed the bla2 gene, whereas SC4-C24 showed intermediate resistance with no genes to confirm resistance in its genome. Mobile genetic elements such as plasmids have been shown to carry β-lactam resistance genes (Partridge et al. 2018). If such plasmid-mediated resistance genes were to be present in C. perfringens and expressed, it could translate to phenotypical resistance to β-lactam antibiotics. Alternatively, this low level of resistance might also suggest an increase in the CLSI-defined susceptibility breakpoint for ampicillin in C. perfringens. Noticeably, all three strains also showed phenotypic resistance to tetracycline, clindamycin and metronidazole. The presence of genes such as ermQ and tetB(P) could explain the phenotypical resistance to clindamycin and tetracycline, respectively. However, no genes associated with metronidazole were identified in any of the assembled genomes. Previous studies have reported phenotypic resistance to metronidazole in C. perfringens (Tansuphasiri et al. 2005; Akhi et al. 2015). However, the understanding of its resistance mechanism in Clostridium species is still unclear (Lynch et al. 2013). It is also important to add that although the assembly delivered almost near genomes, the lack of genetic confirmation for phenotypic resistance may be due to incomplete gene reconstruction of the region where the gene occurs.

The genes that are associated with virulence in Clostridium perfringens strains were analysed using the VFDB as a reference. A total of 35 different virulence genes were successfully annotated with this database, of which 79% are typically associated with Clostridium perfringens strains. These include genes that are responsible for adherence, such as fibronectin-binding protein, chaperonin GroEL and type IV pili, as well as the two members of the double-component VirR/VirS regulon (sensor histidine kinase and DNA-binding response regulator). There are also 19 genes present in all three strains that are capable of producing eight different toxins, namely α-clostripain, phospholipase C, collagenase, hyaluronidase, perfringolysin O, sialidase and haemolysin (Table 2). Interestingly, the remaining 21% of virulence genes are usually associated with other Gram-positive bacteria, namely Listeria ivanovii, Streptococcus pyogenes, Mycoplasma penetrans and Enterococcus faecium. This might imply their potential influence or contribution to the expanding pathogenicity in these three strains as bile-resistance and adherence factors (Table 2). Among the three C. perfringens strains described in this study, SC4-C13 and SC4-C17 only possess 94.7% of the above-mentioned virulence factors, while for SC4-C24 it is 100%.

Clostridium perfringens can produce a variety of more than 20 different extracellular toxins and hydrolytic enzymes, giving it the ability to cause various histotoxic infections in humans and animals (Kiu & Hall 2018). However, not all strains are able to produce all these toxins. This is mainly due to some toxins being strain-specific (Kiu et al. 2017). The C. perfringens strains in this study possess genes able to code a vast array of toxins and enzymes. However, only a combination of several typing toxins is traditionally used to determine the toxinotype of a species. All three strains were identified as toxinotype A due to the presence of the plc gene, which encodes for α-toxin and the absence of the other typing toxin genes. Furthermore, C. perfringens toxinotype A strains are known to be human pathogens, causing diseases such gas gangrene (clostridial myonecrosis), necrotic enteritis, as well as mild diarrhoea (Brynestad & Granum 2002). The presence of perfringolysin O (θ-toxin) has also shown to have a synergistic effect with α-toxin in the pathology of gas gangrene. Even on its own, this cholesterol-dependent cytolysin can form pores on cell membranes and lyse red blood cells, further highlighting the significant role perfringolysin O has in disease development (Awad et al. 2001; Kiu & Hall 2018). However, these two toxins, along with the κ-toxin gene (colA), are tightly regulated by specific regulatory systems in C. perfringens (Ohtani & Shimizu 2016). The VirS/VirR two-component system is one of the most important systems and consists of two genes: the response regulator (virR) and the sensor histidine kinase (virS; Ma et al. 2011). Still, this system assumes responsibility for coordinating the pathogenicity of C. perfringens type A strains.

The hydrolytic enzymes produced by C. perfringens also play an important role during C. perfringens-mediated tissue infection. Sialidase is encoded by several genes that contribute to internal virulence by improving the adherence of C. perfringens to host cells, enhancing the production, binding and activity of certain toxins that are responsible for intestinal infections, and contributes to intestinal colonisation by C. perfringens (Li et al. 2016). Mu-toxins (hyaluronidase), on the other hand, are enzymes that facilitate the breakdown of hyaluronate substrates, improving contact between C. perfringens and the specific cell wall receptors. It has also been found to promote the spread of α-toxin and, in doing so, potentiate its cytolytic activity (Hynes & Walton 2000).

As in most pathogens, adherence plays an essential role in the pathogenesis of Clostridium perfringens. Several adherence genes were identified in the genomes of C. perfringens, each of which encodes for a different approach of adherence and promotes colonisation. Fibronectin-binding protein (fbpA), for instance, makes it possible for C. perfringens to bind to fibronectin found in all human tissues and organs, including intestinal epithelial cells, and in doing so, invade the host cells (Katayama et al. 2009). Whereas type IV pili are elongated, flexible filaments extend from the bacterial cell surface and are implicated in attaching to and invading host cells, attachment to abiotic surfaces, biofilm formation and bacteriophage susceptibility (Rodgers et al. 2011).

Genomic islands (GIs) are regions of the bacterial genome that are acquired through HGT (Dobrindt et al. 2004; Klein et al. 2018). Antibiotic resistance islands (ARIs) and pathogenicity islands (PAIs) are two subgroups of GIs that aid and contribute to the pathogenesis of an organism (Yoon et al. 2015). Figure 2(a) shows the distribution of ARIs and PAIs present in the three Clostridium perfringens strains. All strains analysed possessed PAIs and ARIs homologous to those found in most enteropathogens. A total of 65 ARIs were identified in all three strains that were related to genes from Staphylococcus aureus (24/65), Acinetobacter baumannii (19/65), Pseudomonas aeruginosa (10/65), Campylobacter coli (3/65), Shigella flexneri (3/65), Klebsiella pneumoniae (2/65), Proteus mirabilis (2/65), Escherichia coli O157H7 (1/65) and Staphylococcus hominis (1/65) (Table S3). This included genes which encode for putative enzymes, such as acyl-CoA dehydrogenase, DNA topoisomerase III and enoyl-CoA hydratase/isomerase, as well as regulatory proteins such as KDP operon transcriptional regulatory protein (KdpE), propionate catabolism operon regulator (prpR), putative transcriptional and response (ArmR) regulators, to name a few. Genes encoding for proteins and enzymes contributing to antibiotic resistance (TetA class A, truncated tetracycline-resistant protein and streptothricin acetyltransferase) and tolerance to heavy metals (mercuric reductase and cadmium efflux system proteins B and C) were also identified. Heavy metals such as mercury and cadmium are commonly used for antimicrobial metal compounds in healthcare (Pal et al. 2017). However, the presence of genes encoding for resistance to antibiotics and heavy metal might suggest possible adverse implications for infectious cases of C. perfringens. Iron plays a key role in the growth and survival of many pathogens in any niche. It is utilised for important biological processes such as DNA synthesis, generating energy and protection against reactive oxygen species (Choo et al. 2016). All the C. perfringens strains possess ARIs related to Shigella flexneri that encodes for proteins responsible for the binding and transport of ferric iron (III) (FecC, FecD and FecE; Table S3). Genes associated with iron acquisition have previously been linked to GIs in Clostridium perfringens (Myers et al. 2006).

The distribution of all PAIs identified among the C. perfringens strains are visualised in Figure 2(b). The various PAIs shared among the three strains ranged from hypothetical or unnamed proteins related to Yersinia pestis, Escherichia coli, Edwardsiella tarda and Streptomyces lividans (Table S4). Other examples include PAIs encoding for Epd, integrase, Orf17, Orf29, Orf34, Orf41, Orf59, Pgk and TktA from Photorhabdus luminescens; dihydrofolate reductase type Ib (plasmid encoded), dihydropteroate synthase, GlcD protein, GspE type II secretion protein, haemolysin B, phosphoglycerate transport protein, putative ABC transporter ATP-binding and membrane proteins, putative lysil-tRNA synthetase (LysU) and short-chain dehydrogenase from Escherichia coli; ATPases with chaperone activity ATP-binding subunit, serine/threonine protein phosphatase and transcriptional regulator from Vibrio alginolyticus and glutamate racemase from Helicobacter pylori. Additionally, sequences similar to ones encoding for HrcN and HrpB6 in Xanthomonas species were also present. These two genes encode for the type III secretion system (TTSS) proteins, a requirement for the pathogenicity of several plant pathogens (Villela et al. 2019).

Among the three Clostridium perfringens strains, SC4-C13 and SC4-17 had the most unique PAIs and ARIs (Tables S3 and S4). Noticeably, both these strains harboured genes encoding for putative transposase. Strain SC4-C13 possessed multiple transposases for IS285 (Y. pseudotuberculosis), ISSau4-like and IS256/Tn4001 (Staphylococcus aureus), whereas SC4-C17 harboured a transposase protein for IS30 (Staphylococcus aureus). Additionally, all strains possessed two or more GIs containing integrase genes from Photorhabdus luminescens, Acinetobacter baumannii or Staphylococcus aureus. The presence of genes associated with transposases and integrases has shown to infer high potential for genetic gain or loss in C. perfringens genomes (Kiu et al. 2017).

A previous study identified over 300 GIs in the genomes of several clinical C. perfringens type A strains, of which almost all were chromosomally encoded (Myers et al. 2006). Less is known about the occurrence of GIs in strains isolated from environmental sources (Klein et al. 2018). Pathogens such as C. perfringens are found in various ecological niches, which also include diverse microbial communities that can contribute to HGT events (Kiu et al. 2017).

The Phage Search Tool Enhanced Release (PHASTER) program identifies and characterises prophage regions into three categories based on the completeness scores, namely intact (>90), questionable (70–90) or incomplete (<70). All three Clostridium perfringens strains used in this study contained two or more prophage regions in their genome (Table 3). Strain SC4-C24 were shown to harbour three prophage regions, of which one was confirmed to be intact, and the other two showed incomplete and questionable phage-related sequences. Another intact prophage was identified in strain SC4-C17, along with an additional incomplete region. Strain SC4-C13 showed putative fragments of two phages, one of which scored relatively high (90). Each prophage region varies in size, the number of CDS and GC content. The two intact prophages in SC4-C17 and SC4-C24 strains have sizes of 49.6 and 42.2 Kbp, respectively, with both containing the second-highest number of CDS (57). Functional annotation of the genes present in these prophage regions showed protein sequences encoding for phage-related elements such as integrase, terminase, portal, capsid, head, tail and/or transposase (Table S5). Phages play a very important role in the evolution of bacterial genomes (Ramisetty & Sudhakari 2019). Their interactions with each other have impacted both their survival and persistence. This is because prophages have several advantages to offer in the way of inversion, deletion and insertion via HGT of genetic material (Darmon & Leach 2014). This is because prophages are known to be ‘hotspots’ for carrying significant genes such as those involved in virulence, antibiotic resistance and metabolic pathways (Ramisetty & Sudhakari 2019). This led to further investigation into functional protein sequences in the identified prophage regions of the three C. perfringens genomes. However, only the presence of hypothetical proteins was discovered in these regions and could suggest that these genes encode for unknown functions from remote sources.

The average nucleotide identity for all Clostridium perfringens genomes was successfully calculated using OrthoANI (Figure 3). ANI analyses of the three strains revealed them to be almost identical to each other as indicated by the ANI of >99% between the strains. Although the differences between the three strains appear negligible, based on the dendrogram, strains SC4-13 and SC4-C17 share relatively more similarities between them, as opposed to strain SC4-24 (Figure 3). When comparing the three C. perfringens genomes with that of the reference genomes (ATCC 13124, Str 13 and FORC_025), another high similarity was observed among the genomes, with our genomes sharing ANI values of 98.49–98.63% with the reference genomes. The dendrogram also showed that the three strains were clustered more closely to FORC­_025, followed by the ATCC 13124 strain which was derived from clinical environments (Myers et al. 2006; Kiu et al. 2017), whereas only strain 13 was originally isolated from a natural environment (soil; Shimizu et al. 2002).

Orthologous genes originate from a single gene present in the last common ancestor through a series of speciation events and usually still have the same biological function in the present-day organism (Xu et al. 2019). Therefore, for further analysis of the three C. perfringens strains, a multi-species comparison of the shared core orthologous gene clusters was performed without (Figure 4(a)) and with reference genomes (Figure 4(b)). Results show a large number of (2,779) core orthologous genes shared between the three C. perfringens strains (Figure 4(a)). These clusters were all revealed to be associated with biological processes, molecular functions and cellular components within their core genomes (Table S6). The orthologous analysis also supported the higher similarity between SC4-C13 and SC4-C17 based on them sharing over 130 orthologous clusters, while SC4-C24 has the highest amount of unique orthologous clusters of the three strains. When comparing the core orthologous clusters of the three strains with those of the reference strains, they were shown to share 2,300 clusters (Figure 4(b)). Altogether, only eight orthologous clusters were identified as unique among all the genomes. Interestingly, these clusters were found in strains SC4-C13 and SC4-C24. Strain SC4-C13 had two of the eight unique clusters, one of which was responsible for glycosyltransferase and the other for the biological response (movement, secretion, enzyme production or gene expression) to external cold stimuli of the cell. Strain SC4-C24 possessed the other six unique clusters, three clusters were responsible for the initiation of DNA-template transcription, hydrogen peroxide catabolic and antigen biosynthetic processes, while the other three clusters did not have assigned functions. Although only two of the three C. perfringens strains showed unique orthologous clusters when compared with the reference genomes, this may be attributed to the possibility of isolates niche adaptation (Braga et al. 2016).

The current study has shown that Clostridium perfringens strains from a water environment posses similar genomic features to published clinical C. perfringens. The in silico analysis of the three genomes revealed putative genes encoding for virulence factors normally responsible for pathogen adherence to host cells, production of toxins and extracellular enzymes and the presence of key system regulators managing the pathogenicity of C. perfringens. The strains were also classified as toxinotype A, which commonly consists of human pathogens only. Furthermore, the presence of ARGs encoding for several classes of antibiotics and multidrug resistance efflux pumps could contribute to the emergence of pathogenic multidrug-resistant C. perfringens strains in this environment. Interestingly, the presence of several genes normally associated with other pathogenic genera, as well as a few intact prophages and genomic islands (ARIs and PAIs), might suggest the transfer of genetic material between bacteria through HGT within the water reservoir, thus indirectly implicating the presence of toxic substances in the water (cadmium, copper, mercury etc.), which are associated with HGT in bacteria. Moreover, expression analyses should be the next step of the investigation. Furthermore, comparative genomics showed great similarities between the environmental C. perfringens strains with those of clinical origin, therefore highlighting the importance of monitoring these bacteria for effective epidemiological surveillance purposes. Genomic studies on Clostridium perfringens in natural environments are currently very limited and should be explored further.

This work is based on research supported in part by the Water Research Commission (WRC) of South Africa (contract no. K5/2347//3), National Research Foundation (NRF) of South Africa (UID no. 109207), some NWU funding. We also thank the Centre for High Performance Computing (programme no. CBBI0890) in Cape Town (South Africa) for access to computing resources. Views expressed are those of the authors and not of the funders.

The Supplementary Material for this paper is available online at https://dx.doi.org/10.2166/wh.2020.029.