A Missense LRRK2 Variant Is a Risk Factor for Excessive Inflammatory Responses in Leprosy

1Leprosy is a chronic dermato-neurological infectious disease caused by M. leprae. Leprosy irrespective of its clinical presentation is curable by multi-drug therapy. The current global effort in leprosy control is focused on early recognition of leprosy cases and prevention of permanent disabilities. A common complication of leprosy are excessive pro-inflammatory episodes termed type-1 reactions (T1R). If untreated, T1R can lead to irreversible nerve function impairment due to a pathological cellular immune response directed against host peripheral nerve cells. Up to 50% of all leprosy cases can undergo T1R with the incidence varying according to endemic settings and criteria for case definition.

2Why only a proportion of leprosy patients undergo T1R is not known. However, clinical and environmental factors have been associated with T1R outcome. Individuals categorized as borderline in the Ridley and Jopling clinical spectrum of leprosy are at increased risk to develop T1R while patients with tuberculoid or lepromatous polar leprosy forms rarely develop T1R. Positive bacillary index, PCR detection of M. leprae, and increased age at leprosy diagnosis are other factors associated with T1R-risk. A number of studies have shown a consistent upregulation of pro-inflammatory cytokines, i.e. TNF, IFNγ and the chemokine IP10, in the blood of T1R patients. In a prospective study, the transcriptional profile of leprosy cases destined for T1R displayed a distinct signature from leprosy patients that remained T1R-free. A dysregulated balance between innate pro- and anti-inflammatory responses emerged from this study as a key factor in T1R outcome.

3A number of genes have been shown to be associated with T1Rincluding TNFSF15 and the pathogen recognition genes TLR1, TLR2 and NOD2. All of these genes had also been found to be associated with leprosy per se in other studies. Since leprosy patients are usually not stratified by their T1R status, it is possible that some of the leprosy per se associations were caused by T1R patient subgroups. For example, the TNFSF15 gene had initially been shown to be associated with leprosy per se by a genome wide association study (GWAS) in a Chinese population. However, in a Vietnamese sample stratified by T1R status the association signal could be unambiguously assigned to the T1R group. Of the genes reported by the GWAS, the TNFSF15 and LRRK2 were the only leprosy susceptibility genes not validated for association with leprosy per se in a Vietnamese population. Like TNFSF15, LRRK2 is a gene with an uncertain role in leprosy per se susceptibility. Several groups have evaluated the association of the LRRK2 gene with leprosy susceptibility but results were inconsistent. Given that T1R affects different proportions of leprosy cases according to the studied population, we wondered if inconsistencies in LRRK2 association with leprosy per se were due to different proportions of T1R in each setting. Here, we evaluated a possible role for LRRK2 in T1R-affected families and contrasted the results with T1R-free families. We identified a set of 18 SNPs in LRRK2 preferentially associated with T1R. These variants overlapped with previous associations reported for Crohn’s Disease (CD), Ulcerative Colitis (UC) and Inflammatory Bowel Disease (IBD).

4For the LRRK2 study, a total of 1372 individuals were selected. These individuals were divided in two family-based samples. The first set of families contained 229 leprosy affected offspring that underwent T1R (T1R-affected) and their respective parents. The T1R-free subset was matched to the T1R-affected subset by leprosy clinical subtype of the offspring (Fig 1). Consequently, the second set of families included 229 leprosy affected offspring and their parents in which the offspring had no signs of leprosy reaction (T1R-free). There was no difference in gender and age at leprosy onset regarding T1R outcome between both subsets (S1 Table). The subjects included in the eQTL analyses were part of a study evaluating the transcriptional profile of leprosy patients prior to T1R onset. Briefly, 53 newly diagnosed leprosy cases in the borderline spectrum (19 BT, 30 BB and 4 BL) were enrolled. A blood sample was collected from each subject within 3 months of leprosy diagnosis and none of the subjects suffered T1R at enrolment.

5The samples used in the current study were selected from our records at the Dermato Venerology Hospital, Ho Chi Ming City, Vietnam, as described previously. Written informed consent was obtained from all subjects enrolled in the study and all subjects were anonymized. This study was approved by the regulatory authorities in Ho Chi Minh City, Vietnam, and the Research Ethics Board at the Research Institute of the McGill University Health Centre, Montreal, Canada. The investigation have been conducted according to the principles expressed in the Declaration of Helsinki.

6Genotypes for156 SNPs mapping to a 500 kb window overlapping the LRRK2 and MUC19 genes were obtained via the 660W-quad v1 Illumina array. The variants selected covered 89% of the SNPs with a MAF > 5% at a r2 > 0.5 for the Vietnamese (KHV) and Chinese (CHB) and 84% for Caucasians (CEU) populations from the 1000 genomes project. All genotypes passed standard quality-control presenting call rates greater than 98%, less than 2 Mendelian Errors (ME) and were in Hardy-Weinberg Equilibrium (HWE) with P > 0.05 in 763 leprosy unaffected parents from both T1R-affecte and T1R-free subsets. LRRK2 expression levels were obtained using Illumina HumanHT12 v4 BeadChips as previously described. Briefly, whole blood from the 53 leprosy patients was divided in two aliquots. One aliquot was stimulated with M. leprae sonicate for 24hrs to 30 hrs while the second aliquot was incubated for the same time interval in the absence of M. leprae sonicate (non-stimulated). Total RNA was extracted from all aliquots and used for LRRK2 quantification.

7Family based SNP and haplotype association tests were performed using Transmission Disequilibrium Test (TDT) as implemented in FBAT 2.0.4. Association testing was carried out under the same genetic model in T1R-affected and T1R-free families and the P values for the best genetic model (additive or dominant) were displayed (Fig 1). Due to the highly correlated nature of the genotyped SNPs, we did not perform a correction for multiple testing. Subsequently, a formal heterogeneity test was performed to evaluate preferential association of genetic variants with T1R (Fig 1) by using a modified version of the FBAT statistics (FBATHet) as described by Gaschignard, J. et al. A multivariate analysis was performed to test for independence of T1R-specific associations. For each SNP bin (r2 > 0.5), SNP with the most significant evidence for association was included in the multivariate model. Multivariate analyses were done by stepwise conditional logistic regression (SAS v.9.3). Logistic regression was also used to estimate the odds ratio for each individual SNP in the T1R-affected subset. Briefly, the TDT evaluates the non-random transmission of alleles from heterozygote parents to affected offspring. We used the non-transmitted alleles from the TDT to create up to three unaffected pseudo-sibs per family, one for each possible genotype. We compared the original T1R-affected offspring with T1R-unaffected pseudo sibs in a matched case-control design as described in. Under the additive model, the TDT and the conditional logistic regression statistics result in the same P values. The recombination rate in centimorgan by mega base (cM/Mb) according to the 1000 genomes and the LD structure of the LRRK2/MUC19 locus were obtained with the R packages Locuszoom v.1.3 and snp.ploter v.0.5.1, respectively. The minor allele frequencies (MAF) and HWE for each SNP were estimated using Haploview 4.2. For the eQTL analyses, the correlation between genotypes and gene expression levels was performed with a simple linear regression under both stimulated and non-stimulated conditions using R 3.2.0.

8Genotype associations of LRRK2 and PD were obtained from the PDgene database (www.pdgene.org). LRRK2 variants associated with PD were obtained from a GWAS meta-analysis of 14 studies with a total population sample of 12,771 PD cases and 93,386 controls. Genotype associations of LRRK2 variants with IBD and CD were obtained from the IBDgenetics database (www.ibdgenetics.org). The IBDgenetics population sample consisted of 42,950 IBD cases (22,575 CD and 20,417 UC patients) and 53,536 controls.

9A genome-wide association (GWAS) in a Chinese sample provided the most comprehensive study of LRRK2 in leprosy. Hence, we first analyzed six LRRK2 SNPs (rs1873613, rs10878220, rs1491938, rs12820920, rs11174812 and rs11173979) described in association with leprosy per se by a GWAS in a Chinese population. None of the six variants presented significant evidence of association with leprosy in the T1R-free subset (Table 1). In contrast, SNPs rs10878220, rs1491938 and rs12820920 that belonged to the same r2 < 0.5 SNP bin were nominally associated with disease in the T1R-affected subset (Table 1). These SNPs were preferentially associated with T1R when a formal heterogeneity test was performed (Table 1). Moreover, the direction of association was the same as previously reported for leprosy per se. While SNP rs10878220 is located in the core promoter region of the LRRK2 gene, SNPs rs1491938 and rs12820920 are LRRK2 intronic variants (Table 1). The remaining three SNPs that did not show association with T1R are located 66kb to 367kp upstream of LRRK2 transcription start site (Table 1).

10To further evaluate the association of LRRK2 with T1R, an additional 150 SNP were selected from the LRRK2 gene region. Of all 156 SNPs evaluated for association in the T1R-affected subset, 34 showed P < 0.05 including the initial three SNPs associated with T1R (Fig 2). In the T1R-free family subset, the SNP rs7972711 situated near the 3’ end of the MUC19 gene and rs10878434, rs7303525 and rs11564172 located near 3`end of LRRK2 provided evidence of association (P < 0.05; Fig 2). None of the 34 SNPs associated with T1R showed evidence of association with leprosy in the T1R-free families (Fig 2). When formally tested for heterogeneity of association, 18 out of the 34 SNPs were preferentially associated with T1R (Fig 2 and Table 1). These 18 SNPs preferentially associated with T1R belong to two extended SNP bins (r2 > 0.5) and two single SNP bins (Table 1). Notably, although not nominally significant these 18 SNPs showed the opposite allelic enriched in T1R-free subset relative to the T1R-affected subset (Table 1). The bin tagged by the missense M2397T variant (rs3761863) situated in the WD40 domain of LRRK2 displayed the strongest preferential association with T1R (P = 0.003; odds ratio (OR) = 1.49; 95% confidence interval (CI) = 1.12–1.97 and P Het = 0.008 for M2397 allele under an additive model) in the same direction as previously reported for CD (Table 2).

11To test if the r2 > 0.5 SNP bins were independently associated with T1R we performed a multivariate analysis including the most significantly associated SNP for each bin (Table 2). The multivariate analysis identified the missense M2397T variant as the main association with T1R (P multi = 0.01; Table 2). However, a trend for independent association was observed for the SNP bin tagged by rs1031996 (P multi = 0.05) and the single SNP bin rs1463739 (P multi = 0.13; Table 2) while the association of T1R and rs1427271 disappeared.

12We further investigated the combined effect of M2397T and rs1031996 by conducting a haplotype analysis. We found that the M2397 allele was mostly observed in the presence of the rs1031996 C allele (48% of the T1R-affected offspring) and this haplotype displayed a strong risk effect (Table 3). Importantly, the trend towards protection observed on the T2397 background was not modulated by rs1031996 alleles. Hence, the haplotype analysis confirmed that the main effect of the LRRK2 gene on T1R susceptibility was driven by the M2397T variant. No haplotype association was observed when leprosy per se was considered as phenotype.

13To investigate if SNP alleles associated with T1R were correlated with LRRK2 transcriptional levels we performed an eQTL analysis. Of the 18 SNPs preferentially associated with T1R, nine variants belonging to the same SNP bin as the missense M2397T variant (r2 = 0.5) were eQTLs for LRRK2 in non-stimulated whole blood of 53 individuals (Table 1, S1 Fig). When a tighter LD threshold (r2 = 0.8) was considered, the eQTL variants were separated from the bin tagged by M2397T. This observation suggested that the modest eQTL effect of M2397T was due to the linkage disequilibrium with the causal eQTL. The strongest eQTL effect was observed for SNP rs2404580 with the T1R-risk “T” allele being associated with higher LRRK2 expression in unstimulated cells (P = 5.1E-05; Fig 3). Following stimulation with M. leprae sonicate, an abrogation of the eQTL effect was observed for all nine SNPs (Figs 3 and S1). Clinical subtypes of leprosy had no detectable impact on the eQTL effect of LRRK2 genotypes.

14We identified an amino acid change M2397T in the WD40 domain of LRRK2 and a bin tagged by the variant the variant rs1031996 as being associated specifically with T1R. LRRK2 is a protein that exerts a diverse set of functions. LRRK2 mediates catalytic processes through its enzymatic ROC/COR domain; facilitates signal transduction through a MAPK domain and interacts with other proteins through three scaffold domains, an ankylin repeat (ANK), a leucine rich repeat (LRR) and a WD40 repeat domain. Scaffold domains are also responsible for protein conformation and stability. Although, LRRK2 displays multiple functions, the association of T1R with a variant in the WD40 domain of LRRK2 suggests that protein conformation and/or stability are key factors in T1R. Indeed, the M2397T variant in LRRK2 has previously been shown to impact on LRRK2 protein turnover. The half-life of LRRK2 with the T1R-risk allele M2397 had been estimated at approximately 8 hrs which is substantially shorter than the estimated 18 hrs half-life of the T2397 allele of LRRK2 (Fig 4). Cytoplasmic LRRK2 forms a complex that arrests nuclear factor of activated T-cells (NFAT) in the cytoplasm. A consequence of LRRK2 deficit in the cytoplasm is the translocation of NFAT to the nucleus, which strongly induces the transcription of pro-inflammatory cytokines (Fig 4). Thus, the association of the M2397 allele with T1R-risk is in agreement with an exacerbated pro-inflammatory response in T1R cases (Fig 4). While the M2397T variant is a strong candidate for functional impact in T1R, the association of the bin tagged by rs1031996 requires further investigations.

15A set of eQTLs for LRRK2 was observed in non-stimulated cell from leprosy patients. Due to LD the M2397 allele and the eQTL alleles that correlated with increased LRRK2 expression were preferentially observed on the same haplotype. This suggests that the effect of M2397 on faster LRRK2 turn-over is mitigated by higher levels of LRRK2 message. Interestingly, in the presence of mycobacterial antigen the compensatory effect of higher LRRK2 message for the more rapid turnover of the M2397 protein is strongly abrogated. This implies that the genetic effect of the M2397 allele will be more pronounced in the presence of mycobacterial antigen while it may be largely abolished without such antigen exposure. While these in-vitro findings in whole blood will need to be validated employing purified cell types and in-vivo studies in leprosy lesion, the results obtained provide an example of how environmental stimuli can modulate germline encoded genetic risk factors (Figs 3 and 4).

16LRRK2 variants associated with T1R overlapped previous reported LRRK2 associations with CD and PD. In PD, rare coding variants in the enzymatic and kinase domains of LRRK2 were shown to be causally linked to PD. In addition, common LRRK2 variants were shown in association with PD. Two of these variants, rs1491932 and rs7970326, were observed in association with T1R and borderline evidence T1R specificity for the alleles opposite to PD. Common variants may tag rare variants with stronger effects. However, given our sample size we were unable to evaluate the role of rare variants in T1R. Of 18 SNPs preferentially associated with T1R, 17 were nominally associated with IBD with the same risk allele observed for T1R. The only exception was rs1463739 a SNP located outside of LRRK2 in the 3’ region of MUC19 (Table 2 and Fig 2). When considering IBD subtypes, 14 T1R-risk variants were associated with risk of CD and 11 with risk of UC. The M2397 allele is a risk factor for T1R, CD and UC, suggesting that a faster LRRK2 turnover leads to an increased pro-inflammatory response that is common to these diseases. T1R and CD may share susceptibility to mycobacterial species as common etiology. In PD, pathogen involvement is controversial although some studies suggested that Helicobacter pylori and prions might play a role in disease susceptibility. Further studies will be needed to understand the precise role of LRRK2 in the pathogenesis of these three diseases.