A comprehensive reference catalog of human skin DNA virome reveals novel viral diversity and microenvironmental influences
Zhiming Li, Shenghui Li, Chongyin Han, Yuxiao Chen, Hefu Zhen, Yuzhe Sun, Xiaofeng Zhou, Yanmei Chen, Yan Zheng, Lianyi Han, Jean Krutmann, Chao Nie, Jiucun Wang, Jingjing Xia
TL;DR
This study creates a detailed catalog of viruses on human skin, revealing new viral diversity and how viruses vary in different skin areas.
Contribution
The study provides a comprehensive catalog of previously unknown skin viruses and their associations with skin microenvironments.
Findings
20,927 viral sequences were identified, with 90.85% being previously unknown.
Distinct viral communities were found in oily, dry, and moist skin areas.
Bacteriophages were linked to bacteria like Pseudomonas and Staphylococcus.
Abstract
Human skin serves as a dynamic habitat for a diverse microbiome, including a complex array of viruses whose diversity and roles are not fully understood. A total of 2,760 skin metagenomes from 6 published skin studies were collected. A skin virome catalog was constructed using standard methods in the viromics field. Viral characteristics were identified through cross-cohort meta-analysis and used to characterize viral features across different skin environments. We identified 20,927 viral sequences, which clustered into 2,873 viral operational taxonomic units (vOTUs), uncovering a substantial breadth of viral diversity on human skin. The results also highlight significant differences in viral communities that are associated with varying skin microenvironments. The oily skin is enriched in Papillomaviridae, the dry skin area is enriched in Autographiviridae and Inoviridae, and the moist…
Genes, proteins, chemicals, diseases, species, mutations and cell lines named across the full text — each resolved to its canonical identifier and authoritative record.
Click any figure to enlarge with its caption.
Fig 1
Fig 2
Fig 3
Fig 4Peer Reviews
No public reviews on file for this paper yet. If you reviewed it on a platform where reviews are public (OpenReview, ICLR, NeurIPS, ICML), you can paste yours below so the community can read it here.
Videos
No videos yet. Explain this paper in a talk, walkthrough, or lecture? Add one.
Taxonomy
TopicsPlant Virus Research Studies · Bacteriophages and microbial interactions · Plant and Fungal Interactions Research
INTRODUCTION
The skin, the largest organ of the human body, not only acts as a protective barrier against chemical injury, ultraviolet radiation, and pathogens but also harbors a complex ecosystem of microorganisms including bacteria, fungi, archaea, and viruses (1, 2). This microbial community, known collectively as the skin microbiome, plays a crucial role in modulating immune responses, protecting against pathogen colonization and maintaining skin health (3, 4). While considerable research has focused on the bacterial or fungal component of the skin microbiome, the diversity and function of the viral component remains less explored. Existing research indicates that viruses account for approximately 1.51% of the microbial population on facial skin samples (5). There is also an observed phenomenon known as “viral bloom” in specific samples, where viruses can comprise up to 96% of the microbial community at certain skin sites (2). The most prevalent viruses on the skin surface are bacteriophages, including those that target Cutibacterium and Staphylococcus. Additionally, potential human viral pathogens such as human papillomavirus and Merkel cell polyomavirus are also present (2).
Metagenomic sequencing and bioinformatics have provided tools to delve deeper into the virome’s complexity, without the limitations of traditional culture-based techniques. Viruses, particularly bacteriophages, which are viruses that infect bacteria, are now recognized as pivotal players in shaping microbial community dynamics. They can influence bacterial diversity and density through lytic cycles—where viruses replicate within and lyse their bacterial hosts—and lysogenic cycles—where viral genomes integrate into the bacterial genome, potentially altering host bacterial behavior and pathogenicity (6–8).
The exploration of the skin virome is crucial for several reasons. First, viruses can influence the health and disease states of the skin by impacting microbial balance and host immune responses (9–11). For example, bacteriophages can control the population of pathogenic bacteria, thereby preventing or exacerbating skin conditions such as acne, atopic dermatitis, and wound infections (12, 13). Second, understanding the skin virome might open up therapeutic avenues, such as phage therapy, to combat antibiotic-resistant bacterial strains (13). Finally, the skin virome can serve as a diagnostic marker for certain skin diseases and potentially for systemic disorders, given the skin’s interaction with other organs (14–16).
To enhance our understanding of the skin virome, we have developed a non-redundant viral reference catalog based on human skin samples collected worldwide. This project involved the analysis of 2,760 metagenomic data sets from the United States, China, Singapore, and Italy, covering a variety of microenvironmental conditions. Our goal in establishing this catalog was to map the diversity of the skin virome, investigating its interaction with host and environmental factors.
MATERIALS AND METHODS
Data collection
We conducted a comprehensive review of studies focused on skin microbiomes, excluding those that relied solely on ITS or 16S sequencing due to their limitations in capturing viral genomic data, which is essential for our cataloging objectives. Our inclusion criteria were limited to studies utilizing metagenomic sequencing techniques without enrichment for virus-like particles. We did not filter studies based on variations in sequencing platforms, extraction methods, or library preparation techniques. We collected and downloaded a total of 2,760 skin metagenomic data sets from the NCBI Sequence Read Archive (SRA) databases (https://www.ncbi.nlm.nih.gov/sra) and CNGB (https://db.cngb.org/). To provide a clear demographic context for our data, we summarized the age, gender, country, and specific skin sites of the subjects in Table S1, with additional metadata detailed in Table S2. In the collected data sets, 97.72% (2,697) are generated from healthy individuals, 1.41% (39) are from persons with atopic dermatitis, and 0.87% (24) are from persons with psoriasis.
Data processing
The raw reads were filtered using SOAPnuke (v2.1.9, with default parameters) (17, 18). Human reads were further removed by aligning the filtered reads with the human genome hg38 using bowtie2 (v2.3.5.1, with default parameters) (19). The remaining reads of each data set were assembled de novo using MEGAHIT (v1.2.9) (20) based on various k-mer sizes (k = 21, 33, 55, 77, 99).
Bacterial profiles
The construction of the bacterial profile followed previously described methods (21, 22). High-quality reads were mapped to 3,547 putative prokaryotic species (23) using Bowtie2. Subsequently, the “jgi_summarize_bam_contig_depth” tool (24) with default parameters was used to calculate the coverage depth of contigs for these 3,547 prokaryotic species. To account for differences in sequencing depth among samples, we estimated the relative abundance of each contig using normalized coverage depth. Finally, the median contig relative abundance within each species was used as the species relative abundance.
Virus identification and decontamination
To maximize the recovery of viral diversity, especially small-genome viruses that are common in human-associated viromes, we included contigs larger than 2,000 nt for downstream analysis. This threshold reflects the biological reality that many viruses, such as single-stranded DNA viruses, have genomes below 5 kb. Similar thresholds have been used in recent large-scale studies (25, 26). Contigs larger than 2,000 nt were preliminarily screened through VIBRANT (v1.2.1) (27) and DeepVirFinder (v1.0) (28) (P-value < 0.01 and score > 0.90), identifying 2,480,992 contigs as potential viral sequences. To purify these viral sequences, following previous studies (29), we utilized hmmsearch (v3.1) (30) with default options to search for bacterial, archaeal, and fungal universal single-copy orthologs (BUSCO-v5.4.4) (31) within the viral sequences and removed them. Subsequently, we used BLASTN to align the resulting viral sequences against the human genome to remove sequences potentially originating from human genomic regions. Finally, we used geNomad (32) (v1.8.1, with default parameters) to further filter the sequences (virus score of at least 0.7) and employed CheckV (33) (v1.0.3, with default parameters) to assess the quality and completeness of single-contig virus genomes. Through these steps, we identified 4,199 complete viral genomes, 5,249 high-quality viral genomes with >90% completeness, and 11,479 medium-quality viral genomes with 50%–90% completeness. The contigs were clustered through the following steps: (i) Use BLASTN alignment to calculate the similarity and coverage of contigs. (ii) The HipMCL algorithm (34) was used to perform sequence clustering based on 95% similarity and 85% alignment coverage, resulting in 2,873 potential viral operational taxonomic units (vOTUs), yielding 558 complete viruses, 841 high-quality viruses, and 1,474 medium-quality viruses.
Viral taxonomy
We annotated vOTUs based on the annotation results from geNomad (32) (v1.8.1, with default parameters. This classification method follows the taxonomy in ICTV’s VMR number 19. In the SVD catalog, 98.75% (2,837 out of 2,873) of the viral operational taxonomic units (vOTUs) were classified.
Host prediction
Host prediction for vOTUs was carried out using iPHoP (35) (version 1.3.3, default parameters). The database was constructed using a published collection of skin microbiome genomes (36). The construction process strictly followed the procedure outlined in iPHoP under “Adding bacterial and/or archaeal MAGs to the host database.
Functional annotation
In the Skin Virome Database (SVD), protein-coding sequences of viral operational taxonomic units (vOTUs) were predicted and annotated using Prodigal and EggNOG (v 5.0) (37) and VOGDB (38).
Phylogenetic analysis
Phylogenetic analysis was conducted on 2,873 vOTUs from the Skin Virome Database (SVD) using a phylogenetic approach based on amino acid sequence similarity. We utilized ViPTreeGen (v1.1.3) (39) to generate a viral proteomic tree of the SVD. The proteomic tree was then visualized using iTOL (v6.9) (40).
Quantification of vOTUs
To explore the biogeographical characteristics of the skin virome, we analyzed metagenomic data collected from oily, moist, and dry skin areas, aligning it against all vOTUs in the Skin Virome Database (SVD), covering a total of 2,334 metagenomic data sets. We used bowtie2 with a >95% identity threshold to align high-quality reads with all vOTUs in the SVD. Subsequently, we conducted a sequence-based vOTU abundance analysis using jgi_summarize_bam_contig_depth (default parameters) (24). Furthermore, a vOTU was considered present in a sample only when at least 75% of its contiguous sequence length was covered by mapped reads (coverage breadth ≥75%). The mapping of reads to the SVD generated counts, forming a coverage depth/abundance matrix. Considering that different samples might have varying sequencing depths, we used a normalized coverage depth matrix to estimate the abundance of vOTUs. This approach allows us to fairly compare the relative abundance of vOTUs across different samples, thereby revealing the variations in viral composition across different skin types.
Phage-host association analysis
To explore the associations between phage vOTUs and hosts within the community, we assessed the Spearman correlation between the relative abundances of phage vOTUs and microbial species. If a phage vOTU is predicted to infect a single genus, we calculate the correlation between that phage vOTU and its predicted host. If a phage vOTU is predicted to infect multiple genera, we calculate the correlation for each predicted host. If a phage-host pair (with zero abundance) does not exist in the sample, that sample is excluded from the correlation analysis.
Alpha diversity
Alpha diversity is estimated based on the relative abundances of phage vOTUs. The Shannon diversity index is calculated by calling the diversity function in R (V4.3.1) using the vegan package (41), with the parameter set to index = shannon. Comparisons of alpha diversity across different microenvironments and sites are conducted using the Wilcoxon rank-sum test, and a P-value of less than 0.05 indicates a significant difference.
Statistical analysis
Using the ade4 package (42) within the R platform, PCoA and PCA are performed on samples from different microenvironments and sites. Permutational multivariate analysis of variance (PERMANOVA) is performed using the adonis function from the vegan package. When using PERMANOVA, the skin microenvironment and skin site characteristics are analyzed after adjusting for age and gender. Statistical significance is verified using the wilcox.test and kruskal.test functions. P-values are adjusted using the p.adjust function with the parameter method = fdr. Adjusted P-values less than 0.05 are considered to be statistically significant.
RESULTS
The construction of the human skin viral catalog
We collected a total of 2,760 publicly available metagenomic data from human skin, establishing a skin metagenomic data set that covers the United States, China, Singapore, and Italy (Table S1). This data set comprises samples from 21 different skin sites, categorized into 3 skin microenvironments: oily (1,849 samples), dry (244 samples), and moist (241 samples). After uniform processing (reads preprocessing and host removal, with an average host proportion of 42.86%), these samples represented 34.11 Tb of high-quality non-human metagenomic data and were used to generate a total of 250.56 million long contigs (Table S2). Using an integrated pipeline (see Materials and Methods), approximately 0.96% (n = 2,422,445) of the contigs were detected as high confidence viral sequences (Fig. S1). Completeness prediction results suggested that there were 4,199 complete viral genomes, 5,249 high-quality viral genomes, and 11,479 medium-quality viral genomes (Fig. S1; Table S3). The viral sequences were clustered according to MIUViG standards, with >95% nucleotide similarity and 85% alignment coverage (43), forming a skin virome database (SVD) consisting of 2,873 viral operational taxonomic units (vOTUs) (Fig. S1). This database included 558 complete viral genomes, 841 high-quality viral genomes, and 1,474 medium-quality viral genomes (Fig. 1a; Table S4). The length of the vOTU representative sequences ranged from 2,010 bp to 433,960 bp, with a median length of 28,261 bp (Table S4).
Overview of viral operational taxonomic units (vOTUs) of skin. (a) Proteomic tree of 2,873 high-quality vOTUs. This tree was generated using ViPTreeGen. The outer ring displays the host classification for each vOTU: the innermost circle represents the viral family-level classification assignments. (b) Comparison of vOTUs identified in this study with those in databases. The Venn diagram illustrates the overlap between vOTUs identified in this study and those in public databases. The clustering of phage genomes has a threshold of >95% identity and >85% length coverage. (c) The bar chart displays the classification of the top 14 skin vOTUs. The host classification groups of the vOTUs are shown in different colors. (d) The ubiquity of skin surface viruses in the samples is depicted through multiple graphical representations. The blue bar chart (corresponding to the left Y-axis) shows the number of viral vOTU that appear in a given number of samples (X-axis). The dashed lines indicate the quantity of vOTUs and their corresponding total relative abundance for sample sizes of 2,208 and 1,380. Pie charts display the proportion of vOTUs corresponding to each family when the sample size is 2,208 and 1,380. A linear curve (corresponding to the right Y-axis) quantifies the cumulative relative abundance of these vOTUs in the samples.
The vOTUs of the SVD were compared with several existing viral catalogs, including five human gut virome catalogs (29, 44–47), two human-associated comprehensive viral catalogs (Cenote Human Virome Database (CHVD) (48) and virushostdb (49), one oral virome database (OVD) (50), a Metagenome-Assembled Genome Inventory for Children (MAGIC) (26), and available viral genomes in the RefSeq database. Of these, using a 95% similarity and 85% coverage threshold, only 263 vOTUs in the SVD matched with viral sequences in other catalogs, whereas 2,610 vOTUs did not match any known genomes in the databases (Fig. 1b; Fig. S2).
The 2,873 vOTUs encoded 119,826 protein-coding genes collectively. Functional annotation was performed using both the eggnog databases (37) and VOGDB (38). A total of 53,127 genes were annotated by both databases, among which 44,651 genes were annotated by eggNOG and 37,922 genes were annotated by VOGDB (Fig. S3a). The principal coordinate analysis (PCoA) of the vOTUs functional profiles revealed differences in the clustering of vOTUs based on identified functions (Fig. S3b; PERMANOVA test, P < 0.01). Different viral families show variations in functions related to DNA repair, DNA replication, peptidases and inhibitors, and replication and repair processes (Fig. S4). For example, the Autographiviridae are functionally enriched in DNA repair and recombination proteins, while the Herelleviridae were functionally enriched in replication and repair, and peptidases and inhibitors (Fig. S4).
Taxonomic landscape and host range of skin viruses
The taxonomic annotation results from geNomad (32) were utilized in the SVD catalog, with 96.66% (2,777 out of 2,873) of the vOTUs being classified. Notably, the majority, 76.51% (2,198), were assigned to unidentified family within the class Caudoviricetes. The viral families identified within the SVD include Papillomaviridae (n = 397), Inoviridae (n = 65), and Rountreeviridae (n = 51) (Fig. 1a through c). Papillomaviruses are known member of the human skin virome and represent a group of ancient viruses that had largely co-evolved with their hosts, primarily colonizing the human skin surface and hair follicles (51).
We found that at least 6 vOTUs, which belonged to the class Caudoviricetes, appeared in at least 2,208 samples (representing 80%) and accounted for 14.4% (median 8%) of the total viral abundance (Fig. 1d). In at least 1,380 (50%) of the samples, 78 vOTUs were identified, representing 31.5% (median 27.6%) of the total abundance, primarily including 72 members of unidentified families within the order Caudoviricetes, 6 unclassified vOTUs (Fig. 1d). Those commensal viral species might play a conserved role in human skin.
We used iPHoP to predict skin bacteria that would be infected by vOTUs. Among 1,908 (66.41%) vOTUs were predicted to infect bacteria on the skin surface (Table S4). The predicted infected bacteria include members of the phylum Firmicutes (n = 933), followed by Proteobacteria (n = 606) and Actinobacteriota (n = 297) (Fig. S5ab; Table S4). At the family level, the bacteria predicted to be infected include Staphylococcaceae (n = 664), Moraxellaceae (n = 206), Mycobacteriaceae (n = 184), Enterobacteriaceae (n = 150), and Peptoniphilaceae (n = 100) (Fig. S5c). At the genus level, the predicted infected skin bacteria include Staphylococcus (n = 658), Corynebacterium (n = 154), Acinetobacter (n = 105), Moraxella (n = 81), and Streptococcus (n = 54) (Fig. S5d). In addition to these abundant taxa, the phages were also predicted to infect Enterococcus (n = 14), Klebsiella (n = 43), and Pseudomonas (n = 50) (Table S4).
Further analysis revealed that unclassified phages within Caudoviricetes likely infect major bacteria in the human skin microbiome, including members of the Propionibacteriaceae, Mycobacteriaceae, Staphylococcaceae, Peptoniphilaceae, and Moraxellaceae (Fig. 2a). Caudoviricetes phages belonging to specific families such as Autographiviridae, Herelleviridae, Rountreeviridae, Crassvirales, and Schitoviridae are connected to Pseudomonadaceae, Caulobacteraceae, Enterobacteriaceae, Staphylococcaceae, Mycobacteriaceae, and Streptococcaceae (Fig. 2a). Among the unclassified phages from the class Caudoviricetes, phages predicted to infect multiple bacteria are expected to infect common skin bacteria such as Staphylococcus (n = 613), Corynebacterium (n = 149), Acinetobacter (n = 81), and Moraxella (n = 80) (Fig. 2b). Within the Rountreeviridae family of Caudoviricetes, phages predicted to infect a single genus typically infect Staphylococcus (n = 42), Confluentibacter (n = 1) and Streptococcus (n = 1). (Fig. 2c). In the family of Inoviridae, phages predicted to infect a single genus are expected to infect Acinetobacter (n = 22), Neisseria (n = 7), Pseudomonas (n = 6), and Stenotrophomonas (n = 3) (Fig. 2c).
General overview of bacteriophages vOTUs on the skin surface. (a) A Sankey diagram depicting skin vOTU at the family level and their corresponding hosts determined by iPHoP matching. The width of the links represents the number of relationships, with wider links indicating a greater number. (b and c) The distribution of phages vOTUs within the Caudoviricetes (b), Rountreeviridae and Inoviridae (c), showing the specific hosts they infect.
Close interactions between the skin phages and bacteriome
To explore the interaction between phages and bacteria in the skin, phages and bacterial profiles collected from 2,760 samples were compared. Results showed that the α-diversity (Shannon diversity) and β-diversity (Bray-Curtis distance) were significantly positively correlated (rs = 0.607 and 0.645, respectively; Fig. 3a and b), suggesting a close phages-bacteriome interaction.
The close interactions between the skin virome and bacteriome. Comparison of α diversity (a) and average β diversity (b) between the human skin virome and bacteriome. Each circle in b represents the average Bray-Curtis distance compared with other samples. The regression line is shown in red. (c) Network of associations between vOTUs and their multidrug-resistant bacterial hosts. Triangles represent vOTUs, while squares represent bacteria. The calculation of correlation was performed using Spearman’s rank correlation, and the presence of a gray line indicates a significant correlation between the two variables (FDR < 0.05).
We then analyzed the one-to-one correlations between the relative abundances of phages and their predicted host at the genus level across 2,760 samples. We found a positive correlation with an average rs of 0.357 between them, consistent with previous findings in the gut (rs = 0.18) that phages and their host bacterial species coexist rather than exclude each other in human skin (52). More importantly, bacteria with the potential for multidrug resistance, such as Pseudomonas, Klebsiella, and Staphylococcus,1, 53 were present on the human skin. We found that numerous phages vOTUs were significantly correlated with them, particularly members of the Caudoviricetes including the Rountreeviridae and Herelleviridae (Fig. 3c). These bacteria, of high clinical concern (54), can cause diseases in skin, blood, lungs, gastrointestinal tract, and other parts of the body (55, 56).
The biogeographic characteristics of viruses on the human skin
Human skin is characterized by different microenvironments (oily, dry and moist) and it has been shown that specific skin phenotypes are associated with distinct bacterial diversity (57–60). In this study, three of the six datasets with samples collected from different microenvironments were analyzed. It was observed that in CNP0003934, the viral diversity (Shannon) in dry skin was significantly higher than in oily skin similar to bacterial diversity. However, the viral diversity in oily and moist skin exhibited inconsistent trends in two data set (the SRP002480 and SRP057859); moist skin had higher diversity than oily skin in SRP002480, while in SRP057859, oily skin exhibited higher diversity than moist skin (Fig. 4a). Additionally, variations in diversity between different sites within the same skin microenvironment were observed (Fig. S6; Table S5).
Biogeographic characteristics of the skin virome. (a) Differences in the Shannon index at the vOTU level across different skin microenvironment within each data set. (b). Principal component analysis (PCA) based on the skin vOTU level. Different colors represent different cohorts, and shapes represent different microenvironment. Effect size (R²) and statistical significance are obtained through PERMANOVA (adonis). (c and d) Effect size on the skin virome of different sites within the same microenvironment (c) and across various cohorts (d) and across different microenvironment. In the figure, * indicates P < 0.05, and ** indicates P < 0.01. (e) Composition at the family level of the virome in different skin microenvironment in each data set. (f) Composition at the family level of the virome for different sites within the same skin microenvironment in each data set. (g) Comparison of the relative abundance of viruses at the family level in different microenvironment within each data set. An absolute z-score above 2 is considered statistically significant. Fh, forehead; Ck, cheek; AI, Alar crease; Ea, external auditory canal; Ra, retroarticular crease; Oc, occiput; Ba, back; Mb, manubrium; Na, nare; Ac, antecubital fossa; Id, interdigital web; Pc, popliteal fossa; Ic, inguinal crease; Vf, volar forearm; Hp, hypothenar palm; Tw, toe webspace; Tn, toenail; Ph, plantar heel; Gb, Glabella. CNP0003934: Yi_2024, CNP0000635: Zhiming_2021, SRP002480: Julia_2016, SRP056821: KernRei_2016, SRP057859: Adrian_2017, CNP0002131: Zhiming_2023.
PCoA revealed significant differences in the virome composition among the six human cohorts from different populations (PERMANOVA R² = 10.2%, P < 0.001; Fig. 4b). Skin microenvironments and anatomic sites had a significant impact on the overall virome composition (PERMANOVA R² = 2.1% and R² = 8.5%, P < 0.001). Integrating data from various studies, sites, and microenvironments, we discovered that the most significant influence on the skin virome originated from the differences between studies, which accounted for 6.14% of the variation, followed by the impact of site at 3.14%, and the smallest was from the microenvironment at 0.08% (Table S6). We subsequently quantified the impact of different skin microenvironment and sites on the skin virome within each study and found that both different microenvironment and anatomical sites within the same microenvironment were significantly associated with virome composition (Fig. 4c and d). Furthermore, as individual heterogeneity was shown to be closely related to viral taxonomic variation, we evaluated the influence of age and sex on the skin virome across all the data. After adjusting for the effects of different studies, the impacts of age and sex were 0.08% and 0.06%, respectively (Table S6). Additionally, adjusting for the host’s age and sex did not significantly change the influence of different skin microenvironment and skin sites on the skin virome.
A compositional analysis at the viral family level across different skin microenvironments was conducted and found that unclassified vOTUs accounted for 10.35% of the viral sequences (Fig. 4e and f). Across all data sets, the virome of skin from different microenvironment and sites was dominated by members of the Caudoviricetes that have not been identified at the family level, as well as Papillomaviridae and Inoviridae (Fig. 4e and f). Oily areas were particularly enriched with the most abundant viruses on the skin, Papillomaviridae, while dry areas had a higher concentration of Autographiviridae and Inoviridae. Moist areas were enriched with viruses from the Herelleviridae (Fig. 4g; Table S7). Taken together, the viral composition on the skin surface exhibits microenvironment-specific and anatomic site-specific characteristics.
DISCUSSION
By constructing a comprehensive, non-redundant viral reference catalog from human skin samples across diverse microenvironments and geographic locations, our study has addressed several gaps in the current understanding of the skin virome. We discovered a significant number of vOTUs that do not match any known sequences, constituting approximately 90.85% (n = 2,610) of our catalog. This suggests that viral diversity is more abundant and varied than previously acknowledged. Our analysis also sheds light on the intricate relationships these viruses maintain with their host microenvironments and the broader skin microbiome, revealing how specific microenvironmental and anatomical sites might influence viral roles and adaptations.
A major challenge in our study stemmed from the metagenomic sequencing data, which, unlike more targeted approaches, did not undergo a VLP enrichment step nor included a host depletion step. The absence of VLP enrichment meant that our samples were not specifically concentrated for viruses, thereby including a substantial amount of non-viral genetic material. This lack of enrichment complicates the detection and characterization of viral entities, particularly those that are less abundant or that exist in complex biological matrices where they are overshadowed by more dominant microbial forms. Furthermore, the lack of a host depletion step prior to sequencing resulted in a high background noise of human DNA (61, 62). This background significantly complicates the bioinformatic processing required to distinguish viral DNA from the human DNA, which is predominant in skin samples (2, 5). Consequently, identifying and cataloging viral sequences required the implementation of stringent filtering criteria to ensure the accuracy and reliability of our viral sequence identifications.
Despite these efforts, the limitations inherent in our approach underscore the need for more refined methodologies in future studies. Employing VLP enrichment and host DNA depletion in the initial stages of sample preparation would likely allow for a more comprehensive and detailed exploration of the skin virome. These steps would not only increase the yield of viral sequences but also improve the signal-to-noise ratio, facilitating the detection of rare and novel viruses that could be crucial for understanding the full scope of viral impact on skin health and disease.
Additionally, our results highlight the potential impact of microenvironmental factors and sites on the composition of the virome. The variation in viral communities across different microenvironment and regions illustrates how these viruses may adapt to specific local microenvironments. This specificity to anatomical sites and microenvironment indicates that localized treatment approaches may be crucial, as the skin virome appears to vary with its surrounding conditions.
Conclusion
This extensive exploration of the human skin virome provides a deeper understanding of the viral components of the skin microbiome and highlights their potential role in influencing skin health and disease. As we continue to unravel the complex interactions within the skin microbiome, these findings not only pave the way for innovative therapeutic approaches but also enhance our understanding of skin biology. The way forward is to integrate these findings into broader biological and clinical research to fully realize the potential of the skin virome in promoting health and treating disease.
The reference list from the paper itself. Each links out to its DOI / PubMed record.
- 1Li Z, Xia J, Jiang L, Tan Y, An Y, Zhu X, Ruan J, Chen Z, Zhen H, Ma Y, Jie Z, Xiao L, Yang H, Wang J, Kristiansen K, Xu X, Jin L, Nie C, Krutmann J, Liu X, Wang J. 2021 a. Characterization of the human skin resistome and identification of two microbiota cutotypes. Microbiome 9:47. doi:10.1186/s 40168-020-00995-733597039 PMC 7890624 · doi ↗ · pubmed ↗
- 2Oh J, Byrd AL, Deming C, Conlan S, NISC Comparative Sequencing Program, Kong HH, Segre JA. 2014. Biogeography and individuality shape function in the human skin metagenome. Nature 514:59–64. doi:10.1038/nature 1378625279917 PMC 4185404 · doi ↗ · pubmed ↗
- 3Liu J, Yan R, Zhong Q, Ngo S, Bangayan NJ, Nguyen L, Lui T, Liu M, Erfe MC, Craft N, Tomida S, Li H. 2015. The diversity and host interactions of propionibacterium acnes bacteriophages on human skin. ISME J 9:2078–2093. doi:10.1038/ismej.2015.4725848871 PMC 4542041 · doi ↗ · pubmed ↗
- 4Belkaid Y, Tamoutounour S. 2016. The influence of skin microorganisms on cutaneous immunity. Nat Rev Immunol 16:353–366. doi:10.1038/nri.2016.4827231051 · doi ↗ · pubmed ↗
- 5Li Z, Xia J, Jiang L, Tan Y, An Y, Zhu X, Ruan J, Chen Z, Zhen H, Ma Y, Jie Z, Xiao L, Yang H, Wang J, Kristiansen K, Xu X, Jin L, Nie C, Krutmann J, Liu X, Wang J. 2021 b. Characterization of the human skin resistome and identification of two microbiota cutotypes. Microbiome 9:1–18. doi:10.1186/s 40168-020-00995-733597039 PMC 7890624 · doi ↗ · pubmed ↗
- 6Byrd AL, Belkaid Y, Segre JA. 2018. The human skin microbiome. Nat Rev Microbiol 16:143–155. doi:10.1038/nrmicro.2017.15729332945 · doi ↗ · pubmed ↗
- 7Chevallereau A, Pons BJ, van Houte S, Westra ER. 2022. Interactions between bacterial and phage communities in natural environments. Nat Rev Microbiol 20:49–62. doi:10.1038/s 41579-021-00602-y 34373631 · doi ↗ · pubmed ↗
- 8Erez Z, Steinberger-Levy I, Shamir M, Doron S, Stokar-Avihail A, Peleg Y, Melamed S, Leavitt A, Savidor A, Albeck S, Amitai G, Sorek R. 2017. Communication between viruses guides lysis-lysogeny decisions. Nature 541:488–493. doi:10.1038/nature 2104928099413 PMC 5378303 · doi ↗ · pubmed ↗