Comparison of Skin Transcriptome between Responder and Non-Responder Vitiligo Lesions to Cell Transplantation: A Clinical Trial Study.

Objective: Autologous transplantation of epidermal cells has been used increasingly to treat vitiligo patients and is a simple, safe, and relatively efficient method. However, the outcome is not always satisfactory, and some patients show less or no response to this treatment. This study was evaluated to identify genes expressed differently among responders and non-responders to cell transplantation to find potential markers that could predict 'patients' responses to this type of cell therapy. Materials and
Methods: Eleven stable vitiligo patients who received autologous epidermal cell transplantation were included in this clinical trial study. Before cell transplantation, skin samples were obtained from the recipient's vitiligo lesions. After epidermal cell transplantation, patients were followed for at least six months to assess the response to epidermal cell injection. RNA sequencing was used to determine potential gene expression profile differences between responder and non-responder vitiligo patients.
Results: The RNA sequencing results showed differences in expression levels of 470 genes between the skin specimens of responder versus non-responder patients. There were 269 up-regulated genes and 201 down-regulated genes. Upregulated genes were involved in processes, such as Fatty Acid Omega Oxidation. Down-regulated genes were related to PPAR signaling pathway, and estrogen signaling pathway. Among the most differentially expressed genes (DEGs) with the most altered RNA expression levels in responders versus non-responder patients, we selected three genes (up-regulated genes KRTAP10-11 and down-regulated genes IP6K2 and C9) as potential biomarkers, which are involved in associated pathways.
Conclusion: Based on our findings, it is estimated that proposed genes might predict the response of vitiligo patients to cell therapy. However, further studies are required to clarify the role of these genes in pathogenesis and to characterize gene expression in a larger number of vitiligo patients in the context of epidermal cell transplantation therapy (registration number: IRCT201508201031N16). Copyright© by Royan Institute. All rights reserved.

Introduction

Vitiligo is a relatively common skin disease that appears as hypopigmented or depigmented skin lesions affecting 0.5- 1% of the population worldwide. Although vitiligo is not life-threatening, it can lead to social avoidance and reduce the ‘patients’ quality of life (1, 2). Various surgical and nonsurgical therapies are currently available to treat this disease. Surgical methods such as skin grafting are not suitable for large skin lesions and may not always be applicable for areas such as the lips, eyelids, and genital areas which are commonly affected. Complications include heterogeneity of texture and color with the surrounding recipient skin and infection of the graft site (3).

On the other hand, medical treatments are often less effective than surgical procedures and are associated with relatively high recurrence rates following their discontinuation (4). Phototherapy is another treatment option. However, there is concern on whether ultraviolet (UV) exposure predisposes patients to skin cancer (5-7), but current scientific literature shows conflicting results (8). In recent decades, cell therapy has been introduced as a novel therapeutic choice for patients with stable vitiligo, particularly those resistant to other available treatments (9). However, response to treatment is not seen in all patients and is not predictable. For instance, in Orouji et al. (10) study, a total of 1060 patches in 300 stable vitiligo patients were treated with intra-lesional epidermal cell suspension; however, 22.3% of the patients did not respond to the treatment after 30 months of follow-up.

Considering the risk of local adverse effects on the donor site, such as hypopigmentation and scar formation, the treatment cost, and the psychological problems resulting from poor treatment outcomes, it is necessary to find markers and develop reliable methods to predict the response rate to cell transplantation in vitiligo patients.

Here, we compared the skin transcriptome between responder and non-responder patients treated with epidermal cell suspension to find molecular biomarkers that could potentially predict response to cell transplantation.

Materials and Methods

Case selection

Vitiligo patients were recruited with informed consent from the Dermatology Clinic at Royan Institute between February 2017 and September 2017. Inclusion criteria included focal or generalized vitiligo patients aged 18-50 years old, with clinical stability of at least one year before recruitment to the study, and at least one vitiligo lesion on the trunk or limbs. Pregnant or lactating women, patients with active infectious diseases or undergoing any kind of cell therapy, and patients with a history of UV or laser therapy, or previous treatment with immunosuppressive or cytotoxic medications within six months before entry to the study, were excluded. Vitiligo lesions on the trunk or limbs were selected to obtain biopsy specimens, and none of the biopsies were taken from facial lesions due to aesthetic concerns.

This clinical trial study was approved by the Ethics Committee of the Royan Institute. The trial was registered with the United States National Library of Medicine Clinical Trials.gov (NCT00631865). The clinical trial was also registered with the Iranian Registry of Clinical Trials (IRCT201508201031N16). The Research Ethics Committee of Royan Institute approved skin biopsies for RNA sequencing (IR.ACECR.ROYAN.REC.1396.196).

Skin sampling, cell preparation, and transplantation

Skin sampling, cellular suspension separation techniques, and cell transplantation techniques were adapted from Khodadadi et al. (11) study. Briefly, the skin sample was obtained from the thigh-buttock junction as a donation site, with about one-third to the one-seventh surface area. The skin pieces were incubated with dispase II until the epidermis separated from the dermis. Subsequently, epidermal cells were separated from each other by trypsin/EDTA. The cell suspension was injected intralesionally by a trained dermatologist one day after the skin biopsy.

Obtaining skin biopsy

Selected vitiligo lesions were disinfected and anesthetized by a local injection of 0.2-0.3 ml of 1% lidocaine by the dermatologist. Before cell transplantation, two 2.5 mm diameter punch biopsy samples were obtained from the vitiligo area. The biopsy specimens were immersed inside liquid nitrogen for 30 seconds, and then transferred to the -80˚C freezer until the RNA was extracted.

Assessment of repigmentation

The percentage of repigmentation in the treated lesions was assessed subjectively by a dermatologist before and at 2, 4, and 6 months after the treatment. Grades 0, I, II, III and IV were assigned to patches with 0, 1-24%, 25-49%, 50-74% and 75-100% of repigmentation. Repigmentation was assessed based on lesion color and size changes in photographs. Patches with repigmentation of ≥ 25% were considered as “responder” patches (10).

RNA extraction

Total RNA was extracted from skin biopsies stored at the -80˚C freezer using the TRIzol reagent (Invitrogen, 15596-018) and RNeasy micro kit Cat No. /ID: 74004. RNA quantity and quality were measured using Agilent Technologies 2100 Bioanalyzer (or 2200 TapeStation).

Whole transcriptome RNA sequencing and data analysis

Total RNA was extracted from 11 samples for library preparation, and the library was constructed using the TruSeq RNA Access Library Prep Kit. Sequencing was done using the NovaSeq Sequencer. The quality of produced data was determined using the FastQC V-0.11.8 software. The reads were filtered based on sequencing quality with the Trimmomatic V-0.36. Trimmed reads were then mapped to reference genome GRCh V-38 using HISAT2 V-2.1.0. Finally, the number of reads aligning to each gene was determined using htseq-count V-0.11.2 based on gene annotation file GTF V-91. Normalization and differential expression analysis were done using the DESeq2 package in R statistical software. We used P<0.05, 1.5 ≤ fold change ≤ -1.5. Pathway enrichment analysis of the DEGs was performed by Enrichr database (12). RNA sequencing data has been published in NCBI GEO, accession number PRJNA633437.

Results

Patients characteristics

Skin biopsies were obtained from 11 patients with stable vitiligo (5 females and 6 males) with a mean age of 28.63 years. None of the patients reported any family history of autoimmune diseases such as vitiligo and thyroid diseases.

Clinical outcomes

The follow-up visits were performed at 2, 4, and 6 months after cell transplantation. Six months after cell therapy, based on clinical assessment of treated patches by one dermatologist, 6 patients were considered responders to treatment with repigmentation of 25-49% (grade II) while 5 patients were considered unresponsive with repigmentation of 1-24% (grade 1) (Table 1). Figure 1 shows the extent of repigmentation changes in two responder and non-responder patients.

Table 1

Characteristics of patients

Patient ID	Age (Y)	Sex	Recipient site	Response	Recipient size (cm2)	Injected cells count /cm²
P15	25	Male	Back	R	5	69333
P18	26	Male	Elbow	R	30	85714
P20	32	Male	Hand	R	70	93333
P25	34	Female	Leg	R	5	63063
P26	21	Female	Back	R	9	220000
P27	32	Female	Hand	R	12	104545
P1	18	Female	Elbow	N	3	200000
P9	25	Female	Leg	N	80	34574
P10	26	Male	Thigh	N	80	121951
P17	28	Male	Leg	N	30	107142
P24	48	Male	Leg	N	70	136363

R; Responder and N; Non-responder

Fig 1

Vitiligo lesions. A. Responder and B. Non-responder. m; Months post-transplantation.

Transcriptome profiling and identification of biomarker(s)

We examined differentially expressed genes (DEGs) between responder (n=6) and non-responder (n=5) vitiligo patients. There were 470 DEGs, including 269 up-regulated and 201 down-regulated genes in responder versus non-responder patients (Fig .2A). Up-regulated genes were associated with processes, such as Fatty Acid Omega Oxidation (Fig .2B1). Down-regulated genes were involved in PPAR signaling pathway, and estrogen signaling pathway (Fig .2B2).

Fig 2

Differential expression of genes in responder versus non-responder patients. A. Hierarchical clustering of a heatmap for (DEGs). Different colors indicate up-regulated (red); down-regulated (blue) genes. B. Pathway enrichment analysis results by Enrichr. B1. Up-regulated; B2. Down-regulated genes in responder versus non-responder patients. R; Responder, N; Non-responder and P; Patient.

Ten up-regulated and ten down-regulated DEGs with the most altered RNA expression levels are presented in Tables 2 and 3. Among these genes, we selected three candidate genes, i.e, KRTAP10-11, IP6K2, and C9 that were involved in related pathways. KRTAP10-11 is an upregulated gene involved in keratinization, while IP6K2 is a down-regulated gene that contributes to the regulation of apoptotic processes, and C9 is another down-regulated gene that regulates C9 in the complement cascade.

Discussion

Biomarkers have been one of the most common tools for predicting prognosis or response to treatment in different diseases. RNA biomarkers have several advantages in contrast to other biomarkers such as DNA or protein, and they transmit genetic and regulatory information compared to DNA biomarkers. Compared to protein biomarkers, they reflect cellular states and are more sensitive and specific. Moreover, the cost of finding RNA biomarkers is much lower than protein because specific antibodies are required to detect each protein individually (13).

Transcriptomics technology is often used for biomarker discovery. Among transcriptome analysis methods, RNA sequencing assesses the quantified gene expression measurements at the whole transcriptome level. Moreover, RNA-sequencing provides novel transcript/ gene discovery, alternative splicing, gene fusion, and higher sensitivity in expression level (14, 15).

Numerous studies have been conducted on the potential of using RNA-sequencing technology to discover biomarkers for various biological conditions. In 2017, Wright et al. (16) identified a biomarker panel that can predict response to tumor necrosis factor inhibitors (TNFi) therapy in rheumatoid arthritis using RNA sequencing. A study in 2018 by Moreno-Torres et al. (17) conducted transcriptome profiling for response prediction in treating multiple sclerosis patients with fingolimod. They showed that evaluating differential gene expression before treatment might be useful as a biomarker.

No biomarker has been discovered to predict the response to autologous epidermal cells transplantation in vitiligo patients.

Sequencing results of our study showed that the expression of 470 genes (269 up-regulated, 201 downregulated) were different in responder vitiligo patients compared to non-responders.

Up-regulated genes in responders were involved in fatty acid omega oxidation. Fatty acids regulate tyrosinase synthesis and degradation, and thus regulate pigmentation (18).

As for pathways related to down-regulated genes in responders, PPARs belong to the nuclear hormone receptors subfamily (19). We found decreased PPARs signaling pathway responders’ expression and increased expression of non-responders. PPARα, as an E3 ubiquitin ligase induces Bcl-2 (anti-apoptotic protein) ubiquitination and leads to apoptotic death (20). A previous study has shown that levels of Bcl-2 expression are lower in lesional skin compared to normally pigmented skin in vitiligo patients (21). Genes related to the estrogen signaling pathway are also down-regulated in responders. Generation of H2 O2 by estrogens can cause DNA damage in peripheral blood lymphocytes of vitiligo patients (22). Studies have shown that there are increased levels of H2 O2 in the skin and blood cells of vitiligo patients, thus resulting in accumulation of free radical mediated melanocyte degeneration (23-25).

Also, according to most DEGs and associated terms, three genes were chosen. KRTAP10-11 was up-regulated, and IP6K2 and C9 were down-regulated in responder versus non-responder patients.

The sequencing results demonstrated increased expression level of KRTAP10-11 in responders compared to non-responder patients (FC=2.54, P<0.05). KRTAP10-11 is a gene that encodes keratin-associated protein 10-11 and is involved in keratinization (26). Keratin 10 is a differentiated keratinocyte marker and is expressed in the differentiating suprabasal layers of the epidermis (27). Thus, up-regulation of this gene in responder patients could show its role in keratinization and distribution of received melanin from melanocytes into the epidermis. IP6K2 was also found up-regulated in non-responders compared to responders (FC=-2.1, P<0.05). The role of inositol pyrophosphates (IPs) has been identified in many biological processes (28-30). In response to stress, IP6K2 binds to the tumor suppressor p53 and modulates cell death. p53 inhibits the expression of proarrest target genes or activates pro-apoptotic target genes. IP6K2 selectively modulates inhibiting the induction of p53 proarrest pathway, thus augmenting the p53 response to apoptosis (31). There is no study to date investigating the relation of IP6K2 to vitiligo disease. It seems that the increased expression level of this gene in non-responders may have a role in melanocyte apoptosis.

We showed an increase in C9 expression in nonresponders compared with responder patients (FC=-2, P<0.05). C9 is an essential member of the complement system membrane attack complex (MAC), which creates pores on the cell membrane of target pathogens resulting in their destruction (32). Furthermore, in vitiligo, antipigment cell antibodies can induce melanocyte damage by complement activation and (ADCC) (33). Although no study has previously reported the effects of C9 expression on vitiligo, the increased activity of this gene in nonresponder patients could be attributed to the role of this gene in melanocytes destruction.

Regarding the unknown role of these genes in the pathogenesis of vitiligo, further studies are required to reveal the underlying mechanisms of the gene expression differences.

Conclusion

We proposed candidate biomarkers KRTAP10-11, IP6K2, and C9 that potentially may be used as diagnostic tools in the pathogenesis and prediction of vitiligo patients responses to epidermal cell transplantation. These preliminary data are promising; however, further biomarker screening with larger sample sizes and analyses of candidate genes protein expression patterns are required to predict the response of vitiligo patients to cell transplantation.

Table 2

Ten up-regulated differentially expressed genes (DEGs) with the most altered RNA expression levels in responder patients versus non-responder patients

Gene symbol	Gene title	P value	Fold change	Gene ontology (molecular function)
ZNF486	Zinc finger protein 486	0.00002	2.85	DNA binding; metal ion binding
KRTAP10-13P	Keratin associated protein 10-13	0.0001	2.60	Pseudogene
KRTAP10-11	Keratin-associated protein 10-11	0.0001	2.54	Keratinization
ZDHHC20	Palmitoyltransferase ZDHHC20	0.0002	2.47	Palmitoyltransferase activity; protein-cysteine S-palmitoyltransferase activity; zinc ion binding
PGAM1P4	Phosphoglycerate mutase 1 pseudogene 4	0.0003	2.46	Unknown
AC002454.1		0.0006	2.34	Unknown
AC005550.3		0.0007	2.32	Unknown
RIOX1	Ribosomal oxygenase 1	0.0008	2.30	Histone demethylase activity (H3-K36 specific); histone demethylase activity (H3-K4 specific); iron ion binding; oxidoreductase activity, acting on paired donors, with incorporation or reduction of molecular oxygen, 2-oxoglutarate as one donor, and incorporation of one atom each of oxygen into both donors
C2CD6	C2 Calcium Dependent Domain Containing 6	0.0006	2.23	Unknown
GBX1	Homeobox protein GBX-1	0.0008	2.23	RNA polymerase II transcription factor activity, sequence-specific DNA binding; sequence-specific DNA binding

Table 3

Ten down-regulated differentially expressed genes (DEGs) with the most altered RNA expression levels in responder patients versus non-responder patients

Gene symbol	Gene title	P value	Fold change	Gene ontology (molecular function)
AL691447.3		0.0001	2.62	Unknown
AC022182.3		0.0006	2.30	Unknown
WAPL	Wings apart-like protein homolog	0.0003	2.25	Cell division; meiotic chromosome segregation; mitotic sister chromatid cohesion; negative regulation of chromatin binding; negative regulation of DNA replication; negative regulation of sister chromatid cohesion; positive regulation of fibroblast proliferation; protein localization to chromatin; regulation of chromosome condensation; regulation of cohesion loading; response to toxic substance; viral process
HOXA11	Homeobox protein Hox-A11	0.0006	2.20	RNA polymerase II transcription factor activity, sequence-specific DNA binding; sequence-specific DNA binding
IP6K2	Inositol hexakisphosphate kinase 2	0.0014	2.10	ATP binding; inositol-1,3,4,5,6-pentakisphosphate kinase activity; inositol 5-diphosphate pentakisphosphate 5-kinase activity; inositol diphosphate tetrakisphosphate kinase activity; inositol heptakisphosphate kinase activity; inositol hexakisphosphate 1-kinase activity; inositol hexakisphosphate 3-kinase activity; inositol hexakisphosphate 5-kinase activity; inositol hexakisphosphate kinase activity
ZNF652	Zinc finger protein 652	0.001	2.10	DNA binding; metal ion binding; RNA polymerase II transcription factor activity, sequence-specific DNA binding
MISP3	Uncharacterized protein MISP3	0.002	2.08	Unknown
ANKRD28	Ankyrin Repeat Domain 28	0.001	2.05	Unknown
C9	Complement C9	0.0006	2	Complement system
PCDHA3	Protocadherin alpha-3	0.0026	1.96	Calcium ion binding