Whole Exome Sequencing in Psoriasis Patients Contributes to Studies of Acitretin Treatment Difference

Psoriasis vulgaris is an immune-mediated inflammatory skin disease. Although acitretin is a widely used synthetic retinoid for moderate to severe psoriasis, little is known about patients’ genetics in response to this drug. In this study, 179 patients were enrolled in either the discovery set (13 patients) or replication set (166 patients). The discovery set was sequenced by whole exome sequencing and sequential validation was conducted in the replication set by MassArray assays. Four SNPs (single nucleotide polymorphisms) (rs1105223T>C in CRB2, rs11086065A>G in ANKLE1, rs3821414T>C in ARHGEF3, rs1802073 T>G in SFRP4) were found to be significantly associated with acitretin response in either co-dominant or dominant models via multivariable logistic regression analysis, while CRB2 rs1105223CC (OR = 4.10, 95% CI = 1.46–11.5, p = 0.007) and ANKLE1 rs11086065AG/GG (OR = 2.76, 95% CI = 1.42–5.37, p = 0.003) were associated with no response to acitretin after 8-week treatment. Meanwhile, ARHGEF3 rs3821414CT/CC (OR = 0.25, 95% CI = 0.10–0.68, p = 0.006) and SFRP4 rs1802073GG/GT (OR = 2.40, 95% CI, 1.23–4.70, p = 0.011) were associated with a higher response rate. Four new genetic variations with potential influences on the response to acitretin were found in this study which may serve as genetic markers for acitretin in psoriasis patients.


Introduction
Psoriasis is a common immunologically mediated inflammatory skin disease, characterized by abnormal T-cell activation and inadequate keratinocyte differentiation [1][2][3]. It affects 2%-3% of the world population, and its morbidity is still increasing in recent years [4]. The pathogenesis of psoriasis is complex and involves genetic, environmental, immunological, and even, neurologic factors [5][6][7]. Although great efforts have been made in elucidating the pathogenesis of the disease, the full mechanism is not completely understood.
Acitretin is a synthetic retinoid belonging to the family of retinoid analogs (RA) drugs, and is widely used in moderate to severe psoriasis patients. Acitretin was speculated to regulate the differentiation, proliferation and apoptosis of human epidermal keratinocytes. It is thought to function through interfering with the expression of epidermal growth factor genes [8]. Also, it is reported that acitretin exerted an important influence on Th1 and Th17 cells during the treatment of psoriasis vulgaris; it reduced Th1 and Th17 cell infiltration and attenuated their cytokines in the skin [9]. There is also evidence that acitretin has immunomodulatory properties by inhibiting dermal microvascular endothelial cells and neutrophil migration [10]. The responsiveness of acitretin in psoriasis is notoriously variable. As reported, the response rate of acitretin, defined as 75%, improved the psoriasis area severity index (PASI) (PASI 75) after 12 weeks of treatment by 46%-52% [11]. There is little research focused on the metabolic process of acitretin in vivo, and the mechanism of pharmacology in psoriasis is still unclear.
Few pharmacogenomic studies have focused on this agent so far. Polymorphisms of the apolipoprotein E gene (APOE) and the vascular endothelial growth factor gene (VEGF) have been evaluated as predictors of response to psoriasis patients treated with acitretin. The results revealed that ApoE protein variants did not have any utility as pharmacogenetic markers for predicting patients' response to acitretin [12], and the genetic variant (−460T>C) of VEGF was associated with the response to acitretin in psoriatic patients [13]. However, these studies only focused on isolated polymorphisms within single genes relevant to acitretin metabolism and considerable variation exists across all genes involved in acitretin metabolism. Utilizing whole exome sequencing, we have investigated predictors of outcome to acitretin therapy across each relevant gene in the largest patient cohort studied to date.

Clinical Features of the Psoriatic Patients in Difference Phase
Among the discovery and verification phases, the baseline age, gender and BMI (Body Mass Index) of the patients had no significant difference between the effective and ineffective groups (Table 1). There were also no significant differences with respect to BMI (23.96 ± 3.95 vs. 22.99 ± 3.75, p = 0.340), age (48 ± 16 vs. 42 ± 13, p = 0.233) and gender (p = 0.216) between the discovery and verification phases.

Whole Exome Sequencing Analysis
To reconcile the clinical findings with molecular data in psoriatic patients, a total of 13 patients were selected for whole exome sequencing analysis, and we obtained 38,190 variants in this study. After identification of all the variant calls, Fisher's exact test and the Cochran-Armitage trend as well as different genetic models (dominant, recessive and general) analyses were used to identify variants that were significantly associated with drug response (p < 0.05). In fact, 1790 variants were found to be associated with drug efficacy, the most significant SNPs rs2241984 (MaxSig. p = 9.04 × 10 −5 , Fisher. p = 1.82 × 10 −4 ) as shown in Figure 1; and the top 20 statistics of Pathway Enrichment were shown in Figure 2. Moreover, according to the literature reports, mutation location and significance, 34 positive SNPs were selected and validated by MassArray in independent samples.

Univariate Analysis of Thirty-Four Positive SNPs
A total of 166 psoriatic patients were recruited in the verification phase. For quality control, only SNPs with a frequency above 5% and with a genotyping rate ≥95% were included in the final statistical analysis and 34 SNPs were all qualified. All SNPs were agreed with the Hardy-Weinberg equilibrium except for rs2303022, rs2376558, rs47 and rs76310711 variations. We then analyzed associations between the 34 selected SNPs and drug efficacy ( Table 2). As shown in Table 2, we found that crumbs 2 (CRB2) rs1105223T>C, ankyrin repeat and LEM domain containing 1 (ANKLE1) rs11086065 A>G, Rho guanine nucleotide exchange factor 3 (ARHGEF3) rs3821414 T>C and secreted frizzled related protein 4 (SFRP4) rs1802073 G>T, were significantly associated with drug efficacy in either co-dominant or dominant models.

Multivariable Logistic Regression Analysis of CRB2 rs1105223T>C, ANKLE1 rs11086065 A>G, ARHGEF3 rs3821414 T>C and SFRP4 rs1802073 G>T
In order to further verify the effect of four SNPs, multivariable logistic regression analysis with adjustment for confounders including age, gender and body mass index (BMI) was used. For CRB2 rs1105223T>C variation, 63 patients carried the CRB2 rs1105223TT genotype, 64 patients carried the CRB2 rs1105223CT genotype, and 29 patients carried the CRB2 rs1105223CC genotype. CRB2 rs1105223CC was associated with the ineffective response (OR = 4.098, 95% CI = 1.461-11.493, p = 0.007) compared to the TT genotype, and CRB2 rs1105223TT/CT was also associated with the drug efficacy compared to the CC genotype (OR = 0.588, 95% CI, 0.363-0.955, p = 0.032).

Univariate Analysis of Thirty-Four Positive SNPs
A total of 166 psoriatic patients were recruited in the verification phase. For quality control, only SNPs with a frequency above 5% and with a genotyping rate ≥95% were included in the final statistical analysis and 34 SNPs were all qualified. All SNPs were agreed with the Hardy-Weinberg equilibrium except for rs2303022, rs2376558, rs47 and rs76310711 variations. We then analyzed associations between the 34 selected SNPs and drug efficacy ( Table 2). As shown in Table 2, we found that crumbs 2 (CRB2) rs1105223T>C, ankyrin repeat and LEM domain containing 1 (ANKLE1) rs11086065 A>G, Rho guanine nucleotide exchange factor 3 (ARHGEF3) rs3821414 T>C and secreted frizzled related protein 4 (SFRP4) rs1802073 G>T, were significantly associated with drug efficacy in either co-dominant or dominant models. In order to further verify the effect of four SNPs, multivariable logistic regression analysis with adjustment for confounders including age, gender and body mass index (BMI) was used. For CRB2 rs1105223T>C variation, 63 patients carried the CRB2 rs1105223TT genotype, 64 patients carried the CRB2 rs1105223CT genotype, and 29 patients carried the CRB2 rs1105223CC genotype. CRB2 rs1105223CC was associated with the ineffective response (OR = 4.098, 95% CI = 1.461-11.493, p = 0.007) compared to the TT genotype, and CRB2 rs1105223TT/CT was also associated with the drug efficacy compared to the CC genotype (OR = 0.588, 95% CI, 0.363-0.955, p = 0.032).

Discussion
In this study, we performed whole exome sequencing for 13 psoriasis patients, who experienced either high or extremely low efficacy. Although p values obtained at the discovery stage were individually rather weak, they provide us some biological knowledge for reference. To gain further insight into the potential influence of these genetic markers on the outcome of acitretin treatment in patients with psoriasis, we picked 34 SNPs for validation in an independent set of 166 patients by Sequenom MassArray. The 34 SNPs (as shown in the Table 2) were primarily chosen from the pathway enrichment ( Figure 2) and the results of the whole exome sequencing with p value <0.05. Several studies have previously proved the association between candidate genes involved in metabolic pathways of acitretin and the pathogenic mechanism of psoriasis, such as CSMD1, CCHCR1, GLI1, SFRP4 etc. [14][15][16][17]. In this study, we identified four SNPs that might be associated with the response to acitretin. The four genetic variants-rs1802073 G>T in SFRP4, rs1105223T>C in CRB2, rs3821414T>C in ARHGEF3, rs11086065A>G in ANKLE1-could be validated as predictive markers for the response to acitretin in psoriasis.
Genetic-associated studies identified dozens of psoriasis associated genes and the signaling pathways [18,19] include Notch signaling and Wnt signaling [20,21]. Notch signaling is associated with normally differentiated human epidermis, confirming its involvement in keratinocyte differentiation [22]. Wnt signaling participates in cell proliferation, adhesion and differentiation, suggesting that this pathway might be involved in psoriasis pathogenesis [23,24]. Acitretin is a retinoic acid derivative, although the mechanisms of acitretin to treat psoriasis is unclear, it is possible that the inhibition of Wnt signaling and the activation of Notch signaling [25].
The genetic variant rs1802073 in SFRP4 is a missense variant (a type of nonsynonymous substitution); rs1802073 at position g.7:37947164G>T results in a proline to threonine change at position 320 in the SFRP4 protein. SFRP4 is one of the secreted frizzled-related protein family members, thought to be a negative regulator of the Wnt signaling pathway [26,27]. Recently, it has been shown that the expression of SFRP4 was diminished in the lesional skin of patients with psoriasis [27]. SFRP4 directly inhibits the excessive keratinocyte proliferation evoked, and decreases the severity of the psoriasiform skin phenotype, including decreased acanthosis and reduced leukocyte infiltration [17]. Recently, Green et al indicated that SFRP4 could be a direct target gene of RARs, and RAR agonist induced the significant upregulation of SFRP4 [28]. Furthermore, acitretin is a member of the RAR agonist family, so we speculated that acitretin may alter the transcriptional regulation of SFRP4, such as encoding soluble Wnt antagonists and also ligands and receptors of the Notch pathway [25,29].
The rs1105223 in CRB2 is a missense variant; it is at position g.9:126128211T>C, resulting in a methionine to threonine shift at position 145 in the CRB2 protein. CRB2 is known to contain 15 extracellular EGF-like domains and three extracellular laminin G-like domains; it encodes the extracellular tenth EGF-like domain, and acts as an inhibitory binding protein to influence Notch signaling [30]. The rs11086065 in ANKLE1 is a missense variant at position G.19: 17284194A>G, resulting in a glutamine to arginine shift at position 452 in the ANKLE1 protein. Few references could be found about the function of ANKLE1, but because it contains a GIY-YIG-type (conserved N-terminal catalytic domains connected by linkers to C-terminal DNA-binding domains) nuclease domain, there is likely to be a potential role for ANKLE1 in DNA damage response [31]. A recent study presented a potential linkage of SNPs in the human ANKLE1 gene, showing an association between a function of ANKLE1 in multiple autoimmune syndromes [32] and the increased risk of certain cancers [33]. The rs3821414 in ARHGEF3 is a 3_prime_UTR_variant; this genetic variant does not result in amino acid change. One of the functions of ARHGEF3 is that it modulates differentiation through the activation of RhoA [34]. Although the three SNPs were associated with the ineffective/effective response to acitretin, we cannot conclude whether these SNPs influence the response to acitretin.
In our study, the genetic variants (rs1802073, rs1105223, rs11086065) introduce amino acid changes and may affect protein function, so we used bioinformatics approaches to analyze the effect of these SNPs on protein structure and function, such as PolyPhen, SNPeffect, SIFT and GTEx Protal.
No significant eQTLs were found for four SNPs (rs1802073, rs1105223, rs11086065, rs3821414) in tissue Whole Blood by using GTEx Protal. Furthermore, PolyPhen, SIFT and SNPeffect can predict the function of protein-coding: rs1802073T>G in SFRP4 was judged to be possibly damaging and tolerated; rs1105223T>C in CRB2 was judged to be benign and tolerated; rs11086065A>G in ANKLE1 was judged to be benign and deleterious, respectively. Therefore, it is generally believed that these SNPs may influence the response to acitretin, and further study is needed on the mechanism.

Patients
This study was approved by the Ethic Committee of XiangYa Hospital, the registration number of Chinese Clinical Trial Registry online is ChiCTR-OCH-14004518, and the registration number of ClinicalTrials.gov Protocol Registration and Results System (PRS) is NCT02715960. From April 2014 to July 2016, a total of 179 psoriatic patients (13 patients in the discovery set and 166 for the replication set) with moderate to severe psoriasis were recruited in the Department of Dermatology, Xiang Ya Hospital, Central South University. In the discovery set, five patients were defined as the responders with 100% improvement of PASI from the baseline post and eight patients were non-responders with −20% to −200% improvement of PASI (mean more serious) from the baseline post at 12 weeks of treatment. In the replication set, 100 patients (60.2%) were responders who achieved PASI75 at week 12, and 66 patients (39.8%) were non-responders who failed PASI75. The demographic, photography data as well as clinical data of psoriatic patients at each visit were collected. All patients were treated by a combination of 30 mg/day acitretin (Huapont Pharm., Chongqing, China) and a topical cream, calcipotriol (Bright Future Pharmaceutical Laboratories Ltd., Hongkong) for 12 weeks before the first follow-up.
The psoriatic patients who participated in this study were from the Southern Han Chinese population with the inclusion criteria as follows: (i) patients fulfilled the diagnostic criteria for psoriasis vulgaris; (ii) patients were aged ≥18 years; (iii) no medicines were received in the last four weeks before administration of the study agent. The exclusion criteria were as follows: (i) had other dermatological problems or any other diseases, and required pharmacological treatment; (ii) pregnant women, breastfeeding mothers, or women who were trying to become pregnant within the next 3 years. The study was approved by the institutional review board of Xiang Ya Hospital, Changsha, China. The written informed consent was obtained from each patient before participating in this study. The severity of psoriasis and the drug efficacy were both assessed by PASI. Patients with a PASI score greater than 10 are defined as moderate to severe psoriasis [35]. In conformity with guidelines, PASI75 at week 12 after the initiation of treatment was adopted as the index of response to acitretin in this study, which corresponds to a 75% improvement from baseline. Patients who reached a reduction of less than 75% from the baseline PASI at week 12 were considered as nonresponders.

DNA Extraction
Genomic DNA was extracted from the venous blood samples (5 mL) using a FlexiGene DNA Kit for mammalian blood (Qiagen, Hilden, Germany) according to manufacturer's instructions. Purity/concentration was determined using a Bio-spec Nano Spectrophotometer (Shimadzu Corporation, Japan). All the blood samples were stored at −80 • C until used.

Whole Exome Sequencing
Each captured library was then loaded on the Illumina Hiseq2000 platform (Illumina, San Diego, CA, USA). Whole exome sequencing was performed using Illumina base-calling Software 1.7 for base calling with default parameters. Briefly, DNA was cut into as well as purified 200-300-bp fragments. DNA was then amplified by PCR; clusters of PCR colonies were then sequenced on the Illumina Hiseq2000 platform (Illumina, USA); the sequence of each individual was generated from paired-end 90-bp reads.
The raw sequence data were aligned to the GRCh37 human reference genome using Burrows-Wheeler Aligner (BWA v0.7.7-r411) [36]. PCR duplicates were marked using the Mark Duplicates program in Picard-tools-1.115 tool set. GATK v3.2-2 and Samtools were used for the identification of INDEL, base quality score recalibration (BQSR) and the SNVs (single nucleotide variants) respectively. All variants were annotated using the Annovar program. The Mapping Quality Rank Sum Test (MQRankSum, u-based z-approximation from the Mann-Whitney Rank Sum Test for mapping qualities, only for heterozygous calls) and Read Pos Rank Sum (u-based z-approximation score from the Mann-Whitney Rank Sum Test for the distance from the end of the read for reads with the alternate allele, only for heterozygous calls) were conducted. The basic association tests implemented are the Cochran-Armitage trend test, Fisher's exact test, and different genetic models (dominant, recessive and general). Fisher's exact test was used to perform a standard case/control association analysis to generate significance. Those significance mutations are selected for annotation when the p value <0.05.

Sequenom MassArray Analysis
Thirty-four SNPs (single nucleotide polymorphisms) were analyzed using Sequenom MassArray (Sequenom, San Diego, CA, USA) according to the standard protocol recommended by the manufacturer [37].

Data Statistics and Analysis
The entire analysis was performed in the SPSS 23.0 statistical package (IBM SPSS, Chicago, IL, USA). The allele frequencies in different subgroups were tested by the chi-square analysis method. The chi-square test was used to determine whether genotype distribution of the SNPs agreed with Hardy-Weinberg equilibrium and to compare the distribution of categorical variables between drug response groups. The Fisher's exact test was used when data were spare. Comparisons of continuous variables between or among genotype groups were performed using nonparametric tests (Mann-Whitney U test). Adjusted odds ratio (ORs) and 95% confidence intervals (CIs) were used to describe drug outcomes. OR and 95% CI were calculated by limited backward-LR (likelihood ratio) logistic regression analysis with adjustment by clinical variables. A two-tailed p value less than 0.05 was regarded to be significant. The post hoc power of the sample size in χ-square analysis was operated with G. power (versions 3.1.9.2): the effect size was 1.024 and 0.417; α error was 0.05, and df was 5; and the sample size was 45 and 105 in the discovery and verification phases respectively, thus the power presented from 0.999 to 1.000.

Conclusions
In conclusion, four SNPs (rs1802073T>G in SFRP4, rs1105223T>C in CRB2, rs11086065A>G in ANKLE1, rs3821414T>C in ARHGEF3) were found to be associated with acitretin response via whole exome sequence and sequential validation, and there is accumulating evidence regarding the functional effects of these SNPs, especially rs1802073T>G in SFRP4. In the future, with additional work and validation, these variants will prove to be useful as markers for targeting therapies; they will be applicable more precisely and safely to individual patients, to optimize the treatment of psoriasis and minimize unnecessary expenditure. experiments; Yijing He contributed to the content of the manuscript; Xiang Chen designed the experiments, carried out and summarized data analysis, contributed to the text of the manuscript.

Conflicts of Interest:
The authors declare no conflict of interest.