Charting host-microbe co-metabolism in skin aging and application to metagenomics data.

During aging of human skin, a number of intrinsic and extrinsic factors cause the alteration of the skin's structure, function and cutaneous physiology. Many studies have investigated the influence of the skin microbiome on these alterations, but the molecular mechanisms that dictate the interplay between these factors and the skin microbiome are still not fully understood. To obtain more insight into the connection between the skin microbiome and the human physiological processes involved in skin aging, we performed a systematic study on interconnected pathways of human and bacterial metabolic processes that are known to play a role in skin aging. The bacterial genes in these pathways were subsequently used to create Hidden Markov Models (HMMs), which were applied to screen for presence of defined functionalities in both genomic and metagenomic datasets of skin-associated bacteria. These models were further applied on 16S rRNA gene sequencing data from skin microbiota samples derived from female volunteers of two different age groups (25-28 years ('young') and 59-68 years ('old')). The results show that the main bacterial pathways associated with aging skin are those involved in the production of pigmentation intermediates, fatty acids and ceramides. This study furthermore provides evidence for a relation between skin aging and bacterial enzymes involved in protein glycation. Taken together, the results and insights described in this paper provide new leads for intervening with bacterial processes that are associated with aging of human skin.

Introduction

In the past decade, there has been increased appreciation of the influence of the host microbiome composition and functionalities on human cellular processes, such as gut and skin health and immunity, but also visible aging of the skin. During aging, the human skin changes both macroscopically and microscopically. Wrinkle formation, reduced elasticity and wound healing, as well as a reduced barrier function are hallmarks of this process. These phenotypic changes can be caused by intrinsic factors, such as pigmentation, anatomical variations and hormonal changes, and extrinsic factors, such as temperature, life style and exposure to smoke and sunlight. The molecular processes in human cells that are triggered by these signals range from immune related processes to molecular pathways that affect the structural integrity of the skin and the potential for rejuvenation [reviewed in 1].

The human skin is inhabited by a large number and variety of microorganisms, including bacteria, fungi and viruses. The impact of these organisms on skin health and skin aging becomes more and more evident with the maturation of 16S rRNA gene and metagenomics sequencing techniques that allow for the assessment of changes in the skin microbiome upon aging, onset of disease or therapeutic intervention. Recently, a number of studies have investigated the changes in the skin microbiome associated with skin aging and therapeutic intervention [2-14]. In general, the findings of these studies showed that the skin microbiome consists of a number of microbial genera that are consistently found in all humans (notably species from Staphylococcus, Cutibacterium, Corynebacterium and Acinetobacter) but that the exact composition is influenced by, amongst others, body site, gender, geographic location and age. In the study by Kim et al (2019), an overrepresentation of Alistipes, Prevotella, Porphyromonas, Sphingobacterium, Lactobacillus, Aerococcus, Oscillospira and Ruminococcus was found in the younger group (25–35 years old) of healthy female volunteers and an overrepresentation of Micrococcus, Corynebacterium, Dermacoccus, Actinomyces Streptococcus, Lysinibacillus and Bacillus in the older group (56–63 years old) [4]. Dimitriu and co-workers described an increase in Corynebacterium, Neisseriaceae, Chryseobacterium, Prevotella, Veillonella in older skin and also showed a direct correlation of a number of representatives of these bacterial families with wrinkles and the number of pigmented spots in skin [3]. Similarly, Juge and coworkers described higher levels of Corynebacterium and reduced levels of Propionibacterium in older skin (54–69 years old) [9].

Although over- and underrepresented bacterial taxa can be readily associated with younger or older skin, the limiting step and challenge in the analysis of these data is the interpretation of biological relevance in context of the research question or hypothesis. Specifically, little insight is currently available in literature on the interplay between microbial functionalities and human cellular processes (referred to as co-metabolism). In order to do this correctly, profound knowledge of the physiology of bacterial skin residents needs to be combined with knowledge on the physiology of the host that is related to the skin condition of interest.

Currently, analysis of such multidimensional data is performed by mapping to canonical knowledge that is stored in databases, for instance to determine enriched biological processes that are modified upon intervention or during the progression of a disease [15]. Whereas this is a valid approach, the value of the outcome is intrinsically limited by the lack of databases that store cross-domain knowledge on co-metabolism between microbiological and human cellular processes (e.g. KEGG and GMM) [16-18].

To increase our understanding of the influence of microbial functionalities (represented by the skin microbiome) on molecular processes involved in aging of human skin, there is a need for improved computational data analysis workflows that allow for integration of knowledge and data on both aspects. In this paper, we describe a conceptual framework to substantiate functional links between the human skin microbiome (cheek area) and cellular processes involved in age-related skin appearance. Similar workflows have previously been successfully employed [19, 20].

First, we constructed an expert-curated set of bacterial functionalities (genes) linked to skin aging-related intrinsic processes that was obtained through an extensive literature study. We applied bioinformatics approaches to identify these genes in publicly available bacterial genome sequences. To substantiate literature findings and confirm the link between specific microbial taxa to aging skin, skin microbiota compositional data was generated through 16S rRNA gene amplicon sequencing of microbial DNA isolated from cheek skin samples of female subjects with different age-related skin appearance. The functional potential of the defined microbial communities was inferred by combining compositional data with gene content information from reference genomes of skin bacteria. Finally, the described methodology was validated on a publicly available skin-related shotgun metagenomics dataset.

Materials and methods

Literature searches

Literature searches were carried out in June 2019 by running multiple searches on MedLine using the following search term " "skin aging" OR "skin rejuvenation" OR (skin AND "extrinsic aging") OR (skin AND "intrinsic aging") ". This yielded 8616 abstracts that were subsequently evaluated. Evaluation was done by manually screening the abstracts for statements on molecular pathways, metabolites, microorganisms and microbial genes. For a selected number of abstracts containing such statements, full text papers were collected and evaluated. From these full text papers, citations were screened and, where relevant and not yet retrieved in the first search, added to the set of references. Finally, from this set of literature references, gene symbols and/or locus tags encoding microbial functionalities were retrieved to be used as input for the generation of Hidden Markov Models (HMMs) (see below).

Hidden Markov Models

Reference genomes

For the relevant microbial functionalities that were retrieved from literature, reference genomes were collected from the NCBI sequence repository. Selection was based on availability of genome sequences of specific strains mentioned in the literature references. If the genome sequence of the strain mentioned was not available, a sequence from a representative strain from the same species was included in the reference genome set as an alternative. Alternative strains for which the source of isolation was marked as ‘skin’ were given priority over strains from a non-skin related or unknown origin. The final list of reference genomes is provided in S1 Table. A maximum-likelihood phylogenetic tree of selected reference genomes was generated with FastTree [21], based on concatenated multiple sequence alignment of all genes with exactly one copy on each of the genomes, as determined with OrthAgogue [22].

Generation of Hidden Markov Models

Bacterial functionalities with a common evolutionary origin in the reference genome set (see section ‘Reference genomes’) were identified using OrthAgogue. Multiple sequence alignments were generated with Muscle [23] for each orthologous group containing one or more genes deemed relevant for skin aging processes, as described under “Literature searches”. From these multiple sequence alignments, HMMs were constructed using hmmbuild (https://hmmer.org). The HMMs were designed specifically for skin microorganisms (in contrast to already existing signatures in publicly available databases that are based on a generic and broad range of microorganisms), which limits the application of these models to this specific application. For each model, a score threshold (GA) was determined through a heuristic approach and manual curation by taking the average of the lowest-scoring true positive and the highest-scoring true negative, using all proteins from the orthologous group as true positives and all other proteins from the reference genome set as true negatives. These scores were included in the individual HMM files as gathering threshold (GA). A workflow description for generating the HMMs, including the python script used for this study, are provided in a Github repository.

Relative abundance of bacterial functionalities associated with skin aging

To assess the abundance of bacterial genes in the identified pathways through application of the generated HMM’s (see section ‘Generation of Hidden Markov Models’), publicly available metagenomics datasets of a skin-derived and gut-derived microbial population was collected from literature [24, 25]. From these studies, samples were chosen based on the provided metadata of the subjects and/or sampling site that were included in these studies. The skin samples were chosen based on healthy volunteers collected through a similar sampling technique (swabs) on a similar sampling site (facial skin) as applied in this study (see section below). The gut samples were chosen based on healthy individuals who did not recently use antibiotics. Assembly of these reads was done using Spades (version 3.11.1) with the “—meta” flag enabled. Genes were called using Prodigal (version 2.6.2) [26]. Read mapping was performed by first making a nucleotide database using Bowtie2 (version 2.2.3) and then mapping the reads to the called genes using HUMANn2 (version 0.9.9) with the “—bypass-nucleotide-index", “—bypass-translated-search" and “—nucleotide-database" flags enabled. Bacterial functionalities associated with cellular processes involved in skin aging were identified in the bacterial genomes and shotgun metagenomics datasets by scanning protein-coding sequences (either as provided genome annotations, or as described above) using hmmsearch (version 3.1b2) (https://hmmer.org) with HMMs and score thresholds as described above. Overview tables with per-sample counts (either of genome sequences or shotgun metagenomics samples) of each HMM were generated using Python (https://python.org) and HUMANn2. Model level counts were generated by summing up the counts for the individual members in the model. By using two additional HUMANn2 functions (humann2_join_tables and humann2_renorm_table) the relative abundance was calculated between all samples, to diminish the differences in read number between samples and between the two different datasets. Finally, the relative abundance was normalized for each individual module by dividing all gut and skin samples within one module by the sample with the highest value. Through this approach, all samples within one module are relative to each other, between 0 and 1.

Sample collection and phenotypic assessment

To substantiate literature findings, skin microbiota compositional data was generated through 16S rRNA gene sequencing of microbial DNA isolated from skin swab samples (cheek area) of 2 x 25 healthy female volunteers of two different age groups (20–28 years old (‘young’) and 59–68 years old, (‘old’)).

Study design and skin swab collection

A human trial was performed for single sample collection of cheek samples from two groups of female subjects from the general European descent population (n = 25 per group) in Belgium. The subject age was between 20 and 28 years (younger population group), and between 59 and 68 years (older population group). Exclusion criteria were: nodulo-cystic lesions/acne or sebaceous gland condition; eczema, psoriasis, atopy; prescribed or unprescribed use of skin treatment within 1 month prior to inclusion (oral or topical antibiotic, antifungal or topical steroids); smokers or having smoked in the past 2 years; recent history of chronic alcohol consumption defined as more than 15 standard servings per week or more than 3 servings per day; tanning bed usage less than 1 month prior to inclusion; sunbathing 1 month prior to inclusion; habitual exposure to sun or use of a tanning bed; clear UV light effects in the younger subjects that might be related to aging (inspected by the principal investigator during the recruitment procedure); primary immunodeficiency patients known to have dysbiosis in community diversity; use of tanning dihydroxyacetone (DHA) less than 1 month prior to study start; excessive habitual caffeine use (more than 6 small/medium cups of coffee or soda daily); body mass index higher than 30; pregnancy or lactation. At the moment of inclusion, subjects were requested to: refrain from using OTC products for any kind of skin treatment 7 days prior to the sampling; refrain from any facial creams and make-up on the day of the sampling; refrain from swimming in a chlorinated pool, using a hot tub, sauna/steam baths, 48 hours prior to sampling visit; refrain from any facial treatments, facial masks, scrubs/peelings 2 weeks prior to the sampling visit; to follow specific bath/shower procedures (showering with plain water only, soap should be avoided as much as possible, and no scrubbing of the skin with a towel) during 24 hours prior to sampling visit–on the day of sampling showering was avoided completely (also with plain water).

The study was conducted according to the principles of the Declaration of Helsinki latest version Fortaleza, Brazil, October 2013. The candidates were informed verbally on the aim of the study and the study procedures. All participating subject signed the informed consent form. Data was analyzed anonymously. No IRB approval was requested for this study, for two reasons. Firstly, the study was carried out in Belgium, where IRB approval for this type of studies is not required. Secondly, there was no intervention or invasive procedure involved for the collection of the skin swab samples. Participants were requested to fill out a questionnaire (not considered as a psychological burden), and were only subjected to mild skin characteristics measurements.

During site visit information on the following parameters/metadata was collected: skin hydration (Corneometer, Courage&Khazaka); transepidermal water loss (TEWL) (Aquaflux, Biox); skin pH (Metrhom); skin smoothness, scaliness, sebum spots (Visioscan, Courage&Khazaka); standardized images were taken from the face (frontal and left view) with the Visia-CR (Canfield) to assess invisible and visible spots, red features, pores, porphyrins, color, fine lines, wrinkles and roughness. Additional metadata on skin properties was collected during the study by means of a questionnaire.

Samples were collected by means of swab sampling from the intact cheek skin. In short, skin sampling was performed using a standard custom-made sampling template, which allows for consistent sampling of 2 cm2 of the cheek areas (both sides used). Before sampling of each subject, the template was cleaned with ethanol (70%) and air dried. To further minimize sample cross-contamination, a fresh pair of sterile gloves was worn by the person sampling each individual. The area within the template was swabbed with a sterile HydraFlock collection swab (3206H-25; Puritan Diagnostics, USA). The collection swab was soaked in sterile PBS+ solution (PBS of pH 7.0 with 0.5% Tween-20). Next, the swab was dried very briefly on a sterile gauze before the start of the sampling. Samples were taken by direct swabbing of the cheek skin. The shaft of the swab was held parallel to the skin surface and it was rubbed back and forth 10 times applying firm pressure. Immediately after swabbing, each swab was swirled in a 1.5 ml collection tube with 0.25 ml of sterile PBS+ solution with 2.5 μl of 1 mM EDTA (DNase inhibitor). The swab was left in the solution and the samples were stored on ice before freezing and shipment on dry ice for further processing at NIZO (The Netherlands).

16S rRNA gene sequencing and bioinformatics

A detailed description of the DNA isolation, 16S rRNA gene sequencing, bioinformatics and statistical methods is provided in S1 Appendix. In short, sequencing data was analysed with a bioinformatics workflow based on Qiime 1.8 [27] and functional potential of the bacterial communities was inferred using a method analogous to PICRUSt [28]. Correlation of specific microbial taxa to skin age was primarily based on the skin aging score of the samples collected. This skin aging score was determined using eight skin age-related physiological measurements and was calculated as described in S1 Appendix. Raw sequencing data is available from the European Nucleotide Archive (ENA) with accession number PRJEB45035.

Results

To define the influence of microbial functionalities represented by the skin microbiome on molecular processes involved in aging of human skin, a multistep bioinformatics and data analysis approach was applied (Fig 1). First, scientific literature was evaluated to identify shared metabolic biological pathways involved in skin aging between humans and skin microorganisms. Based on the microbial pathways, Hidden Markov Models (HMMs) were constructed that were subsequently used to analyze the reference genomes from skin organisms, metagenomic data and 16S rRNA gene profiling data. The details of each step are described in the Materials and Methods section. Together, conceptual framework substantiates functional links between the human skin microbiome (cheek area) and cellular processes involved in age-related skin appearance as described in the sections below.

Fig 1

Schematic overview of the workflow applied this study.

The description of the multistep approach is linked to the results obtained and described in this manuscript (figure numbers).

Co-metabolic processes involved in skin aging

Through an extensive literature study, we selected molecular processes in human cells that are involved in skin aging for which associations with microbial functionalities in the skin microbiome were described. Host pathways were designated after the literature search and were used for a targeted investigation of samples from subjects included in this study. An overview is provided in Table 1 and detailed descriptions of each process and co-metabolism pathways involved are provided in S2 Appendix.

Table 1

Co-metabolic process involved in skin aging.

Host Process	Rationale	Targets for co-metabolism	Bacterial genes/ functionalities
UV-B induced immune suppression	UVB radiation causes immune suppression, mediated by urocanic acid.	urocanic acid	hutH
Histidine conversion	Histidine and other amino acids act as natural moisturizers on the skin and display antimicrobial activity.	histidine, other amino acids	hut-operon
Protein glycation	Glycation of collagen type 1 and other structural proteins is a major cause of loss of skin elasticity.	collagen, vimentin, elastin	fructokinases
Pigmentation	Production of melanin via the conversion of tyrosine and phenylalanine is a major pathway in pigmentation.	tyrosine, phenylalanine	tyrosinase
Ceramide metabolism	Incorrect ceramide metabolism has been associated with a decreased barrier function.	ceramides, sphingolipids, fatty acids	ceramide metabolism (cerN, sphR).
Fatty acid metabolism	Fatty acids are of major importance as carbon source, antimicrobial compounds and signalling molecules.	Fatty acids, in particular long chain fatty acids	Fatty acid metabolism, fatty acid biosynthesis, beta-oxidation, Fad operon.
Lipoteichoic acid signaling	LTA interacts with TLR2 and influences the recruitment of immune cells.	lipoteichoic acids	LTA biosynthesis (tag operon, LtaS)
Porphyrin synthesis	Increased porphyrin levels are often used as an indicator of skin aging.	porphyrins, heme	porphyrin biosynthesis, deoR Hem-operon
Proteolysis	Bacterial proteases can degrade structural proteins in the skin	collagen, elastin	lasA, lasB, sspA, sspB, sspC
Oxygen radical production and scavenging	Increased production of ROS and reduced capacity to scavenge these radicals.	oxygen radicals, glutathione	catalase, superoxide dismutase

Most important host processes involved in skin aging, including targets for assumed microbiome co-metabolism and associated bacterial genes or functionalities, as further described in S2 Appendix. The targets for co-metabolism can either be metabolites that are produced by human or bacterial cells, or human proteins (e.g. receptors) that can be targeted by the metabolites.

The pathways listed in Table 1 are centred around a single metabolite (e.g. urocanic acid), entire pathways (e.g. histidine conversion, ceramide biosynthesis) or a more general group of molecules (e.g. fructokinases that may play a role in glycation of multiple distinct proteins). It should be noted that this table is not a comprehensive reflection of the entire spectrum of metabolic pathways that play a role in bacterial-host interaction. For example, it is well known that skin bacteria interact with the host immune system by activation of Toll like receptors (TLRs) [29]. However, some important ligands for these TLRs, such as dsDNA, ssDNA or cell wall fragments, cannot be linked to specific bacterial genes or operons and are therefore not included in this overview.

Genomic distribution of genes involved in skin aging

To study the distribution of microbial functionalities defined in Table 1 across a set of relevant skin microbiota-specific microorganisms, genome sequences of such bacteria were obtained from NCBI as described in the Materials and Methods section (S1 Table). The proteins encoded by these genomes were scanned with HMMs representing the aging-related functionalities described above. For each organism the number of genes predicted to encode these functions was identified (Fig 2).

Fig 2

Skin aging module scores in reference genomes.

Different functionalities are shown in columns, the reference genomes are shown in rows. The columns are clustered based on modules within a known pathway (color-coded) and the individual modules from that pathway are shown on the x-axis. The number of genes predicted to encode the specified functions is shown on a scale from white to red, with red showing the highest number.

This analysis shows that subsets of related species are responsible for specific functionalities of interest, instead of equal distribution of the functionalities across the species. For example, species from the Anaerococcus and Lactobacillus group do in general score low on all the modules; they do not contain the full pathways for porphyrins metabolism or that code for proteins that are involved in proteolysis on the skin. This is in line with the fact that porphyrin metabolism has not been described in literature for species of these groups. In contrast, the porphyrins metabolism pathways are fully represented in species of Corynebacterium and Cutibacterium (formerly known as Propionibacterium) species (Fig 2). Another notable feature is the relative lack of genes for fatty-acid beta oxidation in Staphylococcus species, which is in agreement with the large differences previously described for fatty acid uptake and metabolism between Gram-negative and Gram-positive organisms and among Gram-positive organisms [30].

The illustration in Fig 2 highlights the fact that the skin microbiome as a whole performs functions that cannot be performed by a single species and that complex cross-interactions play a crucial role in understanding host-related processes.

Metagenomic distribution of genes involved in skin aging

To assess the abundance of the identified bacterial genes in the pathways, the HMMs were also applied on metagenomic datasets derived from gut [24] and skin samples [25] (Fig 3). The most notable differences observed between both dataset types are those in the modules related to histidine metabolism, porphyrin metabolism and proteolysis. Importantly, none of the modules is consistently fully covered within a single sample, and especially the fatty acid beta-oxidation is not well represented on the skin. This implicates that not all functionalities are yet fully presented in the current collection of HMMs.

Fig 3

Hits for the HMM modules on metagenomic datasets from gut and skin samples.

Different functionalities are shown in columns whereas the publicly available metagenomic datasets from gut and skin samples are shown in rows. The columns are clustered based on modules within a known pathway (color-coded) and the individual modules from that pathway are shown on the x-axis. The relative abundance of each functionality is shown on a scale from white to red, with red showing the highest level of abundance. Abundance was normalized to values between 0 and 1 within each module (columns) by dividing all gut and skin samples within one module by the sample with the highest value.

Distribution of genes in aging skin

To define bacterial microorganisms that are associated with younger or aging skin, we generated a microbiota compositional data set based on 16S rRNA gene sequencing obtained from microbial DNA isolated from skin samples (superficial layer, cheek area) from healthy, female subjects of two different age groups, as described in S1 Appendix. Analysis of this dataset shows a significant link between samples of subjects from different age groups (‘young’ and ‘old’) and bacterial composition on the OTU (Operational Taxonomic Unit) level (Fig 4). Inspection of the redundancy analysis (RDA) for variation in microbiome composition between the samples of both age groups, shows that, amongst others, Corynebacterium, Acinetobacter, Leptotrichia, Veillonella and Chryseobacterium associate with older skin.

Fig 4

Variation in microbiome composition of cheek skin samples from female volunteers of ‘young’ and ‘old’ age groups, based on OTUs.

Samples were taken from the superficial layer of the cheek. Age group (‘young’ and ‘old’) explains 1.1% of the variation in the microbiome. Separation of samples by age group was significant (p = 0.044). Samples from the ‘young’ and ‘old’ age group are indicated by black circles and blue squares, respectively. Arrows are plotted supplementary and represent 20 bacterial genera that are associated most with the age groups. The length and direction of the arrows indicate the relative strength of the association with either group.

To further quantify the correlation between microbiome composition and skin age, a skin aging (SA) score was calculated for each sample based on eight clinical measurements that represent phenotypic appearance of aging skin (namely pores number, roughness, wrinkles number, porphyrins number, red features number, skin color evenness, spots-visible number, and spots-invisible number) as described in S1 Appendix. In short, a high SA score means that there is a high score on the above parameters that are in general associated with an aging skin. This approach is similar to the methodology used by Dimitru and co-workers [3], with the exception that in the SA score more parameters are included. Fig 5A shows the link between the SA score and the microbiome composition for skin swab samples from female volunteers of the ‘young’ age group. Interestingly, this analysis already shows a significant link between the microbiota composition and SA scores. We only included the samples of the ‘young’ group in this analysis in order to remove the large effect of chronological age on the microbiome (Fig 3). This analysis thus shows a relation with the microbiome that is purely based on actual (age-related) appearance of the skin and not the defined chronological age. Interestingly, a higher SA score in this subject group is mainly driven by Propionibacterium and not by the organisms that were shown to be associated with chronological age (Fig 4).

Fig 5

RDA at the OTU and gene level for all cheek samples derived from subjects of the young age group.

A] RDA on OTU level. The SA score explains 6% of the variation in the microbiome. Separation of samples by SA score was statistically significant (p = 0.002). Blue gradient of the sample symbols indicates relative SA score value (dark blue = high SA score). B] RDA on microbial genes involved in skin aging processes for all cheek samples (superficial layer) derived from subjects of the ‘young’ age group. The SA score explains 7.7% of the variation in the predicted skin aging-related genes of the microbiome. Separation of samples by SA score was statistically significant (p = 0.026). Grey arrows represent the 15 predicted genes that show the highest association with low or high SA score.

Next, the functional potential of the defined skin microbiota was predicted by combining data on the presence and absence of bacterial genes linked to skin aging processes using the compositional profiles derived from the 16S rRNA gene sequence data. As shown in Fig 5B, a number of genes are significantly linked to the SA score, These genes mostly belong to i) the fatty acid beta-oxidation pathway (fadI, fadJ and fadB genes) and ceramide pathway (cerN), which are negatively correlated to the SA score, ii) the genes related to fatty acid biosynthesis (fabI, fabH and fabG genes) and iii) genes from histidine metabolism (hutH, hutU and hutL), which are positively correlated to the SA score. In addition, genes involved in response to oxidative stress (sodA, katA, katG) are related to the SA score, but the direction of this response is ambiguous; katG is associated with low SA score and the functionally overlapping katA is associated with a higher SA score.

Human pathways related to skin aging were identified for which metabolites were also present in molecular pathways of skin commensals. However, this does not necessarily imply that all bacterial pathways are also operational on the skin, or active in such a way that they actually influence the process of skin aging. Therefore, based on the counts for the individual genes of all skin aging associated pathways, a pathway score was calculated by adding up all the individual counts for each gene in a pathway. These pathway scores were subsequently compared to the SA score for the individual subjects by creating a linear model for each pathway (Fig 6 and S2 Table). This analysis revealed that only in the ‘young’ subjects there was a significant correlation between a number of pathways and the skin aging score of these subjects. For subjects in the ‘old’ age group, such a relation was not found (S2 Table). The microbial pathways for radical protection, histidine metabolism and fatty acid biosynthesis are significantly more prevalent in skin with a high SA score (i.e. older looking skin of young female subjects, S2 Table). This higher level pathway integration is in agreement with the redundancy plots in Fig 5 that show that the abundance of a number of the individual genes from these pathways were associated with the SA score.

Fig 6

Normalized functional pathway level scores for subjects of the young age group.

All subjects are in columns and were ordered based on SA score, in which a low SA score represents a young-looking skin and a high score corresponds to older-looking skin (from left to right). Functional pathways involved in skin aging (derived from Table 1) are in rows and were normalized to the maximum score per row.

Discussion

This study focused on an approach for biological interpretation of 16S rRNA gene and metagenomic profiling data in relation to the intrinsic and extrinsic processes that define skin aging. This approach is an addition to the large number of tools that currently exist for the analysis of community sequence data, including methods for sequence read mapping, statistical comparisons of multi-dimensional data and functional inference [reviewed in 31]. Our approach addresses the challenge of translating a list of OTUs and associated taxa to a description of functional pathways and molecules that can shed light on the actual mode of action of skin aging processes. Although we have exemplified its use on bacterial genomes, 16S rRNA gene and bacterial community sequencing data, the HMMs presented in this paper can also directly be used on metatranscriptomic data. The use of HMMs in sequence-based approaches has been used before [19, 20] but the construction of a HMM set for multiple skin aging-related pathways is new.

Application of the models described in this manuscript on a set of skin-related reference genomes shows that there is a clade-specific distribution of these functional pathways. Together with the fact that some essential pathways, such as fatty acid conversion, are not found in typical skin organisms such as S. epidermidis and S. warneri. This finding indicates that the set of modules described in this study is not yet complete and extension towards more clade-specific modules is warranted. The process described in this paper uses a manual process for creating these HMMs based on manual curation of a set of scientific articles. Integrating automated text mining methods [32] and automated HMM construction [33] could facilitate this process.

In this study, in which 16S rRNA gene sequencing data was profiled from a cohort of female subjects of ‘young’ and ‘old’ age groups, an association of (amongst others) Corynebacterium, Veillonella and Chryseobacterium with the old age group was found, which is in agreement with earlier reports [3, 4, 9]. However, when the microbiome was correlated with a set of phenotypic read outs aggregated in a skin aging (SA) score based on eight clinical parameters associated with older skin appearance, predominantly a higher relative abundance of Propionibacterium was associated with older-looking skin. This underscores the importance of selecting the most relevant phenotypic variable when running the statistical analysis.

When aggregating to a pathway level, the data from this study demonstrate that skin with a higher SA score is associated with an increased presence of bacterial metabolic pathways involved in radical protection and protein glycation, and underrepresentation of genes related to pigmentation and LTA synthesis. Based on these data, we assume that these microbial pathways play a role in the physiological processes for skin aging. However, inferring a causal relationship, or even inferring that these pathways are also functionally more active than other pathways in these groups, can only be justified when more data are available on the actual levels of metabolites and transcripts of these genes. Moreover, a more sophisticated method for scoring could be developed that takes into account the dependence and redundancy of gene functions in pathways, rather than simply counting the number of instances of these genes.

Taken together, the approach described in this study (Fig 1) should be regarded as a useful approach for generation of hypotheses with respect to the involvement of bacterial pathways in human physiological process, and as such generate leads for dedicated follow-up experiments. This approach can also be applied on data derived from other body sites (e.g. gut, oral cavity, etc), although the set of reference organisms to include should be adapted to reflect relevance to that particular body site. In case of application of the approach on different site on skin, the same reference organisms can be used as described in this manuscript as these are a relevant representation of organisms that can be found at multiple sites on the skin. In all cases, targeted verification in wet-lab studies will be required. These experiments could be directed at measuring the transcript levels of bacterial genes, profiling the metabolite profiles or in depth (intervention) assessment of the effect of host-related processes [34]. Together with an improved and extended collection of HMMs for skin aging-related processes, these experiments will provide the next step in understanding the molecular processes in skin aging and the role of the microbiome in these processes.

Reference genome list.

List of reference genomes used for the generation of Hidden Markov Models.

(XLSX)

Click here for additional data file.

Significant relationship between SA scores and pathway level scores in skin swab samples of ‘young’ and ‘old’ subjects.

(DOCX)

Click here for additional data file.

Materials and methods for phenotypic assessment.

Detailed description of materials and methods applied for the generation of a 16S rRNA gene dataset on the microbiota of cheek skin samples.

(PDF)

Click here for additional data file.

Co-metabolic processes involved in skin aging.

Detailed description of co-metabolic processes in human cells that are involved in skin aging with links to microbial functionalities in the skin microbiome.

(PDF)

Click here for additional data file.

Raw data HMM scripts.

(HMM)

Click here for additional data file.

Raw data HMM input file with gene symbols, locus tags, organisms and functions.

(XLSX)

Click here for additional data file.

Raw data file with metadata and relative abundances of OTUs, taxa and predicted functionalities.

(XLSX)

Click here for additional data file.

25 Mar 2021

PONE-D-20-39781

Charting host-microbe co-metabolism in skin aging and application to metagenomics data

PLOS ONE

Dear Dr. Kortman,

Thank you for submitting your manuscript to PLOS ONE. After careful consideration, we feel that it has merit but does not fully meet PLOS ONE’s publication criteria as it currently stands. Therefore, we invite you to submit a revised version of the manuscript that addresses the points raised during the review process.

Your manuscript has been evaluated by two experts in the field. They both found the study interesting but identified several points for improvement. Hence, I invite you to resubmit a revised version of your manuscript in which all points of the reviewers have been addressed. Of note, please make sure to use the term "16S rRNA gene" consistently. Terms like "16S sequencing" are not correct.

Please submit your revised manuscript by 10 May 2021. If you will need more time than this to complete your revisions, please reply to this message or contact the journal office at plosone@plos.org. When you're ready to submit your revision, log on to https://www.editorialmanager.com/pone/ and select the 'Submissions Needing Revision' folder to locate your manuscript file.

Please include the following items when submitting your revised manuscript:

A rebuttal letter that responds to each point raised by the academic editor and reviewer(s). You should upload this letter as a separate file labeled 'Response to Reviewers'.

A marked-up copy of your manuscript that highlights changes made to the original version. You should upload this as a separate file labeled 'Revised Manuscript with Track Changes'.

An unmarked version of your revised paper without tracked changes. You should upload this as a separate file labeled 'Manuscript'.

If you would like to make changes to your financial disclosure, please include your updated statement in your cover letter. Guidelines for resubmitting your figure files are available below the reviewer comments at the end of this letter.

If applicable, we recommend that you deposit your laboratory protocols in protocols.io to enhance the reproducibility of your results. Protocols.io assigns your protocol its own identifier (DOI) so that it can be cited independently in the future. For instructions see: http://journals.plos.org/plosone/s/submission-guidelines#loc-laboratory-protocols

We look forward to receiving your revised manuscript.

Kind regards,

Erwin G Zoetendal, PhD

Academic Editor

PLOS ONE

Journal Requirements:

When submitting your revision, we need you to address these additional requirements.

1. Please ensure that your manuscript meets PLOS ONE's style requirements, including those for file naming. The PLOS ONE style templates can be found at and

https://journals.plos.org/plosone/s/file?id=ba62/PLOSOne_formatting_sample_title_authors_affiliations.pdf

2. Thank you for stating the following in the Competing Interests section:

'The authors have declared that no competing interests exist.'

We note that one or more of the authors are employed by a commercial company: NIZO Food Research B.V,  Estée Lauder Companies & TenWise B.V.

Please provide an amended Funding Statement declaring this commercial affiliation, as well as a statement regarding the Role of Funders in your study. If the funding organization did not play a role in the study design, data collection and analysis, decision to publish, or preparation of the manuscript and only provided financial support in the form of authors' salaries and/or research materials, please review your statements relating to the author contributions, and ensure you have specifically and accurately indicated the role(s) that these authors had in your study. You can update author roles in the Author Contributions section of the online submission form.

Please also include the following statement within your amended Funding Statement.

“The funder provided support in the form of salaries for authors [insert relevant initials], but did not have any additional role in the study design, data collection and analysis, decision to publish, or preparation of the manuscript. The specific roles of these authors are articulated in the ‘author contributions’ section.”

If your commercial affiliation did play a role in your study, please state and explain this role within your updated Funding Statement.

2. Please also provide an updated Competing Interests Statement declaring this commercial affiliation along with any other relevant declarations relating to employment, consultancy, patents, products in development, or marketed products, etc.

Within your Competing Interests Statement, please confirm that this commercial affiliation does not alter your adherence to all PLOS ONE policies on sharing data and materials by including the following statement: "This does not alter our adherence to  PLOS ONE policies on sharing data and materials.” (as detailed online in our guide for authors http://journals.plos.org/plosone/s/competing-interests) . If this adherence statement is not accurate and  there are restrictions on sharing of data and/or materials, please state these. Please note that we cannot proceed with consideration of your article until this information has been declared.

Please include both an updated Funding Statement and Competing Interests Statement in your cover letter. We will change the online submission form on your behalf.

Please know it is PLOS ONE policy for corresponding authors to declare, on behalf of all authors, all potential competing interests for the purposes of transparency. PLOS defines a competing interest as anything that interferes with, or could reasonably be perceived as interfering with, the full and objective presentation, peer review, editorial decision-making, or publication of research or non-research articles submitted to one of the journals. Competing interests can be financial or non-financial, professional, or personal. Competing interests can arise in relationship to an organization or another person. Please follow this link to our website for more details on competing interests: http://journals.plos.org/plosone/s/competing-interests

3. We note that you have stated that you will provide repository information for your data at acceptance. Should your manuscript be accepted for publication, we will hold it until you provide the relevant accession numbers or DOIs necessary to access your data. If you wish to make changes to your Data Availability statement, please describe these changes in your cover letter and we will update your Data Availability statement to reflect the information you provide.

Additional Editor Comments (if provided):

[Note: HTML markup is below. Please do not edit.]

Reviewers' comments:

Reviewer's Responses to Questions

Comments to the Author

1. Is the manuscript technically sound, and do the data support the conclusions?

The manuscript must describe a technically sound piece of scientific research with data that supports the conclusions. Experiments must have been conducted rigorously, with appropriate controls, replication, and sample sizes. The conclusions must be drawn appropriately based on the data presented.

Reviewer #1: Partly

Reviewer #2: No

**********

2. Has the statistical analysis been performed appropriately and rigorously?

Reviewer #1: Yes

Reviewer #2: No

**********

3. Have the authors made all data underlying the findings in their manuscript fully available?

The PLOS Data policy requires authors to make all data underlying the findings described in their manuscript fully available without restriction, with rare exception (please refer to the Data Availability Statement in the manuscript PDF file). The data should be provided as part of the manuscript or its supporting information, or deposited to a public repository. For example, in addition to summary statistics, the data points behind means, medians and variance measures should be available. If there are restrictions on publicly sharing data—e.g. participant privacy or use of data from a third party—those must be specified.

Reviewer #1: Yes

Reviewer #2: No

**********

4. Is the manuscript presented in an intelligible fashion and written in standard English?

PLOS ONE does not copyedit accepted manuscripts, so the language in submitted articles must be clear, correct, and unambiguous. Any typographical or grammatical errors should be corrected at revision, so please note any specific errors here.

Reviewer #1: Yes

Reviewer #2: Yes

**********

5. Review Comments to the Author

Please use the space provided to explain your answers to the questions above. You may also include additional comments for the author, including concerns about dual publication, research ethics, or publication ethics. (Please upload your review as an attachment if it exceeds 20,000 characters)

Reviewer #1: This work investigates host-microbe co-metabolism (interconnected functional pathways) in skin aging based on metagenomics data. The work is based on Hidden Markov Models (HMM) that were designed to recognize activity of specific functionally relevant pathways from a combination of (host) genomic and (microbiome) metagenomic data sets on skin-associated bacteria. In addition, the models are applied to young and old adult females using predictions of the metagenomic composition, obtained by combining 16S profiles with sequenced microbial genome information. Specific pathway activity, associated with aging, is reported.

This work is based on a combination of well-established methods from biostatistics, bioinformatics and machine learning. The language and reporting are clear, appropriate citations are provided and I do not observe apparent flaws in the methodology. The work contributes a new algorithmic approach to host-associate microbiome research, where the innovation and benchmarking of such tools is necessary and valuable for the progress of the field. The discussion acknowledges some important limitations, such as the current need for substantial manual input during the method construction.

Main limitations of this work from my (methodological) perspective relate to the level of reporting regarding technical detail, and I comment on these below.

* Major

1) More details on the various steps in data processing and analysis could be added. In particular, a schematic figure (or figures) summarizing the approach and experiments would help to more easily understand the overall study setup and put the various pieces in context. The Methods section should then be compared to such overview in order to make sure that the necessary technical details of each step are provided so that the research is reproducible.

2) I did not find mention of the source code availability. It is paramount that the code for the experiments including data processing, method, and comparisons, would be available for further verification, benchmarking, and development.

3) The work would benefit from verification on another body site but I can see that this might be out of scope for this manuscript; it would be informative, however, to discuss in more length the possible limitations and challenges that one might encounter when extending the analysis to other body sites.

* Minor

All data are said to be fully available without restriction. But this is human data. Does this mean that also detailed personal subject data is available without restriction? Does this respect and guarantee anonymity of the volunteers to a sufficient extent?

The figure quality can be substantially improved but I trust this is done in the next rounds.

60: Kim et al. (year)?

66: Veilonella -> Veillonella?

123: "S1 Table" -> Table S1 (?)

233: "analogous to PICRUST" -> why not just state the name of this method itself? I am not sure if we should take any specific toolkit as a standard in this field

391: remove the word "important" and let the readers evaluate the importance

Reviewer #2: The authors of this manuscript have made an extensive literature review to identify bacterial species and corresponding genes in the human skin microbiota that can be associated to aging. Molecular processes in human cells associated with the corresponding functionalities are also identified.

A bioinformatic analysis follows to identify groups of orthologous genes among the selected genes. These sequences are then used to build protein signatures (Hidden Markov Models). A sample collection obtained from cheek skin from younger and older individuals is used to evaluate the performance of the potential of this signatures to characterize processes associated to aging. A gut microbiome dataset is also used to test the specificity of these signatures.

The manuscript is interesting and combines a review of existing work on skin microbiome associated to skin aging with the analysis of a de novo acquired data set. However, some of the results need to be changed or better explained as they currently don’t support all the presented conclusions.

Comments:

The manuscript is, to some extent exploiting. a circular argument: literature is mined to identify genes, bacterial genomes and host metabolic functions associated to skin aging. Then a data set of microbiome samples from old and young females is analysed to verify that this indeed the case.

While this approach has merit, it is difficult to see how this approach could be used to generate new hypothesis. This part in the discussion should be further expanded.

Hidden Markov Models: The authors present the generated collection of HMMs as a useful tool for further research. However this collection is not presented or described in any way (I was expecting to find some link to this collection, but I have not been able to find it). These should be presented and discussed.

Hidden Markov Models: The authors built a collection of bacterial genes (presented in table 1) through literature mining. From this collection of genes and the genome sequences of bacterial strains associated to skin aging they built Hidden Markov Models through multiple sequence alignments, ortholog search and the use of hmmer. While this seems a sensible approach I really don't understand why they don't use the already existing signatures that can be found in databases such as PFAM, Interpro, TIGRFAM and the like. The advantages of using these standardized signatures are obvious as they would provided a more consistent set of identifiers that would enable comparisons between and across datasets . Also the signatures on these databases have been built by considering a more extended collections of genomes should they could be of more use if different species/genera are to be analysed. The authors need to justify why they chose to disregard the already existing databases of standardized signatures.

Validation on skin derived and gut samples: Figure 2 compares the occurrence of these motifs in skin derived and gut derived microbiome samples. However the results don't clearly show these motifs to be more often present on the skin samples or in the gut samples. Analysis of Figure 2 shows that these motifs are also found in the gut microbiome. The authors should explain better to what extent Figure 2 represents a validation of the approach. Also the skin samples in Figure 2 don't seem to be specifically related to aging, so the question is to what extent they expect to have more occurrence of these motifs in these samples.

Figure 5 shows the link between the pathway level scores and the SA score however it seems that is only for the younger individuals. There seems to be some link, however more exhaustive anlaysis should be presented. One of the main questions is if there is a significant difference in these pathway level scores between individuals with low or high SA. Maybe a t-test or a correlation analysis would provide convincing evidence that would complement the plot on Figure 5. Also it is not clear why the analysis for the older individuals is not presented.

I find Appendix S2 to be very informative as it provides a nice review of literature findings. Table S1 provides the list of reference genomes used for this task. It would become more informative if the genes/functions reported in literature that made the authors select each genome were also be reported.

**********

6. PLOS authors have the option to publish the peer review history of their article (what does this mean?). If published, this will include your full peer review and any attached files.

If you choose “no”, your identity will remain anonymous but your review may still be made public.

Do you want your identity to be public for this peer review? For information about this choice, including consent withdrawal, please see our Privacy Policy.

Reviewer #1: Yes: Leo Lahti

Reviewer #2: No

[NOTE: If reviewer comments were submitted as an attachment file, they will be attached to this email and accessible via the submission site. Please log into your account, locate the manuscript record, and check for the action link "View Attachments". If this link does not appear, there are no attachment files.]

While revising your submission, please upload your figure files to the Preflight Analysis and Conversion Engine (PACE) digital diagnostic tool, https://pacev2.apexcovantage.com/. PACE helps ensure that figures meet PLOS requirements. To use PACE, you must first register as a user. Registration is free. Then, login and navigate to the UPLOAD tab, where you will find detailed instructions on how to use the tool. If you encounter any issues or have any questions when using PACE, please email PLOS at figures@plos.org. Please note that Supporting Information files do not need this step.

3 Jun 2021

Pleae find our responses in the cover letter and response to reviewers files.

Submitted filename: Response to reviewers.docx

Click here for additional data file.

11 Aug 2021

PONE-D-20-39781R1

Charting host-microbe co-metabolism in skin aging and application to metagenomics data

PLOS ONE

Dear Dr. Kortman,

Thank you for submitting your manuscript to PLOS ONE. After careful consideration, we feel that it has merit but does not fully meet PLOS ONE’s publication criteria as it currently stands. Therefore, we invite you to submit a revised version of the manuscript that addresses the points raised during the review process.

My apologies for the delayed decision, which is mainly due to the holiday period. Your manuscript has been evaluated by the same experts. As you can see in their comments, the manuscript was found significantly improved, but lacked important (background) information that hampered a proper detailed evaluation, such as the source code. Please check the PLOS ONE policy on code sharing and adapt accordingly in a revision: .

Please submit your revised manuscript by 30 September 2021. If you will need more time than this to complete your revisions, please reply to this message or contact the journal office at plosone@plos.org. When you're ready to submit your revision, log on to https://www.editorialmanager.com/pone/ and select the 'Submissions Needing Revision' folder to locate your manuscript file.

Please include the following items when submitting your revised manuscript:

A rebuttal letter that responds to each point raised by the academic editor and reviewer(s). You should upload this letter as a separate file labeled 'Response to Reviewers'.

A marked-up copy of your manuscript that highlights changes made to the original version. You should upload this as a separate file labeled 'Revised Manuscript with Track Changes'.

An unmarked version of your revised paper without tracked changes. You should upload this as a separate file labeled 'Manuscript'.

If you would like to make changes to your financial disclosure, please include your updated statement in your cover letter. Guidelines for resubmitting your figure files are available below the reviewer comments at the end of this letter.

If applicable, we recommend that you deposit your laboratory protocols in protocols.io to enhance the reproducibility of your results. Protocols.io assigns your protocol its own identifier (DOI) so that it can be cited independently in the future. For instructions see: http://journals.plos.org/plosone/s/submission-guidelines#loc-laboratory-protocols. Additionally, PLOS ONE offers an option for publishing peer-reviewed Lab Protocol articles, which describe protocols hosted on protocols.io. Read more information on sharing protocols at https://plos.org/protocols?utm_medium=editorial-email&utm_source=authorletters&utm_campaign=protocols.

We look forward to receiving your revised manuscript.

Kind regards,

Erwin G Zoetendal, PhD

Academic Editor

PLOS ONE

[Note: HTML markup is below. Please do not edit.]

Reviewers' comments:

Reviewer's Responses to Questions

Comments to the Author

1. If the authors have adequately addressed your comments raised in a previous round of review and you feel that this manuscript is now acceptable for publication, you may indicate that here to bypass the “Comments to the Author” section, enter your conflict of interest statement in the “Confidential to Editor” section, and submit your "Accept" recommendation.

Reviewer #1: (No Response)

Reviewer #2: (No Response)

**********

2. Is the manuscript technically sound, and do the data support the conclusions?

The manuscript must describe a technically sound piece of scientific research with data that supports the conclusions. Experiments must have been conducted rigorously, with appropriate controls, replication, and sample sizes. The conclusions must be drawn appropriately based on the data presented.

Reviewer #1: Partly

Reviewer #2: Yes

**********

3. Has the statistical analysis been performed appropriately and rigorously?

Reviewer #1: I Don't Know

Reviewer #2: Yes

**********

4. Have the authors made all data underlying the findings in their manuscript fully available?

The PLOS Data policy requires authors to make all data underlying the findings described in their manuscript fully available without restriction, with rare exception (please refer to the Data Availability Statement in the manuscript PDF file). The data should be provided as part of the manuscript or its supporting information, or deposited to a public repository. For example, in addition to summary statistics, the data points behind means, medians and variance measures should be available. If there are restrictions on publicly sharing data—e.g. participant privacy or use of data from a third party—those must be specified.

Reviewer #1: No

Reviewer #2: No

**********

5. Is the manuscript presented in an intelligible fashion and written in standard English?

PLOS ONE does not copyedit accepted manuscripts, so the language in submitted articles must be clear, correct, and unambiguous. Any typographical or grammatical errors should be corrected at revision, so please note any specific errors here.

Reviewer #1: Yes

Reviewer #2: Yes

**********

6. Review Comments to the Author

Please use the space provided to explain your answers to the questions above. You may also include additional comments for the author, including concerns about dual publication, research ethics, or publication ethics. (Please upload your review as an attachment if it exceeds 20,000 characters)

Reviewer #1: Source code for the experiments has not been provided. Sharing the source code is generally recommended as good practice in bioinformatics research. It is a common misconception that sharing source code equals to delivering a ready-made software product. These are two different concepts, source code for the experiments can be provided "as is", in the form that they were when performing the analyses. If important part of the analysis were conducted manually, and are not documented as source code, the reproducibility and replicability of the analyses is severely limited. The source code should be shared unless there are binding justifications for not doing so. The lack of source code prevents the further evaluation of the implementation reliability, accuracy, robustness, and reproducibility.

Reviewer #2: The authors have addressed all my comments. However an important part of my comment pertains the availability of additional data describing the HMMs and the input files. These files are mentioned in the description of the supporting information but are not available for download. Sometimes submission systems don't easily allow this type of supplementary files. If that is the case, I suggest the authors deposit them in some repository (such as gitlab or github) and include the link.

**********

7. PLOS authors have the option to publish the peer review history of their article (what does this mean?). If published, this will include your full peer review and any attached files.

If you choose “no”, your identity will remain anonymous but your review may still be made public.

Do you want your identity to be public for this peer review? For information about this choice, including consent withdrawal, please see our Privacy Policy.

Reviewer #1: No

Reviewer #2: No

[NOTE: If reviewer comments were submitted as an attachment file, they will be attached to this email and accessible via the submission site. Please log into your account, locate the manuscript record, and check for the action link "View Attachments". If this link does not appear, there are no attachment files.]

While revising your submission, please upload your figure files to the Preflight Analysis and Conversion Engine (PACE) digital diagnostic tool, https://pacev2.apexcovantage.com/. PACE helps ensure that figures meet PLOS requirements. To use PACE, you must first register as a user. Registration is free. Then, login and navigate to the UPLOAD tab, where you will find detailed instructions on how to use the tool. If you encounter any issues or have any questions when using PACE, please email PLOS at figures@plos.org. Please note that Supporting Information files do not need this step.

17 Sep 2021

Please find our response in the files "Response to reviewers" and "Cover letter resubmission".

Submitted filename: Response to reviewers.docx

Click here for additional data file.

11 Oct 2021

Charting host-microbe co-metabolism in skin aging and application to metagenomics data

PONE-D-20-39781R2

Dear Dr. Kortman,

We’re pleased to inform you that your manuscript has been judged scientifically suitable for publication and will be formally accepted for publication once it meets all outstanding technical requirements.

Within one week, you’ll receive an e-mail detailing the required amendments. When these have been addressed, you’ll receive a formal acceptance letter and your manuscript will be scheduled for publication.

N.B. I noticed that the older age group is described as 59-65 under heading “Sample collection and phenotypic assessment” while it is described as 59-68 a few lines later under heading “Study design and skin swab collection”. Make sure to correct this difference in the proof.

An invoice for payment will follow shortly after the formal acceptance. To ensure an efficient process, please log into Editorial Manager at http://www.editorialmanager.com/pone/, click the 'Update My Information' link at the top of the page, and double check that your user information is up-to-date. If you have any billing related questions, please contact our Author Billing department directly at authorbilling@plos.org.

If your institution or institutions have a press office, please notify them about your upcoming paper to help maximize its impact. If they’ll be preparing press materials, please inform our press team as soon as possible -- no later than 48 hours after receiving the formal acceptance. Your manuscript will remain under strict press embargo until 2 pm Eastern Time on the date of publication. For more information, please contact onepress@plos.org.

Kind regards,

Erwin G Zoetendal, PhD

Academic Editor

PLOS ONE

Additional Editor Comments (optional):

Reviewers' comments:

Reviewer's Responses to Questions

Comments to the Author

1. If the authors have adequately addressed your comments raised in a previous round of review and you feel that this manuscript is now acceptable for publication, you may indicate that here to bypass the “Comments to the Author” section, enter your conflict of interest statement in the “Confidential to Editor” section, and submit your "Accept" recommendation.

Reviewer #2: All comments have been addressed

**********

2. Is the manuscript technically sound, and do the data support the conclusions?

The manuscript must describe a technically sound piece of scientific research with data that supports the conclusions. Experiments must have been conducted rigorously, with appropriate controls, replication, and sample sizes. The conclusions must be drawn appropriately based on the data presented.

Reviewer #2: (No Response)

**********

3. Has the statistical analysis been performed appropriately and rigorously?

Reviewer #2: (No Response)

**********

4. Have the authors made all data underlying the findings in their manuscript fully available?

The PLOS Data policy requires authors to make all data underlying the findings described in their manuscript fully available without restriction, with rare exception (please refer to the Data Availability Statement in the manuscript PDF file). The data should be provided as part of the manuscript or its supporting information, or deposited to a public repository. For example, in addition to summary statistics, the data points behind means, medians and variance measures should be available. If there are restrictions on publicly sharing data—e.g. participant privacy or use of data from a third party—those must be specified.

Reviewer #2: Yes

**********

5. Is the manuscript presented in an intelligible fashion and written in standard English?

PLOS ONE does not copyedit accepted manuscripts, so the language in submitted articles must be clear, correct, and unambiguous. Any typographical or grammatical errors should be corrected at revision, so please note any specific errors here.

Reviewer #2: (No Response)

**********

6. Review Comments to the Author

Please use the space provided to explain your answers to the questions above. You may also include additional comments for the author, including concerns about dual publication, research ethics, or publication ethics. (Please upload your review as an attachment if it exceeds 20,000 characters)

Reviewer #2: (No Response)

**********

7. PLOS authors have the option to publish the peer review history of their article (what does this mean?). If published, this will include your full peer review and any attached files.

If you choose “no”, your identity will remain anonymous but your review may still be made public.

Do you want your identity to be public for this peer review? For information about this choice, including consent withdrawal, please see our Privacy Policy.

Reviewer #2: No

18 Oct 2021

PONE-D-20-39781R2

Charting host-microbe co-metabolism in skin aging and application to metagenomics data

Dear Dr. Kortman:

I'm pleased to inform you that your manuscript has been deemed suitable for publication in PLOS ONE. Congratulations! Your manuscript is now with our production department.

If your institution or institutions have a press office, please let them know about your upcoming paper now to help maximize its impact. If they'll be preparing press materials, please inform our press team within the next 48 hours. Your manuscript will remain under strict press embargo until 2 pm Eastern Time on the date of publication. For more information please contact onepress@plos.org.

If we can help with anything else, please email us at plosone@plos.org.

Thank you for submitting your work to PLOS ONE and supporting open access.

Kind regards,

PLOS ONE Editorial Office Staff

on behalf of

Dr. Erwin G Zoetendal

Academic Editor

PLOS ONE