Association Between the Diabetic Foot Ulcer and the Bacterial Colony of the Skin Based on 16S rRNA Gene Sequencing: An Observational Study

Introduction

Diabetes is a growing public health problem with increasing prevalence and mortality that threatens the health of people globally.¹ In 2021, the global prevalence of diabetes in adults aged 20 to 79 years was estimated at 537 million, mainly caused by type 2 diabetes, and the number is expected to reach 783 million by 2045.² Among these, diabetic foot ulcer (DFU) is a common complication with a high prevalence among diabetic patients.³ The International Working Group on Diabetic Foot (IWGDF) defines DFU as: patients with initial diabetes or a history of diabetes have ulcers in the foot, usually with lower limb neuropathy and/or peripheral arterial disease (PAD).⁴ Some studies have shown that 19–34% of diabetic patients may be affected by DFU during their lifetime, while the International Diabetes Federation reports that 9.1–26.1 million people will develop DFU each year.⁵ Moreover, DFU has a high mortality rate. Data show that DFU is associated with 5% of deaths within the first 12 months of diabetes onset, and 42% of deaths within 5 years of diabetes onset.⁶ DFU treatments should follow a multidisciplinary approach, such as local epidemiology-based antibiograms, personalized treatment, regular debridement, periodic wound assessment, and dressing changes, which requires years of experience in treating the condition⁷. Unfortunately, the failure rate for DFU patients to heal within a year of treatment is up to 20%, and the recurrence rate is as high as 40%.⁸ It is obvious that the unidentified mechanisms underlying the genesis and maintenance of DFUs make their treatment difficult. Therefore, additional research into the mechanisms is necessary to enhance the effectiveness of DFU treatment.

DFU is a chronic skin disease. Previous studies have shown that the healing and recurrence of chronic skin wounds are associated with altered bacterial diversity and instability in the microbiota composition of the wound.⁹ Skin microbiota plays an ambivalent role in wound healing, altering the microbiota composition of wounds directly affects their healing potential. On the one hand, some bacterial species can hinder the wound healing process and lead to chronic wounds; on the other hand, some bacterial species can accelerate the wound healing process and suppress the colonization of pathogens.^10,11 There is a significant difference in skin microbial composition between diabetic and non-diabetic patients. Compared to healthy mice, a significant increase in Staphylococcus, Aerococcus, Porphyromonadaceae and Proteobacteria as well as an obvious decline in Streptococcus, Lachnospiraceae, and Acinetobacter were observed in the wounds of diabetic mice.¹² There also exists a number of reports on the microbiota of DFU. For example, Gardner et al found that Staphylococcus, especially S. aureus, was the most common organism in the DFU skin. They also reported that the depth of diabetic ulcer was positively correlated with the amount of anaerobe and negatively correlated with the amount of Staphylococcus, and that ulcer duration was positively correlated with the species diversity, richness and relative abundance of Proteobacteria.¹³ Besides, it has been shown that DFU patients exhibit more pairs of pathogenic microorganisms, such as Candida and Pseudomonas.¹⁴ It is evident that alterations in the microbiota may be a major factor in the development of DFU. However, the wide variety of microorganisms and the multiple microbial properties usually present in DFU result in potentially diverse results of bacterial identification. Therefore, additional information is required to clarify the association between DFU and skin bacterial colony. The aim of this study was to look at the bacterial colonies in healthy skin, diabetic skin and DFU skin and to screen for organisms that mediate key changes in the bacterial colonies using 16S rRNA gene sequencing. It is anticipated that this study will shed fresh light on the relationship between DFU and microorganisms.

Materials and Methods

Study Subjects

A total of 48 patients with diabetes aged 18 years or older were recruited at Panyu Central Hospital from March 2021 to March 2022, including 22 patients with DFU (T group), 26 diabetic patients without foot ulcers (TW group), and 10 healthy volunteers (H group). Inclusion criteria: (1) patients were diagnosed with diabetes according to the diagnostic criteria of American Diabetes Association;¹⁵ (2) patients with foot ulcers had an explicit ulcer on the foot or a recurrent foot ulcer for at least 6 months. Exclusion criteria: (1) patients received antibiotic therapy one week prior to sample collection; (2) patients were combined with active wound infection, untreated osteomyelitis, gangrene, immunocompromised disease; (3) wounds were unrelated to diabetes or its complications; and (4) patients with tumors.¹⁶ The flow chart of the study process was showed in Figure 1. This study was approved by the ethics committee of the Ethics Committee of Panyu Central Hospital (ethics number: PYRC-2020-008).

Figure 1 A flow chart of study process.

DNA Extraction

The swab was removed and moistened the tip with sterile saline. The skin of the foot ulcer area in T group (including toe, heel, sole, medial, lateral and dorsal side of the foot) was swabbed. At the same time, the skin of the foot was swabbed in the same position in TW group and H group. Then, the swabs with samples were quickly placed into the sampling tube, snapped the break against the wall of the mouth of the 2 mL collection tube, and left the swab head inside the liquid of the sampling tube. At least two swab samples were collected from each patient. The collected samples were immediately frozen in collection tubes at −2 °C~-5 °C. DNA was subsequently extracted from the skin swabs using a Bacterial Genomic DNA Isolation Kit (DP302, TIANGEN) and sequenced using a 16S rDNA amplicon.

16S rRNA Gene Sequencing

The obtained bacterial DNA was used as a template to amplify the V3-V4 region of the bacterial 16S rRNA gene using PCR. The sequences of the primers used were 343F: 5’- TACGGRAGGCAGCAG-3’ and 798R: 5’- AGGGTATCTAATCCT −3’. The obtained PCR products were subjected to 1% agarose gel electrophoresis (120 V, 30 min). The concentration of the collected DNA was detected using Qubit^® 2.0 (Qubit Flex, Invitrogen). The prepared libraries were then sequenced using tags on the MiSeq platform (Illumina Inc., USA), and 250 bp paired-end reads were obtained. Approximately 10 million reads were generated. These sequences were joined into the original longer tags by overlapping their 3’ ends. Then, the chimeric sequences were identified and removed using UCHIME v4.2 software to obtain valid tags. The tags with 97% sequence similarity were clustered in QIIME software (version 1.8.0)¹⁷ using UCLUST¹⁸ to obtain operational taxonomic units (OUTs). Based on Silva (bacteria, http://www.arb-silva.de) and UNITE (fungi, http://unite.ut.ee/index.php; to OTUs for taxonomic annotation) bioclassification databases, the OUTs were classified and annotated. A 0.005% sequence of all sequence numbers was used as a threshold to filter OUTs.¹⁹ The taxonomic information of the species corresponding to each OTU was obtained by comparing the representative sequences of OTUs with the microbial reference database. Later, the composition of each sample was counted at each level (ie, phylum, class, order, family, genus, species). Based on this, species abundance tables at different taxonomic levels were generated using the QIIME software, and the colony structure charts at each taxonomic level using the R package. Species annotation was performed using the RDP classifier (version 2.2, http://sourceforge.net/projects/rdpclassifier/) with a believable threshold of 0.8.²⁰

Data Analysis

SPSS 23.0 software was used for statistical analysis. Measurement data were expressed as mean ± standard deviation (SD). Analysis of variance or Kruskal–Wallis test was used for comparison between multiple groups. The least square difference test was performed to further analyze the differences between the two groups. Independent samples t-test was used for comparison between two groups. The differences were statistically significant when P < 0.05. The α diversity of the samples (including Chao1, Shannon, Simpson, Richness) was analyzed using Mothur software (version v.1.30). The Bray–Curtis and Euclidean algorithms of QIIME and R software were adopted for β-diversity analysis. In addition, statistical analysis of the variation of microbial colonies in different groups of samples was performed, and linear discriminant analysis (LDA) of multilevel species difference was carried out by using linear discriminant analysis effect size (LEfSe).

Results

Composition of the Bacterial Colony

First, we analyzed the bacterial colony structure by OTUs based on the 16s sequencing results of the skin bacterial colony in three groups of subjects. The results showed that there were a total of 3453 amplified sequences in the TW group, 1475 amplified sequences in the T group, and 4173 amplified sequences in the H group. A total of 2812 sequences were identical between the TW group and the H group, and 981 identical sequences were identical between the TW group and the T group, with 944 amplicons common to all three groups (Figure 2A). The bacterial flora of each subject was further analyzed. As shown in Figure 2B–D, the flora structure formed by Staphylococcus of the genus Staphylococcaceae in the phylum Firmicutes was predominant in the H group, and that formed by Enhydrobacter of the genus Moraxellaceae in the phylum Proteobacteria and by Corynebacterium_1 of the genus Corynebacteriaceae in the phylum Actinobacteria were also present in the H group. In the T group, a significant decrease in the abundance was observed in Staphylococcus of the genus Staphylococcaceae in the phylum Firmicutes, while an increase in the abundance was found in Proteus of the genus Enterobacteriaceae and the phylum Proteobacteria, Corynebacterium_1 of the genus Corynebacteriaceae in the phylum Actinobacteria, and bacterium_NLAE-zl-H13 of the genus Bacteroides in the phylum Bacteroidetes. It can be noticed that when Enterococcus and Streptococcus were the dominant organisms in the group, the abundance of the others decreased significantly. In contrast to the H group, the abundance of Corynebacterium_1 was obviously increased in the TW group. More interestingly, the abundance of Staphylococcus was significantly raised in the TW group compared to the other two groups (Figure 3). Overall, the differences are very obvious in the bacterial colony structure between the TW group and the T and H groups.

Figure 2 Composition of bacterial colonies. (A) Venn diagram of sequencing results of the skin bacterial colony in three groups of subjects by OTU analysis; (B–D) assessment of bacterial colony structure in three groups of samples at phylum, family, and genus levels. Abbreviation: OUT, operational taxonomic unit.

Figure 3 Bacterial colony composition. Bacterial colony structure of the top 15 most abundant different OTUs in T, TW, and H groups. Abbreviation: OUT, operational taxonomic unit.

α Diversity Analysis

The α diversity analysis was subsequently performed on the bacterial colony of the samples from each group. It was observed that the chao1 index was much lower in the T group relative to the H and TW groups, but it was much higher in the H group than in the TW group (Figure 4A). As for the Shannon index, the TW group showed significantly lower than the H group, and the T group showed higher than the H and TW groups (Figure 4B). A notable rise in the Simpson index was observed in the T group relative to the H group; however, there was no significant difference in the Simpson index between the H and TW groups and between the T and TW groups (Figure 4C). In addition, the richness index was significantly decreased in the TW group compared to the H group, and compared to the H and TW groups, the richness index was considerably lower in the T group (Figure 4D). Collectively, the bacterial colony in the skin of patients in the T group had significantly lower bacterial counts, bacterial diversity, species richness and evenness than those in the H group and the TW group.

Figure 4 Variation of α-diversity index among the three groups. (A) α diversity analysis of chao1. (B) α diversity analysis of Shannon. (C) α diversity analysis of Simpson. (D) α diversity analysis of richness.

β Diversity Analysis

Meanwhile, we performed β diversity analysis on the sample data. No matter the results of Bray–Curtis-based PCoA analysis or those of Euclidean-based PCoA analysis, it was discovered that the T group could be clearly distinguished from the H group with statistical significance, as well as clearly distinguished from the TW group. Further observation showed that the microorganisms in the TW group were partially dispersed in the T group and the H group (Figure 5A and B). Thus, it can be seen that the bacterial colony structure of the T group is different from that of the H and TW groups.

Figure 5 β diversity of microbiota in the three groups. (A) PCoA analysis based on Bray–Curtis. (B) PCoA analysis based on Euclidean.

LEfSe Analysis

In order to identify the key microbiota in the skin bacterial colony, we used LEfSe to analyze the species composition of microbiota in each group. The results showed that the relative abundance of genus Chryseobacterium and genus Lactobacillaceae Roseomonas was higher in the H group than in the TW and T groups. In the TW group, the relative abundance of genus Deinococcaceae and genus Staphylococcaceae was the highest. The highest relative abundance was found in Prevotella, Vagococcus, Peptostreptococcus, and Leptotrichiaceae in the T group (Figure 6A). Later, linear discriminant analysis (LDA) was then used to reveal microbiota that differed significantly among groups. In the H group, the changes of Actinobacteria, Micrococcales, Corynebacteriales, Corynebacteriaceae, and Corynebacterium were observably larger than those in the other two groups. The TW group displayed more distinct changes of Staphylococcaceae, Staphylococcus, Bacillales, Bacilli and Firmicutes than the other two groups. In terms of the T group, the changes of Bacteroidales, Bacteroidia, Enterobacteriaceae, Enterobacteriales and Prevotella were significantly more variable than the other two groups (Figure 6B).

Figure 6 LEfSe analysis. (A) composition of species that differ in bacterial colonies in T, TW, and H groups revealed by LEfSe analysis. (B) linear discriminant analysis of microbiota that differ significantly between groups.

Relative Abundance Percentages of Different Operational Taxonomic Units with the Most Significant Abundance in the Three Groups

Further observation of the significantly different OUTs in the three groups was conducted (Figure 7). Among them, OUT_1 is the Staphylococcus in the phylum Firmicutes, family Staphylococcaceae, with relatively low abundance in the T group relative to the TW and H groups. Similarly, OUT_2 is the Corynebacterium_1 in the phylum Actinobacteria, family Corynebacteriaceae, with relatively low abundance in the T group relative to the TW and H groups. OUT_10 is the Enhydrobacter in the phylum Proteobacteria, family Moraxellaceae, with lower abundance in the T group than in the TW and H groups. OUT_6 indicated the Enterobacteriaceae in the phylum Proteobacteria, family Enterobacteriaceae, which showed relatively low abundance in the T group relative to the TW and H groups. OUT_3, which represented Pseudomonas in the phylum Proteobacteria, family Pseudomonadaceae, had higher abundance in the T group than the TW and H groups.

Figure 7 Relative abundance percentages of the different operational taxonomic units with the most significant abundance in the T, TW, and H groups.

Discussion

Bacterial infection is the most common cause of delayed healing in a variety of chronic skin wounds. DFU is not an exception.²¹ However, the lack of appropriate diagnostic tools makes it impossible to distinguish whether the bacteria in the ulcer are a change in the status of the original colonized bacteria or originate from an external infection. Few studies have been performed on the change in flora from healthy skin to diabetic skin to DFU skin. Therefore, we revealed this change by 16S rRNA gene sequencing and found that the microbiota composition of DFU skin is entirely different from that of healthy and diabetic skin.

It is important to emphasize that the 16S rRNA gene varies between bacterial species. 16S rRNA gene sequencing was performed on the selected hypervariable regions within the 16S gene, although shorter sequences can provide results comparable to the whole genome. Sadly, no single hypervariable region has yet been identified that can discriminate between all known bacterial species.²² Right now, 16S rRNA gene sequencing is still a relatively more comprehensive detection method. In this study, participants were divided into healthy controls (H group), diabetic group (TW group) and DFU (T group). OTU analysis was performed after 16S sequencing using QIIME software. QIIME software focuses on sequencing and functional analysis of microbial colony diversity including data preprocessing, OUT analysis and sample difference analysis.¹⁷ As a result, 3453 amplified sequences were discovered in the TW group, 1475 amplified sequences in the T group and 4173 amplified sequences in the H group. Further analysis showed that there was a significant change in both the phylum level, family level and genus level in the T and TW groups compared to the H group. In a previous meta-analysis of DFU isolated strain, the most common pathogen was Staphylococcus aureus, Methicillin-resistant Staphylococcus aureus accounted 18.0%, and Pseudomonas spp., Escherichia coli and Enterococcus spp. were also included in it.²³ Jnana et al compared microorganisms in wounds of diabetic and non-diabetic patients. They reported that in the wounds of diabetic patients, the abundance of Alcaligenes, Pseudomonas, Burkholderia and Corynebacterium changed significantly. On the other hand, Maria et al isolated relatively rare bacteria such as Abiotrophia defectiva, Providencia rettgeri, and Providencia stuartii in DFU patients.²⁴ It is evident that there is an extremely high heterogeneity of the microbiome in DFU wounds. In the results of the present study, it can be observed that Staphylococcus was dominant in H and TW groups. Notably, the abundance of the Staphylococcus was significantly reduced, while that of the Enterobacteriaceae, Corynebacterium_1 and bacterium_NLAE -zl-H13 increased in the T group. Such a change in the microbiota was similar to the findings of previous studies.²⁵

There are numerous indicators of α diversity that can be expressed in Mothur software with Shannon, Berger-Parker, Simpson, Q statistic; observed richness, Chao1, ACE, and jackknife.²⁶ The β diversity analysis mainly reflects the degree of variation between samples. Several distance metrics, such as Unifrac, Bray–Curtis, Euclidean, Jaccard index, Yue & Clayton, and Morisita-Horn, have been often employed using QIIME and R software.²⁷ In this study, the diversity of skin microorganisms was found to have gradually decreased from the H group to the T group in accordance with the results of α diversity analysis. Besides, β diversity analysis results revealed that the species diversity composition of the T group was significantly different from that of the TW and H groups, as the TW and H groups were in a state of gathering with small differences in diversity. These results differ from the findings of most of the previous reports. Some studies have suggested that species richness of microbiota shows a negative correlation with wound healing and a positive correlation with the increased time of wound healing.^7,28 Other relevant studies have indicated that the burn wound healing time of Pseudomonas aeruginosa-colonized mice is significantly prolonged. The longer healing time can be explained by that the biofilms of Pseudomonas aeruginosa disrupt host defenses, thereby inducing a steady state of chronic infection that impairs wound healing.²⁹ Similarly, it has also been shown that Pseudomonas aeruginosa significantly delays wound healing by inhibiting the formation of new blood vessels in wounds.³⁰ What’s more, a study published in 2021 concluded that persistent infection with Staphylococcus aureus in porcine wound biofilm resulted in slowing down wound healing by causing deficiencies in granulation tissue collagen.³¹ It is thus clear that the slow wound healing caused by microorganisms may not only be a result of the flora abundance but may also be mediated by the key flora. In the present study, it could be observed that the change in abundance of Pseudomonas spp. was significantly increased in the T group, which could also be responsible for the persistence of foot ulcers.

In order to identify the key flora in the change of skin flora in the presence of healthy skin to DFU skin, we analyzed the relative abundance percentages of the different OTUs with the most significant abundance in the three groups. According to the analysis results, the abundance of Staphylococcus, Enhydrobacter and Corynebacterium _1 gradually decreased as the microbial colony structure changed from the H group to the T group, while Enterobacteriaceae and Pseudomonas showed an up-regulation trend. However, previous studies revealed that Staphylococcus and Pseudomonas are very common flora in diabetes, which increases the probability of ulcer infection. Antibiotics are also the primary line of defense against staphylococci due to the small number of Escherichia coli.³² As for the findings of this study, Staphylococcus in ulcers also changed significantly in abundance, and the number of Escherichia coli increased considerably, which may be due to the change in the dominant microflora.³³ Additionally, Kurup et al reported that Gram-negative bacteria (63.0%) were more prevalent than Gram-positive bacteria (37.0%) in the wounds of DFU patients; Pseudomonas (18.8%) was the most common isolate, followed by Escherichia coli (13.9%).³⁴ Such a phenomenon again emphasizes the heterogeneity and complexity of the bacterial histology of DFU.

Conclusion

In summary, the bacterial colony structure in the DFU skin changed significantly relative to the skin of healthy or diabetic patients. Notably, in the DFU skin, the abundance of Staphylococcus, Enhydrobacter and Corynebacterium_1 was significantly reduced, while that of Escherichia coli and Pseudomonas was remarkably increased. The difference between a healthy skin and diabetic skin with or without ulcer is the balance between normal and pathogenic microbiota. Altering the microbiota composition of wounds may help the treatment of DFU. 16S rRNA gene sequencing is a practical detection method with kinds of computational methods to achieve the variation of microbial colony analysis. However, the bacterial colony structure of DFU skin has obvious individual difference, and further studies are needed to break through the limitation brought on by individual difference.