Melanocytes in the Skin – Comparative Whole Transcriptome Analysis of Main Skin Cell Types

Abstract

Melanocytes possess several functions besides a role in pigment synthesis, but detailed characteristics of the cells are still unclear. We used whole transcriptome sequencing (RNA-Seq) to assess differential gene expression of cultivated normal human melanocytes with respect to keratinocytes, fibroblasts and whole skin. The present results reveal cultivated melanocytes as highly proliferative cells with possible stem cell-like properties. The enhanced readiness to regenerate makes melanocytes the most vulnerable cells in the skin and explains their high risk of developing into malignant melanoma.

About the Authors:

Paula Reemann

* E-mail: [email protected]

Affiliation: Department of Physiology, Institute of Biomedicine and Translational Medicine, University of Tartu, Tartu, Estonia

Ene Reimann

Affiliations Core Facility of Clinical Genomics, Department of Pathophysiology, Department of Biomedicine, Institute of Biomedicine and Translational Medicine, University of Tartu, Tartu, Estonia, The Institute of Veterinary Medicine and Animal Sciences of the Estonian University of Life Sciences, Competence Centre on Reproductive Medicine and Biology, Tartu, Estonia

Sten Ilmjärv

Affiliation: Department of Physiology, Institute of Biomedicine and Translational Medicine, University of Tartu and Quretec Ltd (private limited company), Tartu, Estonia

Orm Porosaar

Affiliation: Department of Pediatric Surgery, Tallinn Children's Hospital, Tallinn, Estonia

Helgi Silm

Affiliation: Clinic of Dermatology, Tartu University Hospital, Department of Dermatology, University of Tartu, Tartu, Estonia

Viljar Jaks

Affiliation: Institute of Molecular and Cell Biology, University of Tartu, Tartu, Estonia

Eero Vasar

Affiliation: Department of Physiology, Institute of Biomedicine and Translational Medicine, University of Tartu and Centre of Translational Medicine, University of Tartu, Tartu, Estonia

Külli Kingo

Affiliation: Clinic of Dermatology, Tartu University Hospital, Department of Dermatology, University of Tartu, Tartu, Estonia

Sulev Kõks

Affiliations The Institute of Veterinary Medicine and Animal Sciences of the Estonian University of Life Sciences, Competence Centre on Reproductive Medicine and Biology, Tartu, Estonia, Department of Pathophysiology, Institute of Biomedicine and Translational Medicine, University of Tartu and Centre of Translational Medicine, University of Tartu, Tartu, Estonia

Introduction

Skin is a highly organized and differentiated structure, which consist of various cell types. Keratinocytes (KC) and fibroblasts (FB) together form the majority of cellular components in the skin (7×105–9×105 KC per mm2 [1] and 4×103 mid-dermis FB to 105 papillary FB per mm3 [2], [3]). Therefore the functional properties of the highly outnumbered melanocytes (MC) have received relatively less attention. The average number of pigment-producing MC depends on the body site and is estimated to be between 500 to 2000 MC per mm2 [1], [4]. Interestingly, differences in ethnic background manifest in the intensity of melanogenesis and the morphology of dendrites, but not in the quantity of MC [5]. Despite their small number, MC have proven to have several roles besides melanogenesis, a well-characterized property of MC. They are able to secrete a wide range of signaling molecules, e.g. proinflammatory cytokines, immunosuppressive molecules, neuromediators etc. [2], [6]–[9]. MC interacts highly with surrounding KC, which have been shown to regulate MC survival, dendricity, melanogenesis, and the expression of cell surface receptors [10].

Numerous gene expression analyses of different skin cell populations have been performed in both physiological and pathological states using an array of detection techniques ranging from quantitative real time polymerase chain reaction (qPCR) and in situ hybridization to high throughput methods such as serial analysis of gene expression and microarrays [11]–[13]. However, all these methods have specific limitations. In contrast, the use of high-throughput RNA-Seq on rRNA-depleted samples allows the detection of nearly all coding and non-coding RNA species in a given sample.

In the present study we outline, to our knowledge for the first time, the differences of MC compared to other main cell types of the skin at the level of complete transcriptome. We used whole skin samples and cultivated primary skin cells, harvested from the same body site of healthy subjects of similar age.

Materials and Methods

Ethics Statement and Patients

All procedures were carried out in accordance with the ethical standards. This study (including written consent form) has been approved by Research Ethics Committee of the University of Tartu (approval number 178/T-19). The patients with no concurrent diseases and signs of infection, were recruited from among elective patients present at the Department of Pediatric Surgery, Tallinn Children's Hospital. A written informed consent was obtained from all parents or caretakers of patients under 18 years. Additionally, separate written informed consent was obtained from all patients aged 8–17 years. Nine pediatric foreskins from healthy volunteers (aged 5 months to 10 year) were used. Collected tissue samples were by-products of circumcise procedures and no additional intervention was caused by our investigation. All samples were coded and information of the donor identity was only available for the physician. All patient related information was stored separately from the samples and the data.

Cell culture

From each tissue sample three skin cell types (keratinocytes, melanocytes and fibroblasts) were harvested. After rinsing in phosphate-buffered saline (PBS w/o Ca, mg, PAA Laboratories), subcutaneous fat was removed from skin pieces and tissues were incubated in dispase II (2.4 U/ml, Sigma-Aldrich) at +4°C overnight. Epidermis was peeled off from the dermis, transferred into 0.05% trypsin/0.02% EDTA (Life Technologies) for 3 min at 37°C. Enzymatic digestion was stopped by adding the trypsin inhibitor (Sigma-Aldrich). MC culture dishes were precoated with gelatin before cell seeding. EpiLife basal medium with human keratinocyte growth supplement (Life Technologies) and melanocyte growth medium M2 with supplement mix (PromoCell) were used to cultivate corresponding cells.

A piece of dermis was used for isolation of FB by migration method as follows. The dermis was cut into 4×4 mm pieces and attached onto a culture dish, covered with 10 ml Dulbecco's Modified of Eagle's Medium (DMEM) (PAA Laboratories) supplemented with10% foetal bovine serum (Sigma-Aldrich). The medium was changed every 2nd day throughout the study. The skin cells were cultivated at low passage number (2–3) to minimize the influence of culturing conditions.

RNA Extraction and Library Preparation and Whole Transcriptome Sequencing

Cultivated skin cells underwent to lysis and RNA extraction process using Trizol® (Invitrogen) as described in [12]. The skin biopsies were homogenized using Precellys® 24 system and previously optimized protocol [12] were used for total RNA extraction, followed by DNAse I (Qiagen) treatment. The purity and concentration of samples was checked with both Qubit spectrophotometer and Nano Drop ND-1000 and the integrity of the RNA (RIN) was evaluated using Agilent 2100 Bioanalyzer. We chose 12 total RNA samples, with the highest RIN (9–10), extracted from 4 KC, 4 MC, 2 FB and 2 whole skin samples for library preparation. Extracted mRNA was enriched using RiboMinusTM Eukaryote kit (Invitrogen) according to manufacturer's instructions. The final quantity of RNA was 10 µg per reaction. The cDNA library with size-selected in the range of 150–250 bp and following bar-coding preparation procedure was based on a protocol provided by Applied Biosystems. Samples were sequenced using SOLiD 5500xl platform with 75 bp forward and 35 reverse primers.

Analysis of RNA-Seq data

Sequencing of cDNA libraries resulted in 24842284 to 44324428 paired reads per sample. For greater mapping quality the initial 75 bp F3 and 35 bp F5 reads were trimmed to 45 and 25 base pairs, respectively. All color-spaced reads were aligned to human reference genome (Ensembl, release 73) using TopHat v2.1.0 [14] that used Bowtie version 1.0.0 [15]. RPKM (reads per kilobase of transcript per million mapped reads) values for gene expression levels were calculated with Cufflinks v2.0.2 [16] and raw counts were retrieved with HTSeq version 0.5.3p9 (http://www-huber.embl.de/users/anders/HTSeq/) using gene annotations of protein coding genes downloaded from Ensembl (release 73). Differential expression was estimated on raw counts with edgeR [17]. All programs were used with their default parameters with TopHat set to not to find novel junctions.

Modeling background regions

To estimate the number of truly expressed genes we modeled intergenic regions using a methodology described in [18]. Models of intergenic regions are expected to reflect the level of background expression (noise), which is taken as the baseline when estimating the number of expressed genes. For each gene, the length of the background region was equal to the gene's longest combined transcript (the sum of all transcribed nucleotides) and it extended upstream from position-1000 relative to the transcription start site. Only background regions that did not overlap with any expressed sequence tags (EST) were used in the analysis. A gene was considered as expressed only if the RPKM value in all samples of the corresponding cell type was above the cutoff (0.95). Conversely, the gene was labeled as not expressed if the RPKM value was below 0.95 in at least one of the samples. EST annotations were downloaded using UCSC Table Browser (http://genome.ucsc.edu/cgi-bin/hgTables).

Differential expression analysis of gene expression

Differential expression was estimated between MC and samples from KC, FB and the whole skin using edgeR [17]. A gene was considered as differentially expressed if the FDR-adjusted p-value was below 0.05 and if the gene was expressed in at least one of the cell types.

To identify a gene as expressed only in MC and not in KC, FB (termed as “uniquely expressed in MC”) it had to meet the following requirements: 1) RPKM >0.95 in MC (gene is labeled as expressed in MC), 2) RPKM <0.95 in KC, FB (gene is labeled as not expressed in KC, FB), 3) differential expression FDR <0.05 (gene is differentially expressed in MC with respect to KC and FB).

Pathway analysis of differentially expressed genes

Gene ontology enrichment analysis of differentially expressed genes was performed using g:GOSt (http://biit.cs.ut.ee/gprofiler/index.cgi) [19]. Bases on the p-values of g:GOSt analysis, GOsummaries package were used to generate wordclouds of gene names (http://cran.r-project.org/web/packages/GOsummaries/index.html). The word sizes in wordclouds are defined by the p-values.

Additionally, multidimensional scaling test for visualizing the level of similarity of individual samples in study groups were performed using edgeR. The results confirm the homogenity and purity of cell populations (S1 Fig.).

Results and Discussion

Overall differences between cultivated MC, KC, FB and whole skin tissue

As expected, the total number of expressed genes was the highest in whole skin samples (10,871 genes), since other cell types besides KC, MC and FB (epithelial cells, Merkel cells, Langerhans cells etc.) are found in a skin biopsy (Table 1). Interestingly, the total number of expressed genes was the lowest in KC (Table 1). Thereat, 7766 genes were commonly expressed in all study groups (whole skin and KC, MC and FB). The list of the detected genes and their RPKM values can be found in the S1 Table. Similarly, the number of genes considered as uniquely expressed was the largest in the whole skin sample (290 genes, Table 1). When comparing specific cell populations, FB had a higher number of uniquely expressed genes compared to MC and KC (277, 122, 138 uniquely expressed genes, respectively) (Table 1). It is likely that the true number of uniquely expressed genes is higher as we applied a relatively strict cutoff criterion (RPKM >0.95) when calling gene expression as present or absent. The genes, uniquely expressed in MC are displayed in Table 2.

[Figure omitted. See PDF.]

Table 1. The number of detected and uniquely expressed genes.

https://doi.org/10.1371/journal.pone.0115717.t001

[Figure omitted. See PDF.]

Table 2. Uniquely expressed genes in MC.

https://doi.org/10.1371/journal.pone.0115717.t002

Based on the differential gene expression analysis, we identified significantly fewer similarities between MC and whole skin gene expression patterns than when comparing KC or FB to the whole skin. In melanocytes, 6231 genes were differentially expressed compared to the whole skin. Of those, 3680 were upregulated in MC and 2551 downregulated in MC with respect to whole skin samples (Fig. 1). The number of differentially expressed genes with respect to whole skin was 4480 in KC and 4454 in FB. This finding can be explained by the relatively small proportion of MC in the total cell number of skin. The entire list of differentially expressed genes can be found in S2 Table.

[Figure omitted. See PDF.]

Figure 1. The number of differentially expressed genes in each study group - melanocyte (MC), keratinocyte (KC), fibroblast (FB) and whole skin (Skin).

Red triangles – upregulated genes, blue triangles – downregulated genes.

https://doi.org/10.1371/journal.pone.0115717.g001

Differential gene expression in MC

Previous studies have mainly pointed out the role of melanocytes as pigment-producing cells in the skin. The gene expression profile of MC has been compared to other normal skin cells but also to pathologic melanoma cells [11], [13]. We confirmed the high expression level of previously identified melanocyte-specific genes such as DCT, TYR, KIT, EDNRB, MITF, and TYRP1, in MC compared to FB and KC (S2 Table), as reported also in a previous microarray study [11]. DCT, TYR and TYRP1 encode enzymes acting in the melanin synthesis pathway [20]. EDNRB and MITF are crucial for melanoblast proliferation and KIT is needed for the differentiation of melanoblasts into TYR-positive MC [21], [22]. We also showed MC1R, PLA1A, NPM2 to be uniquely expressed in MC, but not in KC or FB (Table 2) corroborating previously published data [23]. MC1R, a receptor for melanocyte-stimulating hormones and adrenocorticotropic hormone, is involved in regulating the pigmentation of the skin and hair. PLA1A and NPM2 have both been shown to be essential in melanoma progression [24], [25]. In agreement with previous studies [11], [13] the expression of ABCC2, DNAJA1, GPR143, MLANA, OCA2, QPCT, RRAGD, TBC1D7 and GPR137B was detected at a higher level in MC compared to KC and FB (S2 Table). Majority of these genes are also related to the melanogenesis pathway, controlling the growth and maturation of melanosomes, being involved in melanosome biogenesis, transporting melanin precursor molecules or being candidate genes in melanocytic tumor progression [26]–[30].

In the present study our goal was to detect other potential functions and key genes of MC. Comparing MC, KC and FB, we found several gene groups, which distinguish MC from other skin cells.

The most prominent group of genes differentially expressed in MC compared to FB (Fig. 1) was that of encoding various histone proteins (HIST1H1A, HIST1H1B, HIST1H2AA, HIST1H2AE, HIST1H2AG, HIST1H2AI, HIST1H2BB, HIST1H2BH, HIST1H2BI, HIST1H2BN, HIST1H3D, HIST1H3F, HIST1H3H, HIST1H3J, HIST1H4D, HIST1H4I, HIST1H4L, HIST2H3D, HIST3H2A, HIST3H2BB and HIST4H4) (S2 Table). Another group of genes, the expression of which drastically differed in MC compared to FB, were cell division cycle protein genes (CDC20, CDC25A, CDC25C, CDC6, CDCA2, CDCA5, CDCA8) (S2 Table). Additionally, a set of kinesin family genes (KIF13B, KIF20B, KIF21A, KIF22, KIF24, KIF2C, KIFC1 was differentially expressed in MC compared to FB (S2 Table). We saw also significantly higher expression of calcium-binding proteins S100A1, S100A14, S100A8, S100A9 and S100B genes in MC compared to FB (S2 Table).

Above mentioned gene groups, which are prominently expressed in MC but not in FB suggest that MC are active and intensively proliferative in the culture conditions. For instance, the high prevalence of histone genes in MC reveal to intensified DNA synthesis as histones are responsible for nucleosome structure and proper DNA wrapping [31]. This corresponds well with the high proliferative activity we saw during the cell cultivation process. This conclusion is substantiated by the increased expression of a number of cell cycle regulating genes, which are involved in the regulation of cell cycle at several steps and levels and kinesin genes, which are related to cell movements and intracellular trafficking, including chromosome and centrosome positioning during mitosis [32]. Calcium-binding proteins are responsible for numerous cellular processes, e.g. cell cycle regulation and differentiation, but have also been suggested to have tumor suppressor functions and are highly expressed in cells with stem cell properties [29], [33].

Compared to KC, MC expressed a set of major histocompatibility complex protein genes such as HLA-B, HLA-DMA, HLA-DPB1, HLA-DRA and HLA-F (S2 Table). In addition genes, which are related to viral and bacterial infection defense mechanism - interferon induced protein genes [34] IFI27, IFI35, IFI44, IFI44L, IFI6, IFIT1, IFIT2, IFIT3, IFITM1, IFITM2 and IFITM3, were highly expressed in MC compared to KC (S2 Table). KC have been shown to be the key players in modulating the immunological status of physiological and pathological skin; being the first sensors of harmful agents, they secrete inhibitory and stimulating cytokines, and activate other immune competent cells (e.g. Langerhans cells) [35]. Our data suggest that the role of MC in cutaneous immune system regulation might be more extensive than anticipated so far.

Importantly, we could identify a set of MC-specific genes previously not described by other researchers. In many cases the functions of these genes are poorly characterized and further experiments are needed to identify their precise role in MC (Table 2).

A number of genes that were unique for MC (CA8, CHRNA6, CTTNBP2, EPHA5, FAXC, KCNN2, SCG2, SLITRK2) (Table 2) refer to MC's origin from the neural crest [36]–[44]. Additionally, a number of genes (ANO5, CGREF1, EGFL8, ENTHD1 and ITPR1) involved in calcium-mediated processes were also uniquely expressed in MC (Table 2).

Next we identified a number of genes that were not uniquely expressed in MC but were specifically upregulated in melanocyte culture. Their specific role in melanocytes is unclear, but based on the existing biological data these genes can be divided into following functional classes: tumorogenesis, inflammation and stem cell related genes.

Melanocytes and tumorogenesis.

Genes belonging to pathways involved in tumor progression are more characteristic to MC than to KC and FB. High expression of tumor suppressive genes [45]–[48] and its candidate genes such as ARL9, C10orf90, SLC22A18, DAPK1, BEX1, PYHIN1, IGSF8 could be observed in MC compared to KC and FB (S2 Table). Interestingly, a novel tumor suppressor IGLON family genes (IGLON2, IGLON 3 and IGLON4) were all detectable in MC. IGLON3 also known as LSAMP was prominently expressed in MC compared to KC and FB (S2 Table). IGLON family genes are mainly described as regulators of nerve growth factors but recent studies have shown their strong tumor suppressive capacities [49].

A number of genes, which normally have a role in growth and dividing processes or in apoptosis regulation, can play a role in cancer formation. For example several RING-type zinc finger proteins, which function to avoid uncontrolled proliferation and to be a part of embryonic development, act as cancer development modulators [50]. As an example, certain RING finger proteins RNF144A, RNF157 and RNF187 were specifically upregulated in MC (S2 Table). DNAJA1 is increased in MC (S2 Table), which overexpression reduces cancer cell survival [51].

Interestingly, laminin alpha 1 gene LAMA1 was highly expressed in MC compared to KC and FB (S2 Table). Laminins are integral parts of the extracellular matrix [52]. LAMA1 is present mostly in the early stages in most tissues of the embryo and is not common for adult tissues [53]. These results suggest that MC might have an important role in basement membrane formation and remodeling and might refer to a higher potential of MC to transform into tumorous cells.

The family of sialyltransferases, which comprises a large group of enzymes, responsible for the synthesis of sialylated glycans, regulates immune response including virus binding [54]. Sialylated glycans can be found on the surface of many tumor cells where they counteract the recognition of malignant cells by the immune system [55]. Our analysis identified several sialyltransferase genes ST3GAL4, ST3GAL5, ST3GAL6, ST6GAL1, ST6GALNAC3, ST8SIA1 and ST8SIA6, which were expressed at a higher level in MC compared to KC and FB (S2 Table).

Further, we confirmed a differential expression of genes, which have already been shown to be related to melanoma development. Such as chemokines, which major role is to guide the migration of cells and mediate immune response are important for tumor invasion and metastatic behavior [56]. We also showed that CXCL5, CCL28 and chemokine-like protein FAM19A5 were significantly more expressed in MC compared to KC and FB (S2 Table). Additionally, a few genes, which regulate angiogenesis (semaphorins SEMA4C and SEMA6A, matrix metalloproteinases MMP8 and MMP17 [57]–[59] and are thereby essential for malignant processes, had a higher expression level in MC compared to KC and FB (S2 Table). ABCC2, ABCB5 and ABCB6, which were also highly expressed in MC compared to KC and FB (S2 Table), are the members of the ABC transporter family. ABC proteins have been shown to be involved in multidrug resistance in cancer treatment, but they also promote the pluripotency of embryonic cells and sustain the self-renewal of stem cells [60].

Melanocytes and inflammation.

Susceptibility of cells to malignancies is strongly connected both inflammatory processes, but also their stem-cell-like properties. Inflammation influences cancer development at different levels - predisposing precancerosis, misdirect immune system, initiating invasion process etc. [61].

A glaring example about the relation of immune response and tumorogenesis is a family of tumor necrosis factors (TNFs) and their receptors. Being strongly engaged both in immune system modulation and apoptosis regulation, they trigger infiltration of inflammatory cells into tumorous tissue [62]. The cross-regulation of TNF and interferon regulatory factors have been proposed recently [63]. In line with this, the tumor necrosis factors receptors TNFRSF14, TNFRSF19 and interferon regulatory factors IFI6 and IRF4 were highly expressed in MC compared to KB and FB (S2 Table).

In our previous study [12] a quantitative real-time PCR (QRT-PCR) analysis based predesigned TaqMan Gene Expression Assays for selected interleukin 10 (IL10) family cytokine's genes showed the differential gene expression in cultivated MC relative to KC and FB. The data correlates well with the results obtained in the present study. For example IL20RA and IL20RB had significantly higher expression in KC compared to MC. Also IL22RA1 could be found in KC and not in MC, whereas IL22RA2 gene was expressed in whole skin but not in MC, KC or FB (S1 Table). Among the studied IL10 family cytokines (IL10, IL19, IL20, IL22, IL24, IL26, IL28B, IL29) and their receptors (IL10RA, IL10RB, IL20RA, IL20RB, IL22RA1, IL22RA2, IL28RA), IL24 was the most prominent cytokine in MC, which was hardly detectable in KC and FB (S2 Table). IL10 family of cytokines are responsible for host defense mechanisms and have both have both pro-inflammatory and anti-inflammatory roles [64].

Chronic inflammation is strongly connected to oxidative stress processes [65]. Melanin biosynthesis itself generates a large amount of free radicals [66], therefore it is crucial to have an efficient control system, which can balance the inflammatory process before it damages DNA or destroys the cell. For example, we found FOXO3 transcription factor, which coordinates reduction/oxidation balance in neural stem cells [67] and ATM, which assists cells in recognizing damaged, but is also modulating the antioxidant system, and glutathione peroxidase genes GPX3, GPX7 and GPX8, which regulate intracellular reactive oxygen species balance [68], to be upregulated in MC compared to KC and FB (S2 Table).

Stem cells-like properties.

Cultured MC expressed a wide range of genes, characteristic for stem cells. Evidence shows that several pathways that are important in normal stem cells (BCL2 family genes, Notch, Sonic hedgehog and Wnt signaling pathways), may also act in cancer development [69], [70]. For instance, we saw the highest level of antiapoptotic BCL2 and BCL2A1 and stem cell factor inducer [71] RCAN1 expression level in MC when compared to KC, FB and the whole skin (S2 Table). Interestingly, CD200 was uniquely expressed in MC (Table 2). CD200 has been proposed to be a follicular stem cell marker, but its expression increases also with apoptosis and cancers overexpressing CD200 expand and metastasize more rapidly [72]. And, as mentioned above, the expression of S100 calcium-binding proteins, which are specifically expressed in cells with stem cell properties was increased in MC (S2 Table). Tumor cells and stem cells have similarities in their self-renewal process; they have extensive proliferative potential and stem cells are often a target for malignant genetic transformations [69]. These stem-cell-like properties have brought forth MC as a potential source for induced pluripotent stem cells (iPSCs) [73].

Pathway analysis of melanocytes

Pathway analysis of differentially expressed genes also described MC as active and intensively dividing cell population in cell culture (Fig. 2). We identified a number of pathways prevalent in MC, which were characteristic for ongoing regenerative process and could be related to cell dividing processes: genes regulating mitotic activity and cell cycle, DNA replication and packing, assembling and metabolism of different cellular components (cytoskeleton, structural macromolecules), formation of lysosome and Golgi complex etc. (Fig. 2). Consequently, gene expression profile corresponding to enhanced metabolic activity could be also observed.

[Figure omitted. See PDF.]

Figure 2. Comparative pathway analysis of MC and KC, FB and whole skin.

Red plots indicate pathways, which were prominently expressed in MC. Blue plots mark pathways, which were downregulated in MC (A, B, C) and concomitantly upregulated in the whole skin (A), KC (B) and FB (C), respectively. The word sizes in wordclouds are defined by the p-values.

https://doi.org/10.1371/journal.pone.0115717.g002

Although the gene expression pattern analysis describes the cultured MC as highly proliferative cells we have to consider that the high proliferation rate is characteristic only to MC's cell in culture and does not reflect the actual situation in vivo. The cell culture conditions include less cell-cell contacts and a high level of growth factors reminiscent of active regenerative state (like in case of wound healing). Thus the analysis of gene expression pattern of cultured cells does not reflect the homeostatic state of the cells in a tissue but rather are expected to describe their response to the injury. On the other hand, a number of well-characterized cell type-specific genes could be readily identified from each cell population analyzed suggesting that the low-passage cultured cells have well retained their identity. Since the cell culture has remained as a gold standard in order to obtain sufficient amount of relatively homogenous cell populations for tissue engineering and toxicity testing, knowing the characteristic properties of cells in culture is instrumental for their further use in ex-vivo applications.

Conclusions

In this study we have identified a number of genes and pathways, which are characteristic or unique for MC compared to KC and FB. We also demonstrated the difference between gene expression pattern of MC culture and the whole skin. The data presented provide an insight into the various possible roles of MC in the skin. As expected by the rapid growth in the cell culture, our differential gene expression and pathway analyses described MC as cells with a high proliferative capacity in vitro compared to KC an FB. That might suggest they have preserved the readiness to regenerate and some stem-cells-like properties more than KC and FB. However, these properties make MC the most vulnerable cells in the skin and provide an explanation to their increased susceptibility to harmful environment agents (eg UV exposure) and high incidence rate of malignant melanoma. On the other hand, the increased stem cell-like properties might give MC a good self-renewing potential and also advocate for their use as a potential source for induced pluripotent stem cells for therapeutic purposes.

Supporting Information

[Figure omitted. See PDF.]

S1 Fig.

Multidimensional scaling blot for visualizing the level of similarity of individual samples in KC, MC, FB and whole skin (WS) groups. The function plotMDS in edgeR package evaluates the similarity of MC, KC, FB and whole skin replicate samples that were used in our experiment. It calculates the root-mean-square of top 500 genes with largest absolute log2 fold change between the two samples, termed leading log2-fold-change. From the MDS plot it is quite clear that the distinction between the MC, KC and FB is prominent over the distinction of samples retrieved from the same individual.

https://doi.org/10.1371/journal.pone.0115717.s001

(PDF)

S1 Table.

RPKM values of genes we detected in MC, KC, FB and the whole skin. A gene was considered as expressed only if the RPKM value in all samples of the corresponding cell type was above the cutoff (0.95). Conversely, the gene was labeled as not expressed if the RPKM value was below 0.95 in at least one of the samples.

https://doi.org/10.1371/journal.pone.0115717.s002

(XLSX)

S2 Table.

List of differentially expressed genes comparing MC to KC, FB and whole skin.

https://doi.org/10.1371/journal.pone.0115717.s003

(XLSX)

S3 Table.

List of pathways (based on differentially expressed genes), prominent in MC compared to KC, FB and whole skin.

https://doi.org/10.1371/journal.pone.0115717.s004

(XLSX)

Author Contributions

Conceived and designed the experiments: PR KK SK. Performed the experiments: PR ER. Analyzed the data: SI PR. Contributed reagents/materials/analysis tools: OP EV SK KK HS. Wrote the paper: PR VJ SI.

Citation: Reemann P, Reimann E, Ilmjärv S, Porosaar O, Silm H, Jaks V, et al. (2014) Melanocytes in the Skin – Comparative Whole Transcriptome Analysis of Main Skin Cell Types. PLoS ONE 9(12): e115717. https://doi.org/10.1371/journal.pone.0115717

References

1. Hoath SB, Leahy DG (2003) The organization of human epidermis: functional epidermal units and phi proportionality. J Invest Dermatol 121:1440–1446.

2. Miller CC, Godeau G, Lebreton-DeCoster C, Desmouliere A, Pellat B, et al. (2003) Validation of a morphometric method for evaluating fibroblast numbers in normal and pathologic tissues. Exp Dermatol 12:403–411.

3. Randolph RK, Simon M (1998) Dermal fibroblasts actively metabolize retinoic acid but not retinol. J Invest Dermatol 111:478–484.

4. Thingnes J, Lavelle TJ, Hovig E, Omholt SW (2012) Understanding the Melanocyte Distribution in Human Epidermis: An Agent-Based Computational Model Approach. PLoS One 7.

5. Lin JY, Fisher DE (2007) Melanocyte biology and skin pigmentation. Nature 445:843–850.

6. Chakraborty AK, Funasaka Y, Slominski A, Ermak G, Hwang J, et al. (1996) Production and release of proopiomelanocortin (POMC) derived peptides by human melanocytes and keratinocytes in culture: regulation by ultraviolet B. Biochim Biophys Acta 1313:130–138.

7. Funasaka Y, Chakraborty AK, Hayashi Y, Komoto M, Ohashi A, et al. (1998) Modulation of melanocyte-stimulating hormone receptor expression on normal human melanocytes: evidence for a regulatory role of ultraviolet B, interleukin-1alpha, interleukin-1beta, endothelin-1 and tumour necrosis factor-alpha. Br J Dermatol 139:216–224.

8. Costin GE, Hearing VJ (2007) Human skin pigmentation: melanocytes modulate skin color in response to stress. FASEB J 21:976–994.

9. Slominski A (2009) Neuroendocrine activity of the melanocyte. Exp Dermatol 18:760–763.

10. Goldsmith LA, Fitzpatrick TB, Stephen Katz S, Gilchrest B, Paller A, et al. (2012) Fitzpatrick's dermatology in general medicine. New York: McGraw-Hill Professional. p. p.

11. Lee JS, Kim DH, Choi DK, Kim CD, Ahn GB, et al. (2013) Comparison of Gene Expression Profiles between Keratinocytes, Melanocytes and Fibroblasts. Ann Dermatol 25:36–45.

12. Reemann P, Reimann E, Suutre S, Paavo M, Loite U, et al. (2013) Expression of Class II Cytokine Genes in Children's Skin. Acta Derm Venereol.

13. Hoek K, Rimm DL, Williams KR, Zhao H, Ariyan S, et al. (2004) Expression profiling reveals novel pathways in the transformation of melanocytes to melanomas. Cancer Res 64:5270–5282.

14. Kim D, Pertea G, Trapnell C, Pimentel H, Kelley R, et al. (2013) TopHat2: accurate alignment of transcriptomes in the presence of insertions, deletions and gene fusions. Genome Biol 14:R36.

15. Langmead B, Trapnell C, Pop M, Salzberg SL (2009) Ultrafast and memory-efficient alignment of short DNA sequences to the human genome. Genome Biol 10:R25.

16. Roberts A, Trapnell C, Donaghey J, Rinn JL, Pachter L (2011) Improving RNA-Seq expression estimates by correcting for fragment bias. Genome Biol 12:R22.

17. Robinson MD, McCarthy DJ, Smyth GK (2010) edgeR: a Bioconductor package for differential expression analysis of digital gene expression data. Bioinformatics 26:139–140.

18. Ramskold D, Wang ET, Burge CB, Sandberg R (2009) An abundance of ubiquitously expressed genes revealed by tissue transcriptome sequence data. PLoS Comput Biol 5:e1000598.

19. Reimand J, Arak T, Vilo J (2011) g:Profiler-a web server for functional interpretation of gene lists (2011 update). Nucleic Acids Research 39:W307–W315.

20. Cichorek M, Wachulska M, Stasiewicz A, Tyminska A (2013) Skin melanocytes: biology and development. Postepy Dermatol Alergol 30:30–41.

21. Hou L, Pavan WJ (2008) Transcriptional and signaling regulation in neural crest stem cell-derived melanocyte development: do all roads lead to Mitf? Cell Res 18:1163–1176.

22. Kawakami A, Fisher DE (2011) Key discoveries in melanocyte development. J Invest Dermatol 131:E2–4.

23. Roberts DW, Newton RA, Beaumont KA, Helen Leonard J, Sturm RA (2006) Quantitative analysis of MC1R gene expression in human skin cell cultures. Pigment Cell Res 19:76–89.

24. Liu W, Peng Y, Tobin DJ (2013) A new 12-gene diagnostic biomarker signature of melanoma revealed by integrated microarray analysis. PeerJ 1:e49.

25. Koga Y, Pelizzola M, Cheng E, Krauthammer M, Sznol M, et al. (2009) Genome-wide screen of promoter methylation identifies novel markers in melanoma. Genome Research 19:1462–1470.

26. Du J, Miller AJ, Widlund HR, Horstmann MA, Ramaswamy S, et al. (2003) MLANA/MART1 and SILV/PMEL17/GP100 are transcriptionally regulated by MITF in melanocytes and melanoma. Am J Pathol 163:333–343.

27. Lee ST, Nicholls RD, Jong MT, Fukai K, Spritz RA (1995) Organization and sequence of the human P gene and identification of a new family of transport proteins. Genomics 26:354–363.

28. Gillis JS (2006) Microarray evidence of glutaminyl cyclase gene expression in melanoma: implications for tumor antigen specific immunotherapy. J Transl Med 4:27.

29. Auge JM, Molina R, Filella X, Bosch E, Gonzalez Cao M, et al. (2005) S-100beta and MIA in advanced melanoma in relation to prognostic factors. Anticancer Res 25:1779–1782.

30. de Souza CF, Xander P, Monteiro AC, Silva AG, da Silva DC, et al. (2012) Mining gene expression signature for the detection of pre-malignant melanocytes and early melanomas with risk for metastasis. PLoS One 7:e44800.

31. Zentner GE, Henikoff S (2013) Regulation of nucleosome dynamics by histone modifications. Nat Struct Mol Biol 20:259–266.

32. Verhey KJ, Kaul N, Soppina V (2011) Kinesin Assembly and Movement in Cells. Annual Review of Biophysics, Vol 40 40:267–288.

33. Chen ZX, Wallis K, Fell SM, Sobrado VR, Hemmer MC, et al. (2014) RNA Helicase A Is a Downstream Mediator of KIF1Bbeta Tumor-Suppressor Function in Neuroblastoma. Cancer Discov 4:434–451.

34. Itsui Y, Sakamoto N, Kurosaki M, Kanazawa N, Tanabe Y, et al. (2006) Expressional screening of interferon-stimulated genes for antiviral activity against hepatitis C virus replication. Journal of Viral Hepatitis 13:690–700.

35. Freinkel RK, Woodley D (2001) The biology of the skin. New York: Parthenon Pub. Group. 432 p. p.

36. Kos R, Reedy MV, Johnson RL, Erickson CA (2001) The winged-helix transcription factor FoxD3 is important for establishing the neural crest lineage and repressing melanogenesis in avian embryos. Development 128:1467–1479.

37. Chakrabarty K, Von Oerthel L, Hellemons A, Clotman F, Espana A, et al. (2012) Genome wide expression profiling of the mesodiencephalic region identifies novel factors involved in early and late dopaminergic development. Biol Open 1:693–704.

38. Jiao Y, Yan J, Zhao Y, Donahue LR, Beamer WG, et al. (2005) Carbonic anhydrase-related protein VIII deficiency is associated with a distinctive lifelong gait disorder in waddles mice. Genetics 171:1239–1246.

39. Chen YK, Hsueh YP (2012) Cortactin-binding protein 2 modulates the mobility of cortactin and regulates dendritic spine formation and maintenance. J Neurosci 32:1043–1055.

40. Akaneya Y, Sohya K, Kitamura A, Kimura F, Washburn C, et al. (2010) Ephrin-A5 and EphA5 interaction induces synaptogenesis during early hippocampal development. PLoS One 5:e12486.

41. Hill KK, Bedian V, Juang JL, Hoffmann FM (1995) Genetic Interactions between the Drosophila Abelson (Abl) Tyrosine Kinase and Failed Axon Connections (Fax), a Novel Protein in Axon Bundles. Genetics 141:595–606.

42. Szatanik M, Vibert N, Vassias I, Guenet JL, Eugene D, et al. (2008) Behavioral effects of a deletion in Kcnn2, the gene encoding the SK2 subunit of small-conductance Ca2+-activated K+ channels. Neurogenetics 9:237–248.

43. Hotta K, Hosaka M, Tanabe A, Takeuchi T (2009) Secretogranin II binds to secretogranin III and forms secretory granules with orexin, neuropeptide Y, and POMC. Journal of Endocrinology 202:111–121.

44. Yim YS, Kwon Y, Nam J, Yoon HI, Lee K, et al. (2013) Slitrks control excitatory and inhibitory synapse formation with LAR receptor protein tyrosine phosphatases. Proc Natl Acad Sci U S A 110:4057–4062.

45. Louro R, Nakaya HI, Paquola AC, Martins EA, da Silva AM, et al. (2004) RASL11A, member of a novel small monomeric GTPase gene family, is down-regulated in prostate tumors. Biochem Biophys Res Commun 316:618–627.

46. Zhang X, Zhang Q, Zhang J, Qiu L, Yan S-S, et al. (2010) FATS is a transcriptional target of p53 and associated with antitumor activity. Molecular Cancer 9:244.

47. Ahmed IA, Pusch CM, Hamed T, Rashad H, Idris A, et al. (2010) Epigenetic alterations by methylation of RASSF1A and DAPK1 promoter sequences in mammary carcinoma detected in extracellular tumor DNA. Cancer Genet Cytogenet 199:96–100.

48. Foltz G, Ryu GY, Yoon JG, Nelson T, Fahey J, et al. (2006) Genome-wide analysis of epigenetic silencing identifies BEX1 and BEX2 as candidate tumor suppressor genes in malignant glioma. Cancer Res 66:6665–6674.

49. Chen J, Lui WO, Vos MD, Clark GJ, Takahashi M, et al. (2003) The t(1;3) breakpoint-spanning genes LSAMP and NORE1 are involved in clear cell renal cell carcinomas. Cancer Cell 4:405–413.

50. Marzook H, Li DQ, Nair VS, Mudvari P, Reddy SD, et al. (2012) Metastasis-associated protein 1 drives tumor cell migration and invasion through transcriptional repression of RING finger protein 144A. J Biol Chem 287:5615–5626.

51. Stark JL, Mehla K, Chaika N, Acton TB, Xiao R, et al. (2014) Structure and Function of Human DnaJ Homologue Subfamily A Member 1 (DNAJA1) and Its Relationship to Pancreatic Cancer. Biochemistry 53:1360–1372.

52. Bair EL, Chen ML, McDaniel K, Sekiguchi K, Cress AE, et al. (2005) Membrane type 1 matrix metalloprotease cleaves laminin-10 and promotes prostate cancer cell migration. Neoplasia 7:380–389.

53. Ekblom P, Lonai P, Talts JF (2003) Expression and biological role of laminin-1. Matrix Biol 22:35–47.

54. Kumari K, Gulati S, Smith DF, Gulati U, Cummings RD, et al. (2007) Receptor binding specificity of recent human H3N2 influenza viruses. Virol J 4:42.

55. Videira PA, Correia M, Malagolini N, Crespo HJ, Ligeiro D, et al. (2009) ST3Gal.I sialyltransferase relevance in bladder cancer tissues and cell lines. BMC Cancer 9:357.

56. Wu S, Singh S, Varney ML, Kindle S, Singh RK (2012) Modulation of CXCL-8 expression in human melanoma cells regulates tumor growth, angiogenesis, invasion, and metastasis. Cancer Med 1:306–317.

57. Neufeld G, Sabag AD, Rabinovicz N, Kessler O (2012) Semaphorins in angiogenesis and tumor progression. Cold Spring Harb Perspect Med 2:a006718.

58. Garcia-Areas R, Libreros S, Amat S, Keating P, Carrio R, et al. (2014) Semaphorin7A promotes tumor growth and exerts a pro-angiogenic effect in macrophages of mammary tumor-bearing mice. Front Physiol 5:17.

59. Fang C, Wen G, Zhang L, Lin L, Moore A, et al. (2013) An important role of matrix metalloproteinase-8 in angiogenesis in vitro and in vivo. Cardiovasc Res 99:146–155.

60. Padmanabhan R, Chen KG, Gillet JP, Handley M, Mallon BS, et al. (2012) Regulation and expression of the ATP-binding cassette transporter ABCG2 in human embryonic stem cells. Stem Cells 30:2175–2187.

61. Elinav E, Nowarski R, Thaiss CA, Hu B, Jin C, et al. (2013) Inflammation-induced cancer: crosstalk between tumours, immune cells and microorganisms. Nat Rev Cancer 13:759–771.

62. Wajant H (2009) The role of TNF in cancer. Results Probl Cell Differ 49:1–15.

63. Cantaert T, Baeten D, Tak PP, van Baarsen LG (2010) Type I IFN and TNFalpha cross-regulation in immune-mediated inflammatory disease: basic concepts and clinical relevance. Arthritis Res Ther 12:219.

64. Cao S, Zhang X, Edwards JP, Mosser DM (2006) NF-kappaB1 (p50) homodimers differentially regulate pro- and anti-inflammatory cytokines in macrophages. J Biol Chem 281:26041–26050.

65. Reuter S, Gupta SC, Chaturvedi MM, Aggarwal BB (2010) Oxidative stress, inflammation, and cancer: how are they linked? Free Radic Biol Med 49:1603–1616.

66. Jimenez-Cervantes C, Martinez-Esparza M, Perez C, Daum N, Solano F, et al. (2001) Inhibition of melanogenesis in response to oxidative stress: transient downregulation of melanocyte differentiation markers and possible involvement of microphthalmia transcription factor. J Cell Sci 114:2335–2344.

67. Yeo H, Lyssiotis CA, Zhang YQ, Ying HQ, Asara JM, et al. (2013) FoxO3 coordinates metabolic pathways to maintain redox balance in neural stem cells. Embo Journal 32:2589–2602.

68. Wang K, Zhang T, Dong Q, Nice EC, Huang C, et al. (2013) Redox homeostasis: the linchpin in stem cell self-renewal and differentiation. Cell Death Dis 4:e537.

69. Reya T, Morrison SJ, Clarke MF, Weissman IL (2001) Stem cells, cancer, and cancer stem cells. Nature 414:105–111.

70. White RM, Zon LI (2008) Melanocytes in development, regeneration, and cancer. Cell Stem Cell 3:242–252.

71. Wu Z, Li Y, MacNeil AJ, Junkins RD, Berman JN, et al. Calcineurin-Rcan1 interaction contributes to stem cell factor-mediated mast cell activation. J Immunol 191:5885–5894.

72. Rosenblum MD, Olasz E, Woodliff JE, Johnson BD, Konkol MC, et al. (2004) CD200 is a novel p53-target gene involved in apoptosis-associated immune tolerance. Blood 103:2691–2698.

73. Utikal J, Maherali N, Kulalert W, Hochedlinger K (2009) Sox2 is dispensable for the reprogramming of melanocytes and melanoma cells into induced pluripotent stem cells. J Cell Sci 122:3502–3510.