Characterization of Staphylococcus and Corynebacterium Clusters in the Human Axillary Region

The skin microbial community is regarded as essential for human health and well-being, but likewise plays an important role in the formation of body odor in, for instance, the axillae. Few molecular-based research was done on the axillary microbiome. This study typified the axillary microbiome of a group of 53 healthy subjects. A profound view was obtained of the interpersonal, intrapersonal and temporal diversity of the human axillary microbiota. Denaturing gradient gel electrophoresis (DGGE) and next generation sequencing on 16S rRNA gene region were combined and used as extent to each other. Two important clusters were characterized, where Staphylococcus and Corynebacterium species were the abundant species. Females predominantly clustered within the Staphylococcus cluster (87%, n = 17), whereas males clustered more in the Corynebacterium cluster (39%, n = 36). The axillary microbiota was unique to each individual. Left-right asymmetry occurred in about half of the human population. For the first time, an elaborate study was performed on the dynamics of the axillary microbiome. A relatively stable axillary microbiome was noticed, although a few subjects evolved towards another stable community. The deodorant usage had a proportional linear influence on the species diversity of the axillary microbiome.

background was subtracted, differences in the intensity of the lanes were compensated during normalization and bands and band classes were detected. Clustering was done with Pearson correlation and the unweighted pair group with mathematical averages (UPGMA) dendrogram method. Relevant and non-relevant clusters were separated by the statistical Cluster Cutoff method (BioNumerics Manual 5.10). Similarities and abundances were extracted from the software and statistical analysis was performed using SPSS version 19. Significant cut-off values were indicated in the paper. Triplicate specimens of each sample were loaded on DGGE-gel, pooled and checked. No differences were seen.

Denaturing Gradient Gel Electrophoresis (DGGE) diversity indices
The range-weighted richness (Rr) is a DGGE specific range of values which indicate the richness and genetic diversity of species within the bacterial community. It is correlated with the distribution of the bands in the DGGE pattern and the percentage denaturant gradient of the gel needed to represent the sample's total diversity (within the limits of the technique). This is mathematically expressed as Rr=N 2 ×D g , where N represents the total number of bands in the pattern, and D g the denaturant gradient comprised between the first and the last band of the pattern [8]. The community organization (Co) describes the species abundance distribution in the microbial community and is calculated as the Gini coefficient times 100 [9]. The Gini coefficient (ranging from zero to one) is a single value that describes a specific degree of evenness measuring the normalized area between a given Pareto-Lorenz (PL) curve and the perfect evenness line. The higher the Gini coefficient, the more uneven a community is. The PL evenness distribution curve was constructed based on the DGGE profiles as previously described [8,10]. The community dynamics (Dy) was studied computing the moving window analysis 5 (MWA) plot of consecutive DGGE profiles of the same subject. The microbial community rate of change was conducted using the UPGMA and distance matrices of each DGGE based on the Pearson correlation similarity coefficient to cluster the succeeding samples [8].

Sanger sequencing
The 16S rRNA genes of five isolated pure strains, plated on blood agar plates, was amplified by PCR using the forward primer P63F and the reverse primer P1387R [11]. The PCR program was performed and checked as described in Table S1B. Sanger sequencing was performed on the 16S rRNA amplicons and aligned and compared with sequences from the NCBI database. The closest match of each isolate was identified. Afterwards, sequences of all strains were submitted to GenBank and the submission numbers are presented in Table S2.

Pyrosequencing
Amplicon pyrosequencing was performed on the total DNA extracted from nine specific individual samples. Barcoded amplicons were prepared with the primers 530F-mod [12] and 1061R [13] amplifying a 562 bp DNA fragment flanking the V4, V5 and V6 regions of the 16S rRNA gene [14], extended as amplicon fusion primers with respective primer L adaptor, key sequence and multiplex identifiers (MID) on the forward primer. Amplicons were generated by using FastStart High Fidelity Taq DNA Polymerase kit (Roche) under the conditions mentioned in Table S1B. Amplicons were purified with the High Pure PCR Product Purification Kit (Roche) and pooled as specified by the manufacturer. The purity and quality of the PCR products was verified on agarose gel. Emulsion PCR, emulsion breaking and sequencing were performed applying the GS FLX Titanium chemistry protocols and using a 454 GS FLX pyrosequencer 6 (Roche) as recommended by the manufacturer. For this study, nine amplicons were sequenced in a pool of 17 mixed amplicons on 1/4th of an FLX picotitre plate. Quality filtering of the pyrosequencing reads was performed using the automatic amplicon pipeline of the GS Run Processor (Roche), with a modification of the valley filter (vfScanAll-Flows false instead of TiOnly) to extract sequences. The raw flowgrams were processed and analyzed in an in-house Mothur [15]

Alpha diversity analysis
The alpha diversity was calculated to characterize the diversity of one individual axillary sample. Figure S6 displays the range-weighted richness and the community organization of 43 individual DGGE samples. For the pyrosequenced samples, the analysis was performed on both the subsampled (with normalization at 7135 sequenses) and complete (without normalization) 7 dataset. The Shannon diversity index and the Chao 1 richness estimator were calculated and plotted in Figure S5 and S4, respectively. To assess the completeness of sampling also a rarefaction analysis was performed, and plotted in Figure 5. Detailed alpha diversity characteristics are displayed in Table S3. The Shannon diversity index, Chao1 richness estimator and observed richness of the axillary samples of the nine subjects were plotted in function of their weekly deodorant usage ( Figure 6).

Beta diversity analysis
The beta diversity analysis was calculated to study the difference between the individual microbial axillary communities. A heatmap was generated in Figure 4 of the top 25 OTU's of the subsampled dataset with hierarchical complete linkage clustering based on Bray-Curtis similarities. Next, a heatmap of the Yue and Clayton θ dissimilarity index of the pairwise comparison of the pyrosequenced samples is shown in Figure S1A. Figure S1B represents a heatmap of the similarities of the different samples analyzed by DGGE.

Hypothesis testing
Additional hypothesis tests were conducted in order to identify whether the two cluster are dissimilar or not. The null-hypothesis (H 0 ) is that the community structures are similar. The parsimony method (aka P-test) is a generic test that describes whether two or more communities have the same structure. The results of the parsimony testing procedure are displayed in Table   S4A. A cluster dendrogram is figured in Figure S3A clustered by UPGMA with the Bray-Curtis index for community structure. Nonmetric multidimensional scaling (NMDS) is efficient at identifying underlying gradients and representing relationships based on various types of 8 distance measures. The NMDS plot is displayed in Figure S3B of the pyrosequenced samples, based on the abundance-based Jaccard distance measure. Using molecular variance (amova) it is tested whether the centers of the clouds of the axillary clusters are separated or not. In Table   S4B, it is tested whether the observed cluster separation between the Corynebacterium and Staphylococcus clusters in the NMDS plots is statistically significant.

DGGE versus Sanger sequenced isolates versus 454 pyrosequencing
In this research, DGGE fingerprinting, Sanger sequencing of isolates and 454 pyrosequencing were successfully combined. DGGE is a relatively rapid method, well suited for mixed microbial communities and is an interesting technique to study the community dynamics and diversity.
Firmicutes and Actinobacteria are important phyla on the skin, with very dissimilar GC content.
Firmicutes are known to possess a low GC content, whereas Actinobacteria are known to have a high GC content. As such, DGGE is a very solid technique to differentiate amongst the two. The disadvantages of DGGE were resolved by combining with 454 pyrosequencing and sequencing of isolated bacteria. Although isolation is only possible for a subset of the skin bacteria, some bands on DGGE were identified by combining the DGGE pattern of the pure culture with the pattern of the mixed axillary culture. 454 pyrosequencing made identification on genus level and quantification possible, but was a less flexible and a rather slow method, due to the rigorous analytic and statistical work. Figure 1  Callewaert et al., Figure S1. Heatmap to assess the interpersonal diversity. (A) Heatmap of the Yue and Clayton θ dissimilarity ሺ1 െ ‫ܦ‬ ఏ ೊ ሻ index of the pairwise comparison of 9 pyrosequenced sampled community structures to assess the interpersonal diversity. The higher the index (red), the more similar the communities; (B) Heatmap of the DGGE similarities of 43 sampled community structures. The higher the index (red), the more similar the communities. The nine pyrosequenced samples are indicated with MID (multiplex identifiers).

A B
Callewaert et al., Figure S2. Stacked bar sample-wise taxonomic distribution of the sequences (A) on the order level; (B) on the class level; (C) on the phylum level.

A B
Callewaert et al., Figure S2. Stacked bar sample-wise taxonomic distribution of the sequences (A) on the order level; (B) on the class level; (C) on the phylum level.
13 C Callewaert et al., Figure S3. (A) Cluster dendrogram by UPGMA with the Bray-Curtis index for community structure. (B) NMDS plot based on the abundance-based Jaccard distance measure. Red sample names are from the Corynebacterium cluster, blue sample names Staphylococcus cluster.

A B
Callewaert et al., Figure S4. Chao1 richness estimator of the data on the complete dataset (A) and on the normalized dataset (B).

A B
Callewaert et al., Figure S5. Shannon community diversity index datacurve on the complete dataset (A) and on the normalized dataset (B). The normalized and complete dataset displayed a stable curve for the index, indicating a reliable diversity index.

A B
Callewaert et al., Figure S6. Diversity indices of the DGGE results: Range-weighted richness (Rr) -indication of species richness -and community organization (Co) -indication of species evenness -of 43 DGGE samples.
Callewaert et al., Figure S7. Information and permission form for subjects giving axillary samples. 18

Axillary samples: Information and permission form
We will take 6 samples of your axillary bacteria by means of cotton swabs (3 left and 3 right). Your data will be handled confidentially. Research approved by the Ghent University Hospital Ethical Committee.

Person code Gender
Date & time Age Here used for the pyrosequencing results.