Analysis of receptor tyrosine kinase genetics identifies two novel risk loci in GAS6 and PROS1 in Behçet’s disease

The TAM kinase (Tyro3, Axl, Mer) and its two ligands (Gas6 and protein S) have been shown to play an important regulatory role in the innate immune response. The present study aimed to investigate whether the tag single-nucleotide polymorphisms (tag SNPs) of these 5 protein-coding genes are associated with Behçet’s disease (BD). A two-stage association study was performed in a total of 907 BD patients and 1780 healthy controls. Altogether 32 polymorphisms were tested, using a Sequenom MassARRAY genotyping method in the first stage and a PCR-restriction fragment length polymorphism (PCR-RFLP) assay in the replication phase. Real-time PCR was performed to test the relative mRNA expression level of GAS6 and PROS1 from different SNP genotyped healthy individuals. The frequency of the C allele and CC genotype of rs9577873 in GAS6 (Pc = 4.92 × 10−5, Pc = 1.91 × 10−5, respectively) and A allele and AA genotype of rs4857037 in PROS1 (Pc = 1.85 × 10−6, Pc = 4.52 × 10−7, respectively) were significantly increased in BD. GAS6 expression in CC carriers of rs9577873 was significantly lower than that in CT/TT individuals (P = 0.001). Decreased expression of GAS6 and increased pro-inflammatory cytokines (IL-6 and IFN-γ: P = 4.23 × 10−4, P = 0.011, respectively) in individuals carrying the CC genotype suggest that the TAM-GAS6/PROS1 signal pathway may be involved in the pathogenesis of BD.

cytokine signals via the JAK-STAT pathway has also been identified as a risk factor for BD [10][11][12] . Other kinases may also play a role in BD development and this recently became a subject of interest in our laboratory.
The TAM receptors (TYRO3, AXL and Mer (gene name MERTK)) consist of a subfamily of receptor tyrosine kinases (RTK) 13 that have two vitamin-K related agonists: GAS6 and protein S (ProS, PS; gene name PROS1) [14][15][16] . This small subfamily of RTKs has been implicated as pivotal effectors 17 during immune modulation, including the regulation of innate immune reactions [18][19][20] , phagocytosis of apoptotic cells [21][22][23] , function and homeostasis of immune cells including dendritic cells 24,25 , monocytes/macrophages 26 , T cells 27 and NK cells 28,29 , vascular integrity and even pathogenesis of cancer 30 . GAS6 has a variable potency to activate these 3 receptors (Axl> Tyro3> Mer) while ProS was specific to the latter two receptors, but does not interact with Axl 31 . Abnormal TAM signal components have been detected in the demyelination process of both experimental autoimmune encephalomyelitis (EAE) models 32 and in multiple sclerosis (MS) patients 33 . Several SNPs within the MERTK gene have been identified from a genome-wide association study (GWAS) in MS in Australians and Europeans 33 . A mutation of MERTK has been reported to be associated with retinitis pigmentosa (RP) in the Royal College of Surgeons (RCS) rat and in humans because it compromised the phagocytic function of retinal pigment epithelial (RPE) cells [34][35][36][37] . A defective regulation of TAM receptors or their ligands has also been found in psoriasis 38 , primary Sjögren's syndrome 39 and inflammatory bowel disease (IBD) 40 . Previous studies also revealed that TAM-knockout mice easily develop systemic autoimmunity 41,42 . Appropriate apoptotic processes are essential to the development and maintenance of the immune system 43 . Abnormal clearance of apoptotic cells 22,44 or resistance to cell death is believed to play a role in the pathogenesis and recurrences of autoimmune disorders such as rheumatoid arthritis (RA), systemic lupus erythematosus (SLE) and BD 8,45 . Significantly altered serum concentrations of soluble forms of TAM-GAS6/PS components have been found in BD 46,47 , SLE 48 , RA 49 and IBD 50 patients suggesting a potential role of TAM related signals in their pathogenesis.
Further support for a role of the TAM-GAS6/PS family in immune mediated disease comes from genetic association studies. Gene polymorphisms of the TAM family were associated with several immune-related disorders, including systemic lupus erythematosus 51,52 , type 2 diabetes 53 , diabetic nephropathy and systemic sclerosis 33 . In view of the given roles of these genes in immune disorders, we hypothesized that TYRO3, AXL, MERTK, GAS6 and PROS1 may also be associated with BD. Although abnormal expression of TAM-ligands in several tumor-associated or immune related diseases including BD 46 and IBD 24 have been reported, the relation between TAM-GAS6/PROS1 signal pathway gene polymorphisms with uveitis has not yet been reported in Han Chinese and was therefore the subject of the study presented here.

Results
Clinical feature of enrolled BD cases. Clinical and demographic information of the enrolled BD cases and controls (1154 male/ 626 female) are shown in Supplementary Table S1. Since the patients visited a department of Ophthalmology, the enrolled BD patients all had uveitis. Almost 85% of BD patients were male whereas 64.8% of controls were male.
Genotype and allele frequencies of tested SNPs in cases and controls in the 1st, 2nd stage and combined studies. Thirty-two SNPs were genotyped in 412 BD and 612 healthy controls in the 1st stage study. The results showed significantly increased frequencies of the GAS6/rs9577873 C allele and CC genotype (P c = 3.72 × 10 −2 , OR = 1.662; P c = 2.42 × 10 −2 , OR = 1.772) and PROS1/rs4857037 A allele and AA genotype (P c = 2.99 × 10 −3 , OR = 2.015; P c = 2.37 × 10 −3 , OR = 2.117) (Tables 1 and 2). However, none of the remaining SNPs showed a significant association with ocular BD (Supplementary Table S2).
A stratified analysis was performed to study whether an association of polymorphisms of GAS/rs9577873 and PROS1/rs4857037 was associated with some of the main clinical features of BD. However, no significant association was found for these SNP genotypes/alleles and 5 clinical manifestations in BD (Supplementary  Table S3). Stratification for gender showed that genotype and allele frequency for both GAS6/rs9577873 and PROS1/rs4857037 showed a stronger significant difference in male (GAS6/rs9577873 C allele and CC genotype: P = 3.15 × 10 −6 , OR = 1.638; P = 2.55 × 10 −6 , OR = 1.720; PROS1/rs4857037 A allele and AA genotype: P = 4.16 × 10 −7 , OR = 1.840; P = 1.62 × 10 −7 , OR = 1.958) compared with female patients (Table 3).
Pairwise linkage disequilibrium (LD) and haplotype association analysis were performed using the SHEsis website. Six SNPs in the PROS1 gene (rs12634349-rs4857037-rs7616142-rs6803590-rs8178607-rs13062355) were in linkage disequilibrium with D' ranging from 0.914 to 1.00 and r 2 ranging from 0.009 to 0.707. Twelve SNPs in the GAS6 gene (rs9604488-rs7994900-rs7492052-rs6602910-rs12868833-rs7319547-rs7399860-rs9577924-rs 7323932-rs9604466-rs9577873-rs7399637) were also in linkage disequilibrium with D' ranging from 0.140 to 0.976 and r 2 ranging from 0.005 to 0.682. The global haplotype frequencies were significantly different between the case and control group (P < 0.001). Furthermore, we investigated two kinds of PROS1 haplotypes (AATACA; AATGCG) that were more frequent in the case group than in the normal control groups, whereas the GATACA haplotype and AATACG haplotype were less frequent in the BD group than in the normal group (Supplementary  Tables S4-8). There was no significant difference in the frequency distribution of the other haplotypes in these two groups.
Scientific RepoRts | 6:26662 | DOI: 10.1038/srep26662 The effect of rs9577873 on GAS6 expression and cytokine production. In order to find a possible functional association with GAS6/rs9577873, we carried out a real-time quantitative PCR assay to evaluate the GAS6 expression in PBMCs from known genotyped healthy individuals. We did not use patients for this functional analysis because inflammation and immunosuppressive drug use would confound the data. The results indicated that the GAS6 mRNA level in CC individuals of SNP rs9577873 was significantly lower compared with CT/TT individuals (Fig. 1a, P = 0.001). Several cytokines as well as Th1 and Th17 responses have been shown to play an important role in BD pathogenesis. GAS6 has been reported as an inhibitor of pro-inflammatory cytokines in monocytes/macrophages 26 . Further investigations were therefore performed to study whether various genotypes of rs9577873 altered cytokine production by PBMCs from healthy individuals. An elevated secretion of IL-6 ( Fig. 2b, P = 4.23 × 10 −4 ) and IFN-γ ( Fig. 2d, P = 0.011) by stimulated PBMCs was observed in CC genotype when compared to CT genotype carriers. The different genotype carriers did not show differences concerning IL-1β , TNF-α , IL-17 or IL-10 production by stimulated PBMCs (Fig. 2a,c,e,f).
The effect of rs4857037 on PROS1 expression and cytokine production. We also investigated the functional consequences of polymorphisms of rs4857037 in the PROS1 gene. No significant association in the expression of PROS1 by non-stimulated PBMCs from various rs4857037 genotype carriers could be detected (Fig. 1b). We also investigated whether the expression of PROS1 was affected by different rs4857037 genotypes in LPS-stimulated PBMCs. The PROS1 mRNA level in AA individuals of SNP rs4857037 was significantly higher than in AG/GG individuals when PBMCs had been stimulated with LPS ( Fig. 1b, P = 0.033). Cytokine expression (IL-1β , IL-6, TNF-α , IL-10, IFN-γ and IL-17) by stimulated PBMCs was however not associated with PROS1/ rs4857037 genotype (Supplementary Fig. S1a-f).

Discussion
In the current study, we describe a novel association between receptor tyrosine kinase pathway (TAM-GAS6/ PROS1) genes with ocular BD in a Chinese Han population. The observed association was restricted to the ligands GAS6 and PROS1 and no association could be detected with gene polymorphisms of the TAM receptors-TYRO3,AXL and MERTK. A significantly higher CC genotype and C allele carrier frequency of the GAS6/ rs9577873 gene and AA genotype and A allele carrier frequency of the PROS1/rs4857037 gene was observed in BD. Functional analysis showed a significant down regulation of GAS6 mRNA expression by PBMCs in CC genotype compared with CT or TT carriers. PROS1/rs4857037 AA genotyped individuals showed an increased PROS1 expression compared to AG/GG carriers in LPS-stimulated PBMCs but not in non-stimulated cells. An elevated secretion of IL-6 and IFN-γ by LPS or anti-CD3 combined with anti-CD28 antibodies treated PBMCs was observed in CC genotype compared to CT genotype of GAS6/rs9577873 carriers. Stratification for gender showed that genotype and allele frequency for both GAS6/rs9577873 and PROS1/rs4857037 showed stronger significant differences in male compared with female patients. The observed SNP associations may thus be related via a regulation of gene transcription and modulation of inflammatory cytokine expression thereby leading to a higher risk of developing ocular BD. TAMs belong to a family of tyrosine kinase receptors that have received little attention until recently, and were mostly related to TAM involvement in cancer 17 . Our study confirms earlier studies showing an association between the TAM-GAS6/PS pathway and immune related disorders such as psoriasis 38 , primary Sjögren's syndrome 39 , SLE 48 , RA 49 and IBD 40,50 . Studies on the role of gene polymorphisms in this pathway were limited to vascular and some autoimmune diseases including stroke 54   Vasculitis is a prominent feature of BD and it is possible that the role of TAM-GAS6/PS in this disease may be related to its involvement in the control of vascular integrity 17 . The TAM-GAS6/PS signal pathway has been shown to play a pivotal role in blocking the innate immune response and has also been associated with the regulation of cancer development 17 . The TAM-GAS6/PS system is widely expressed in a variety of immune cells, and has been shown to play a role in the EAE 32 model and in retinitis pigmentosa 35 . TAM-GAS6/PS also plays critical roles in the maintenance of retinal function via its effect on retinal pigment epithelial phagocytosis and altering the expression of microRNAs 37,55 . It has also been shown to play a role in a mouse model of human uveitis known as experimental autoimmune uveoretinitis (EAU), which can be induced by immunization of the experimental animals with retinal autoantigens, such as interphotoreceptor retinoid-binding protein (IRBP). It has been  reported that TAM receptor knockout mice were more susceptible to IRBP immunization supporting its negative control of the inflammatory process 42 . Our findings are in agreement with these previous studies showing that suppression of the TAM signal or knocking-out TAM in mice makes them more prone to autoimmunity 41,42 . Apoptosis is essential to the development and maintenance of the immune system 56 and the TAM-GAS6/ PS system has been shown to mainly exert its role in immunity by its effect on the phagocytosis of apoptotic cells 21,23,25 and thereby maintaining immune homeostasis 14,18 . Changes in the plasma concentration of Axl/Gas6 have been reported in rheumatic diseases such as SLE 48 , BD 46,47 and RA 49 . Significantly increased Axl levels were also observed in patients with BD 47 . This may be caused by a dysfunctional feedback or interaction of the TAM-GAS6/PS system resulting in an elevated Axl level. There is a growing body of evidence suggesting that chronic inflammation in different autoimmune diseases such as BD 45,57 , RA 58 and multiple sclerosis is caused by the apoptosis refractory nature of activated T cells and a continuous stimulation caused by apoptotic debris 21 . The abnormal apoptotic process seems critical for autoimmune disease pathogenesis and may be involved in causing the recurrent and chronic character of the disease. Collagen-induced arthritis (CIA) is used as an animal model for RA and several studies have shown that agonists of the TAM pathway may be used to block inflammation in this model 59 .
Earlier GWAS approaches in BD did not reveal an association with the TAM-GAS6/PS pathway. The reason for this discrepancy may be explained as follows. Although GWAS is a powerful approach to scan and find complicated disease related loci, it is based on the assumption of indirect association mapping using reference linkage SNPs and has a strict P value set at the 10 −8 level. Its result may vary according to ethnicity, sample sources, sample size and GWAS chip coverage in the discovery phase. For example, the Affymetrix 500k chip (Affymetrix Company, Santa Clara, CA, USA) can capture about 65% of the common variants and the Ilumina 317k chips (Illumina Company, San Diego, CA, USA) can capture about 75% of the common variants. GWAS therefore does not cover all possible common variants. Our study used a p value with a lower threshold as used for GWAS and may therefore have picked up loci not detected by GWAS.
Previous findings suggested that TAM-ligands may be associated with BD activity especially in neuro-BD, where the serum GAS6 level was shown to be decreased in BD 47 . Although the consequences of diminished mRNA levels of GAS6 in CC genotype carriers is not clear, it may suggest an involvement of the TAM-GAS6/PS signal pathway during abnormal apoptotic processes 21 in patients with BD. Further experiments investigating apoptotic function and hnRNA expression in relation to TAM polymorphisms are necessary to clarify this subject. Further confirmation of the role of this pathway in BD pathogenesis is necessary and this may involve the measurement of TAM-GAS6/PS related products in the blood of our BD patients. It is also not clear whether the observed association is confined to BD patients with uveitis or whether it can also be observed in other uveitis entities. Since BD is a multisystemic disease, it would also be interesting to investigate whether the observed association is restricted to ocular BD. Future investigations are therefore needed in a large cohort of BD patients recruited via other medical departments including rheumatology and dermatology to address this issue. Our study was performed in a Han Chinese population and verification studies are needed in other ethnic populations. Our study does not exclude the possibility that other SNPs or copy number variants (CNVs) of TAM and its ligand genes are associated with this disease. Additional SNPs around the tagging SNP should be genotyped to identify the underlying haplotype. Further investigations are needed to address these questions and studies concerning specific agonists or analogues targeting this pathway may offer therapeutic opportunities in the near future.
In conclusion, our results show that polymorphisms in the genes encoding the ligands of the receptor tyrosine kinase family, GAS6 and PROS1 confer genetic susceptibility for ocular BD in Han Chinese.

Materials and Methods
Study participant recruitment. A total of 412 patients who fulfilled the criteria for Behçet't Disease according to the International Study Group diagnostic criteria 60 were included in the first phase study. Six hundred and twelve age, geographically and ethnically matched healthy Chinese Han volunteers served as controls. Another set of 495 BD patients and 1168 healthy controls were included in the replication study. Ethical considerations. The experimental protocols and study design were approved by the local ethical research committee of the First Affiliated Hospital of Chongqing Medical University. All experiments were carried out in accordance with the approved guidelines. The ethical standards of the Declaration of Helsinki were followed during all the experimental procedures. All study participants were well informed and signed an informed consent before their enrollment.
Tag SNP selection. The choice of SNPs was mainly based on tagSNPs. Thirty-two tagSNPs involving 5 TAM signal genes were chosen in the present study. After a search in the public database HapMap and HaploView (V4.0; Daly lab at the Broad Institute, Cambridge, MA, USA) and specific analysis for the Han Chinese in Beijing (CHB) population, our candidate tagSNPs were chosen based on a minor allele frequency (MAF) > 0.05 and r 2 was set at 0.8. We chose a total of thirty-two SNPs: two in AXL, one in TYRO3, eleven in MERTK, twelve in GAS6 and six in PROS1.
Genomic DNA preparation and SNP genotyping analysis. Peripheral  a Nanodrop 2000 (Thermo Fisher Scientific, Wilmington, DE, USA), quality checked, standardized and stored at − 20 °C until assayed. The primers used for genotyping were designed by MassARRAY Assay design software. SNP genotyping in the discovery cohort was determined using the Sequenom MassARRAY system platform (Sequenom Inc, San Diego, California, USA) and iPLEX reagents according to the manufacturer's instructions (Agena Bioscience, California, USA). The PCR reaction was performed on the GeneAmp PCR System 9700 instrument (ABI, Foster City, CA, USA). Subjects in the replication phase were genotyped using the PCR-RFLP method. Appropriate primers and restriction enzymes (Table 4) were used to amplify the target DNA sequence and digest the PCR product respectively. After digesting the PCR products with restriction enzymes at 37 °C, the digested products were subjected to electrophoresis using 4% or 5% agarose gels, stained by GoldView (SBS Genetech, Beijing, China) and visualized under a Bio-Rad imaging system. Randomly selected samples (3%) were checked using direct sequencing by Sangon Biotech Company (Sangon Biotech, Shanghai, China) and the results were in complete accordance with the Sequenom MassARRAY system genotyping and PCR-RFLP outcomes.
Cell isolation and culture. Fresh peripheral blood from healthy volunteers was used to isolate peripheral blood mononuclear cells (PBMCs) by Ficoll-Hypaque density gradient centrifugation. A density of 1 × 10 6 PBMCs per well were seeded into 24-well culture plates with complete medium RPMI 1640 (containing 10%FBS, 100 U/mL penicillin and 100 μ g/mL streptomycin) and stimulated with 100 ng/mL lipopolysaccharide (LPS) (Sigma, Missouri, USA) for 24 h to detect IL-1β , TNF-α and IL-6 production. Antigen stimulation was simulated using a cocktail of anti-CD3 and anti-CD28 antibodies (5:1) (Miltenyi Biotec, Palo Alto, CA) for 3 days at 37 °C in a humidified 5% CO 2 incubator where after IL-10, IL-17 and IFN-γ were measured in the culture supernatants.

RNA preparation, reverse transcription and Real-time PCR.
Total RNA was extracted from LPS-stimulated PBMCs and non-stimulated PBMCs using the TRIzol (Invitrogen, San Diego, California, USA) method. Reverse transcription was then performed using a Takara transcriptase kit (Takara, Dalian, China). Relative mRNA expression assays were performed on an ABI 7500 real-time instrument. The expression of GAS6 and PROS1 mRNA relative to β -actin was calculated using the 2 −ΔΔCt method. The GAS6 and PROS1 expression was detected using primers reported elsewhere 39 . All samples were evaluated in duplicate with at least three experimental replicates.
Cytokine ELISAs. The concentration of IL-1β , TNF-α , IL-6, IL-10, IL-17 and IFN-γ in the culture supernatants of LPS or anti-CD3 combined with anti-CD28 antibodies treated PBMCs was measured using the human Duoset ELISA development kit (R&D Systems, Minneapolis, Minnesota, USA). Recombinant proteins were used to generate the standard curve. Test results were evaluated using an ELISA reader (SpectraMax M2e, Molecular Devices, USA) at 450 nm. Statistical analysis. Genotype and allele frequency data were calculated by Typer4.0 software from the MassARRAY system platform or by a direct count from the PCR-RFLP result. The SHEsis website was used to test the Hardy-Weinberg equilibrium (HWE) in controls of all tested SNPs. No HWE deviation was observed. Pairwise linkage disequilibrium (LD) and haplotype association analysis were also performed using the SHEsis website. SPSS (SPSS Inc., Chicago, Illinois, USA) version 17.0 was used to analyze the χ 2 test, P value, odds ratio (OR) as well as 95% confidence intervals (95% CIs) data. To correct for multiple comparisons, the P values were adjusted as corrected P values (P c ) with the Bonferroni correction approach according to the number of analyses performed. Statistical significance level was set at P c < 0.05. Continuous variables were summarized through the mean ± SD. Student t test or the Nonparametric Mann-Whitney U test was selected for independent group comparisons. P values of 0.05 or less (Two-tailed) were considered as statistically significant.