Yue-miao Zhang1,2,3†
Xing-zi Liu1,2,3†
Xu-jie Zhou1,2,3
Li-jun Liu1,2,3Su-fang Shi1,2,3Ping Hou1,2,3Ji-cheng Lv1,2,3
Hong Zhang1,2,3*- 1Renal Division, Department of Medicine, Peking University First Hospital, Beijing, China
- 2Key Laboratory of Renal Disease, Ministry of Health of China, Beijing, China
- 3Key Laboratory of Chronic Kidney Disease Prevention and Treatment (Peking University), Ministry of Education of China, Beijing, China
C-C chemokine receptor 6 (CCR6) is a susceptibility gene of various immune-related diseases, which was suggested to be shared with immunoglobulin A nephropathy (IgAN). In this study, we aimed to identify the functional variants. First, we analyzed the associations of CCR6 common and rare variants detected by multi-platform chips with IgAN susceptibility using imputation and identified 68 significantly associated common variants located in the regulatory region. Among them, rs3093023 showed both statistical significance (rs3093023-A, odds ratio [OR] = 1.15, P = 2.00 × 10−2) and the expression quantitative trait loci (eQTL) effect (P = 1.45 × 10−3). It was independently replicated (rs3093023-A, OR = 1.18, P = 5.56 × 10−3) and the association was reinforced in the meta-analysis (rs3093023-A, OR = 1.17, P = 6.14 × 10−7). Although rs3093023 was in a strong linkage disequilibrium with the reported CCR6 functional variant dinucleotide polymorphism, CCR6DNP, the alleles of rs3093023 (G>A) rather than of CCR6DNP were shown differential nuclear protein binding effect by electrophoretic mobility shift assay. The RegulomeDB and JASPAR databases predicted Pou2f1 as the potential transcription factor, which was negatively associated with CCR6 mRNA (r = −0.60, P = 3.94 × 10−9). At the mRNA level, the eQTL effect of CCR6 was validated (P = 4.39 × 10−2), and CCR6 was positively associated with the expression of CCR4 and IL-17A rather than that of CXCR3 and IFNG. At the protein level, a higher CCR6+ cell ratio was observed in a risk allele dose-dependent manner in lymphocytes (P = 3.57 × 10−2), CD3+ T cells (P = 4.54 × 10−2), and CD4+ T cells (P = 1.32 × 10−2), but not in CD8+ T cells. Clinical-pathological analysis showed that rs3093023 risk allele was significantly associated with diastolic blood pressure, serum creatinine, and high ratio of tubular atrophy/interstitial fibrosis. Overall, the rs3093023 was prioritized as the function variant in CCR6, which may contribute to IgAN susceptibility by regulating Th17 cells.
Introduction
Immunoglobulin A nephropathy (IgAN) is the most common form of primary glomerulonephritis worldwide and is a major cause of end-stage renal disease in young adults. The variable prevalence among different ethnicities (1) and high familial aggregation (2) suggest that genetic components play an important role in IgAN pathogenesis. However, the known susceptibility loci of IgAN together explains no more than 10% of the disease heritability (3–7), suggesting that the genetic loci remain to be investigated.
Genetic studies have demonstrated the association between C-C motif chemokine receptor 6 gene (CCR6) and susceptibility to various immune-related diseases, such as rheumatoid arthritis (8, 9), psoriasis (10), and lupus nephritis (11). CCR6 is the unique receptor for chemokine C-C motif chemokine ligand 20 (CCL20) and is a specific surface marker for interleukin (IL)-17- producing Th17 cells (12). CCR6–CCL20 signaling was reported to play an important role in recruiting Th17 cells to inflammatory sites (13). For example, synoviocytes from mouse and human arthritic joints could produce CCL20 and consequently recruit Th17 cells to the inflamed joints (14). In patients with psoriasis, IL-17-treated epidermal keratinocytes promoted CCL20 production, mediating Th17 cell entry to the lesional skin (15). Furthermore, patients with lupus nephritis have increased levels of Th17 cells (16); in glomerulonephritis, Th17 cells migrate from gut to kidneys via CCR6–CCL20 signaling (17).
Recent studies have highlighted the roles of Th17 cells and CCR6–CCL20 signaling in IgAN (18). Th17 cell levels were increased in patients with IgAN compared with healthy controls (19). IgAN mice infected with hemolytic streptococcus showed higher Th17 cell levels in renal lymphocytes, while the addition of CCL20 antibodies could suppress the increase (20). Furthermore, IgA1-stimulated mesangial cells from patients with IgAN produced more CCL20 and attracted CCR6+ Th17 cells to kidneys (21). In genetic studies, CCR6 was suggested to be a shared gene of IgAN (6); however, the functional variants associated with IgAN susceptibility were still unclear.
In this study, we conducted a two-stage genetic association study to identify CCR6 functional variants associated with IgAN susceptibility, followed by in silico regulatory effect prediction, expression association analysis, and clinical-pathological association analysis.
Participants and Methods
Participants
In the discovery stage, three independent cohorts were enrolled (Figure 1). Cohort 1 consisted of 1,228 patients with IgAN (age: 31.10 ± 10.70 years, male ratio: 54.00%) and 966 healthy controls (age: 31.54 ± 8.43 years, male ratio: 59.13%), cohort 2 consisted of 500 patients with IgAN (age: 36.70 ± 12.50 years, male ratio: 54.70%) and 2,372 healthy controls (age: 34.81 ± 14.9 years, male ratio: 55.32%), and cohort 3 consisted of 640 patients with IgAN (age: 36.50 ± 13.50 years, male ratio: 48.50%) and 320 healthy controls (age: 36.47 ± 13.45 years, male ratio: 56.56%). For replication, we enrolled 1,030 patients with IgAN (age: 35.91 ± 13.07 years, male ratio: 51.55%) and 1,228 healthy controls (age: 38.6 ± 12.97 years, male ratio: 51.47%). All participants were unrelated northern Han Chinese. The study was approved by the Medical Ethics Committee of Peking University First Hospital (No. [2019]320), and the informed consent was provided by the participants.
Figure 1. Flowchart of the study design. In the discovery stage, we conducted a genetic association study between the CCR6 common and rare variants and IgAN using multi-platform chips. The replication stage was aimed at replicating the functional loci discovered in the discovery stage using an independent population. After the two-stage genetic association analysis, we conducted in silico functional prediction, expression association analysis, and clinical-pathological association analysis.
All patients with IgAN were diagnosed based on the dominant deposition of IgA in the glomerular mesangium by immunofluorescence, while patients with secondary causes, such as IgA vasculitis, systemic lupus erythematosus, or liver cirrhosis, were excluded. Their baseline demographic and clinical data at renal biopsy, such as age, sex, systolic blood pressure (SBP), diastolic blood pressure (DBP), serum creatinine (Scr), and 24 h proteinuria, were collected. The estimated glomerular filtration rate (eGFR) was calculated using the Chronic Kidney Disease Epidemiology Collaboration equation (22). Pathological characteristics were scored according to the Oxford classification (23).
Genotyping of the CCR6 Variants
In the discovery stage, variants in a region spanning 50 kb (chromosome 6: 167, 510–167, 560 kb) encompassing CCR6 and its downstream and upstream were derived from Illumina chips (USA). More specifically, to analyze the associations of common variants, genotyping information was obtained from the Illumina Human 610-Quad chip (cohort 1) and the Illumina Human OmniZhongHua-8 chip (cohort 2). To analyze the associations of rare variants, genotyping information was obtained from the Illumina Human Exome chip (cohort 3). In the replication stage, the candidate functional variants were genotyped using TaqMan allele discrimination assays (Applied Biosystems, USA).
In silico Functional Annotation
To identify the candidate functional variants, regulatory features were annotated using the RegulomeDB database (24) (http://www.regulomedb.org), and the expression quantitative trait loci (eQTL) were analyzed using the blood eQTL database (25) (http://www.genenetwork.nl/bloodeqtlbrowser/). The potential transcription factor was predicted using the RegulomeDB and JASPAR (26) (http://jaspar.binf.ku.dk/) databases. The regulatory effects of the functional SNP and its proxies, such as enhancer histone marks, promoter histone marks, and DNase hypersensitivity, in different tissues/cells were annotated using the HaploReg database (27) (http://www.broadinstitute.org/mammals/haploreg/haploreg.php).
Sequencing of CCR6 Dinucleotide Polymorphism (CCR6DNP)
A triallelic dinucleotide polymorphism of CCR6, CCR6DNP, was reported to be a disease-causing variation that can affect the CCR6 expression in rheumatoid arthritis (8). To further investigate the potential functional variant, we determined whether CCR6DNP could tag our candidate functional variants. We randomly selected 361 people (166 cases vs. 195 controls) based on the genotypes of the candidate functional variants. Genomic DNA was extracted from whole blood samples and amplified using forward (5′-CAACCACCTTTGAAAGAGCAG-3′) and reverse (5′-CCCTTGTTCATCCCAACCT-3′) primers (8). The DNA products underwent Sanger sequencing, and the signals were analyzed by Chromas (Technelysium Pty, Australia) (Supplementary Figure 1).
Isolation of Peripheral Blood Mononuclear Cells (PBMCs)
PBMCs were isolated from fresh peripheral blood samples collected from the participants using Ficoll-Paque density gradient centrifugation. Briefly, each sample was centrifuged at 600 g for 30 min to separate plasma. The sample was mixed with phosphate-buffered saline (PBS) in a 1:1 ratio and then slowly transferred onto Ficoll-Paque PLUS (GE Healthcare Bioscience, USA) in a 2:1 ratio. After 30-min centrifugation at 600 g, the upper thin white layer of PBMCs was transferred to a centrifuge tube and resuspended in PBS for further use.
Electrophoretic Mobility Shift Assay (EMSA)
Nuclear protein was extracted from the PBMCs using the NE-PER Nuclear and Cytoplasmic Extraction Reagent Kit (ThermoFisher Scientific, USA). Biotinylated oligonucleotide probes were designed to correspond to the genomic sequences surrounding the candidate functional variants (8) (Table 1). EMSA was performed using a LightShift Chemiluminescent EMSA Kit (Thermo Scientific, USA) according to the manufacturer's instruction. Briefly, the binding reaction mixtures were set up containing 10 μg nuclear extract, 10× binding buffer (100 mM Tris, 500 mM KCl, 10 mM dithiothreitol), 50% glycerol, 100 mM MgCl2, poly(dI:dC) [poly(deoxyinosinic-deoxycytidylic) acid], and 1% Nonidet P-40 of up to 20 μl. The mixtures were incubated on ice for 15 min before the probe was added and were incubated for another 20 min at room temperature. For the competition assay, 1 μl of 100-fold excess amounts of unlabeled cold oligonucleotides against the biotinylated oligonucleotide probes used were incubated with the nuclear extract. Then, the DNA-nuclear protein complexes were run on 5% polyacrylamide gels at 50 V for 30 min and transferred to a nylon membrane at 200 mA for 30 min. After incubating the membrane with an antibody, the membrane was exposed and scanned using an LAS-3000 Imaging System (GE Healthcare Bioscience, USA). The experiment was repeated in triplicate.
Table 1. Oligonucleotide probes for EMSA.
Gene Expression
Total RNA from PBMCs were extracted with TRIzol Reagent (Invitrogen, USA) following the manufacturer's instruction. To ensure that the total RNA was free of DNA, the treatment with DNaseI was applied. The quantifications of CXCR3 (C-X-C motif chemokine receptor 3), CCR4, CCR6, IL-17A (interleukin-17A), INFG (interferon-γ), Pou2f1 (POU class 2 homeobox 1), and Sox_10 (SRY-box transcription factor 10) expression of PBMCs were tested by RNA microarrays, which was performed by Compass Biotechnology Co., Ltd (Beijing, China) using the HumanHT-12 v4.0 Gene Expression BeadChip (Illumina, USA).
Flow Cytometry
The CCR6 expression of lymphocytes and T-cell subtypes were detected by flow cytometry. Briefly, 500 μl of peripheral blood samples were incubated with antibodies against surface markers (anti-CD3 FITC, anti-CD4 APC, anti-CD8 PerCP-Cy5.5, and anti-CCR6 PE) at 4°C in the dark for 30 min. After erythrocyte lysis, the samples were washed with PBS and acquired by the BD FACSCanto II Flow Cytometry System (BD Biosciences, USA). About 10,000 lymphocytes were measured for each experiment. All manipulations were performed within 2 h after fresh peripheral blood was collected to avoid much death cells. The data were analyzed using FlowJo (Tree Star, USA); Supplementary Figure 2 shows the details of the gating strategy.
Statistical Analysis
Quality control for the discovery cohorts was performed using PLINK (http://www.cog-genomics.org/plink2/). We excluded samples with (i) sex inconsistencies, (ii) missing genotype call rate > 5%, and (iii) related or cryptically related individuals and single-nucleotide polymorphisms (SNPs) with (i) the missing call rate > 5%, (ii) minor allele frequency <1%, and (iii) deviation from the Hardy–Weinberg equilibrium. For rare variant quality control, minor allele frequency <1% was not applied. Imputation was performed using IMPUTE2. The reference panel comprised data from 1000 Genomes Project Phase 3 haplotypes (Nov 2014, build 37) (28). Supplementary Table 1 shows the details of the quality control and imputation.
Allelic and genotypic associations were assessed using PLINK to yield odds ratio (OR) with 95% confidence interval (95% CI), and p-value was adjusted for the false discovery rate method (29). Meta-analysis was performed with STATA (Stata Corp, USA) using the Mantel–Haenszel test, and the pooled OR was calculated using a fixed-effects model. Statistical power was calculated using Power and Sample Size Calculation (http://biostat.mc.vanderbilt.edu/PowerSampleSize). Quantitative variables with normal distribution were expressed as means and standard deviations, and Student's t-test was performed. Quantitative variables with non-normal distribution were expressed as medians and interquartile ranges, and Mann–Whitney U test was performed. Qualitative variables were summarized as absolute frequencies and percentages, and chi-square test was performed. The association between two continuous variables was determined using Pearson's correlation and linear regression analyses. The results were analyzed using SPSS 22.0 (SPSS Inc., USA). A two-tailed P < 0.05 was considered statistically significant.
Results
Identification of CCR6 Functional Variants
After adjusting for the false discovery rate, common variant-based meta-analysis of cohort 1 and cohort 2 identified 68 significantly associated regulatory variants, and rs3093018 was the top associated SNP (OR = 0.85, P = 6.88 × 10−4) (Supplementary Figure 3, Supplementary Table 2). The rare variant-based association analysis of cohort 3 failed to detect significantly associated rare variants. However, the exome chip replicated the associations of rs3093023 and rs3093024, which were among the 68 regulatory variants (rs3093023, P = 0.009; rs3093024, P = 0.010) (Supplementary Table 3, Supplementary Figure 4). In silico functional annotation showed that eight SNPs (rs3093023, rs1331299, rs975822, rs6907666, rs968333, rs6908364, rs6456158, and rs3093026) were likely to affect transcription factor binding and gene expression, while seven SNPs (rs3093023, rs968334, rs1571878, rs3093024, rs3093025, rs4709148, and rs10946216) showed the eQTL effects (Supplementary Table 4). Accordingly, rs3093023, with both statistical association and the eQTL effect, was prioritized as the candidate functional SNP. Consistently, linkage disequilibrium (LD) analysis showed that the other SNPs with the eQTL effect were in strong LD (r2 > 0.8) with rs3093023 (Supplementary Table 5).
Replication and Meta-Analysis
We validated the genetic association between rs3093023 and IgAN in an independent replication cohort. rs3093024, which was in almost absolute LD with rs3093023 (r2 = 0.98), was used as the positive control. rs3093023 was associated with IgAN susceptibility in the allele model (OR = 1.18, 95% CI: 1.05–1.33, P = 5.56 × 10−3), the co-dominant model (OR = 1.37, 95% CI: 1.08–1.74, P = 1.00 × 10−2), and the dominant model (OR = 1.26, 95% CI: 1.06–1.51, P = 9.22 × 10−3) (Table 2). Similar results were observed for rs3093024. The meta-analysis reinforced the genetic association of rs3093023 (OR = 1.17, 95% CI: 1.10–1.25, P = 6.14 × 10−7) (Figure 2A). Furthermore, the frequency of the rs3093023 A risk allele was high in east Asians, followed by Europeans and Americans; its frequency was lowest in Africans (Figure 2B), which is consistent with the epidemiological trend of IgAN (30).
Table 2. Allele and genotype frequencies of rs3093023 and rs3093024 in the replication cohort.
Figure 2. Meta-analysis of rs3093023. (A) The combined analysis with the discovery and replication cohorts. (B) Racial distribution of the rs3093023 alleles. OR, odds ratio; 95% CI, 95% confidence interval.
Functional Analysis of rs3093023
The distribution of the CCR6DNP genotypes in 361 people was as follows: TG/TG (19.94%), TG/CG (40.17%), TG/CA (7.48%), CG/CG (25.21%), CG/CA (6.92%), and CA/CA (0.28%). The association analysis showed that the rs3093023 A allele could tag the CCR6DNP TG allele in the entire genotyped cohort (cases, r = 0.87, P = 6.22 × 10−53; controls, r = 0.83, P = 2.13 × 10−51; all individuals, r = 0.85, P = 6.15 × 10−103) (Figure 3A, Supplementary Table 6).
Figure 3. Functional annotation of rs3093023. (A) The association between rs3093023 and CCR6DNP. (B) EMSA. (C) Prediction of transcription factors. (D) Analysis of expression association. EMSA, electrophoretic mobility shift assay.
We performed EMSA to investigate whether rs3093023 or CCR6DNP was the functional variant in IgAN. The oligonucleotide probe containing the rs3093023 allele showed significant mobility shifts, and the protective allele G bound much more protein than the risk allele A (Figure 3B, lanes 1 and 2). Moreover, the excessive unlabeled competitor probes abolished the DNA-nuclear protein complexes (Figure 3B, lanes 3 and 4). However, the oligonucleotide probe containing the CCR6DNP allele showed negligible mobility shifts (Figure 3B, lanes 5–7), and the excessive competitor probes did not abolish the DNA-nuclear protein complexes (Figure 3B, lanes 8–10).
In the RegulomeDB and JASPAR databases, we observed a higher affinity of Pou2f1 and Sox_10 for the regulatory motif for the rs3093023 G allele relative to the rs3093023 A allele (Figure 3C). However, the CCR6 mRNA expression was negatively associated with the Pou2f1 expression (r = −0.60, P = 3.94 × 10−9) rather than with the Sox_10 expression (r = 0.06, P = 0.578) (Figure 3D), suggesting that Pou2f1 may be the key transcription factor regulating the CCR6 expression via rs3093023.
CCR6 Expression Is Increased in Individuals With the Risk Genotype
We further used the HaploReg database to annotate the regulatory effects of rs3093023 and its proxy SNPs (rs968334, rs3093025, rs61056617, rs3093024, rs1854853, rs200505068, rs10946216, rs3093019, rs3093018, rs1571878, and rs3093017) in 127 tissues/cells. The results showed that the regulatory effects of rs3093023 and its proxies were significantly enriched in 23 tissues/cells (p < 0.05) (Supplementary Table 7). Among the 23 tissues/cells, T cells (12/23, 52.15%), mainly T-helper (Th) cells (6/23, 26.09%), accounted for the most in Figure 4, which suggests that rs3093023 may participate in the pathogenesis of IgAN by regulating T cells.
Figure 4. Tissue enrichment analysis of regulatory effects of rs3093023 and its proxies. (A) The significantly enriched tissues/cells (p < 0.05) of the regulatory effects of rs3093023 and its proxies. (B) The distribution of T cell subset enriched for regulatory effects of rs3093023 and its proxies. SNPs, single-nucleotide polymorphisms.
At mRNA level, individuals with the risk genotype A/A had significantly higher CCR6 mRNA levels compared with individuals with the protective genotype G/G, and these levels were dose-dependent on the number of A alleles (A/A vs. A/G vs. G/G: 7.91 ± 0.29 vs. 7.76 ± 0.45 vs. 7.51 ± 0.51, respectively, P = 4.39 × 10−2) (Figure 5A). The CCR6 mRNA levels were also higher in patients with IgAN than in controls (7.84 ± 0.48 vs. 7.58 ± 0.40, P = 2.64 × 10−2) (Figure 5B). Moreover, the CCR6 mRNA expression was positively associated with the expression of CCR4 (r = 0.63, P = 2.9 × 10−10) and IL-17A (r = 0.40, P = 2.57 × 10−4) but not with the expression of CXCR3 (r = 0.09, P = 0.450) and IFNG (r = 0.18, P = 0.112) (Figure 5C), suggesting that rs3093023 may regulate the CCR6 mRNA expression in Th17 cells.
Figure 5. CCR6 expression in individuals with rs3093023 genotypes. (A) The CCR6 mRNA expression in individuals with different rs3093023 genotypes. (B) The CCR6 mRNA expression in IgAN patients and controls. (C) The association between CCR6 mRNA expression and other cell surface markers of Th1 and Th17 cells. (D) The CCR6 protein expression in individuals with different rs3093023 genotypes. (E) The CCR6 protein expression in IgAN patients and controls.
At the protein level, we analyzed the associations between CCR6+ T cell subsets and rs3093023 genotypes. Consistently, higher CCR6+ cell ratios were observed in a risk allele dose-dependent manner (A/A > A/G > G/G) in lymphocytes (26.07 ± 7.23 vs. 24.41 ± 6.57 vs. 23.07 ± 6.54, respectively, P = 3.57 × 10−2), CD3+ T cells (15.96 ± 4.84 vs. 14.37 ± 4.34 vs. 14.21 ± 4.04, respectively, P = 4.54 × 10−2), and CD4+ T cells (11.56 ± 3.64 vs. 10.35 ± 3.56 vs. 9.81 ± 3.11, respectively, P = 1.32 × 10−2), but not in CD8+ T cells (2.98 ± 2.32 vs. 2.64 ± 1.69 vs. 2.51 ± 1.88, respectively, P = 0.388) (Figure 5D). A higher CCR6+ cell ratio was also observed in patients with IgAN for lymphocytes (28.31 ± 7.10 vs. 21.98 ± 5.47, P = 1.36 × 10−12), CD3+ T cells (16.66 ± 4.61 vs. 13.52 ± 3.43, P = 2.02 × 10−8), and CD4+ T cells (12.39 ± 3.46 vs. 9.21 ± 2.58, P = 2.25 × 10−13) (Figure 5E).
Association Between rs3093023 Genotypes With IgAN Phenotypes
DBP and Scr were significantly higher in patients with the risk genotype A/A than in patients with the protective genotype G/G (Table 3), and the increases were dose-dependent on the number of risk allele A (80.73 ± 13.04 mmHg vs. 78.85 ± 11.55 mmHg vs. 78.68 ± 11.64 mmHg, respectively, P = 0.023; 116.03 ± 83.67 μmol/L vs. 109.50 ± 76.23 μmol/L vs. 101.20 ± 52.02 μmol/L, respectively, P = 0.006). Moreover, the incidence of tubular atrophy/interstitial fibrosis increased with the number of risk allele A (47.09% vs. 46.17% vs. 38.03%, P = 0.034).
Table 3. Association of rs3093023 genotypes with IgAN phenotypes.
Discussion
In the discovery stage, we detected 68 common variants in CCR6 associated with IgAN susceptibility. rs3093023 showed both statistical significance and functional evidence, which was independently replicated. EMSA showed differential affinity of nuclear proteins for the rs3093023 alleles (G > A). In silico analysis identified Pou2f1 as a potential transcription factor, which was negatively associated with the CCR6 expression. Moreover, CCR6 transcript- and protein- expression levels increased with the number of risk allele A, and the CCR6 mRNA expression was positively associated with the expression of CCR4 and IL-17A rather than with the expression of CXCR3 and INFG. The clinical–pathological analysis showed that DBP, Scr, and incidence of tubular atrophy/interstitial fibrosis increased in a risk allele dose-dependent manner. Taken together, the results suggest that CCR6 rs3093023 may be involved in IgAN susceptibility by regulating Th17 cells.
CCR6 is a shared susceptibility gene between IgAN and other immune-related diseases (6). However, the associated CCR6 variants have not been reported in IgAN. A rheumatoid arthritis GWAS in Japanese identified CCR6 rs3093024 as the causal variant associated with the CCR6 expression (6). rs3093023 A causes disease susceptibility in patients with systemic sclerosis (9). Furthermore, the functional CCR6DNP affects CCR6 transcription in rheumatoid arthritis, systemic sclerosis, and Crohn's diseases (6, 8, 9). In this study, we first investigated the associations between IgAN and the CCR6 common and rare variants. Consistently, most of the CCR6 susceptibility variants associated with complex immune-related diseases are in the non-coding regions (31). A common variant, rs3093023, was prioritized as a functional variant. Although rs3093023 was in strong LD with the reported functional variant CCR6DNP (6), the allelic differences in nuclear protein binding were only observed for rs3093023, suggesting that rs3093023 may be the functional variant involved in IgAN.
Furthermore, our expression analysis indicated that rs3093023 and its proxy SNPs mainly alters the CCR6 expression in Th cells. Accumulating evidence indicates that changes in CD4+ T cells, including T regulatory (Treg) and effector Th cells, are involved in IgAN pathogenesis (18). Th2 polarization and Th17- and Tfh (T follicular helper)-derived interleukin release may promote aberrant IgA1 galactosylation in IgAN (32, 33). CCR6 is a surface marker for IL-17-producing Th17 cells (12) and IFNγ-producing Th1 cells (34). The expression of CCR6 and CCR4 together can identify Th17 cells, whereas the expression of CCR6 and CXCR3 identifies Th1 cells (34). In the present study, CCR6 transcript and protein expression levels increased with the number of risk allele. More importantly, CCR6 mRNA was positively associated with the expression of CCR4 and IL-17A rather than with the expression of CXCR3 and IFNG. Moreover, the significant association between rs3093023 genotypes and CCR6+ cells was mainly observed in CD4+ T cells. The results indicate that rs3093023 likely modulates the CCR6 expression in Th17 cells. Increased circulatory Th17 cells have been found in patients with IgAN and a murine model of IgAN (19, 35). The levels of the cytokines secreted by Th17 cells, including CCL20, IL-17A, IL-6, and IL-21, were also increased in kidneys of mice with IgAN (20). IgA1 could promote CCL20 production in mesangial cells from patients with IgAN to attract inflammatory CCR6+ Th17 cells to kidneys to accelerate renal injury (21). Notably, B cells have high amounts of CCR6 mRNA and protein (36). CCR6 and the corresponding ligand CCL20 can recruit B cells to renal inflammatory sites (37). However, the present study mainly focused on T cells with the most enrichment of rs3093023, and future studies can evaluate the role of B cells. Taken together, our findings provide further evidence from the genetic point of view that Th17 cells are involved in the pathogenesis of IgAN.
Finally, we observed that DBP, Scr, and the incidence of tubular atrophy/interstitial fibrosis were increased in a risk allele dose-dependent manner, indicating that the functional variant rs3093023 may be associated with IgAN severity. Higher Th17 cell numbers promote IL-17 production, inducing ongoing inflammatory response and renal injury (38). IL-17A activates human proximal renal tubular cells, which subsequently release high amounts of IL-6 and IL-8 (39). In the murine models of nephrotoxic nephritis (40) and ischemia/reperfusion (41), IL-17 deficiency ameliorated glomerular crescent formation and tubulointerstitial injury. In IgAN, serum levels of IL-17A secreted from Th17 cells were increased, and were positively associated with proteinuria (19). Immunohistochemistry has shown that IL-17A is expressed at the renal tubule site in most patients with IgAN (19). Compared with patients without IL-17A expression, patients with IL-17A expression have higher Scr, increased proteinuria, and more severe tubulointerstitial damage (19), which are strongly associated with progression to renal failure (42).
Overall, rs3093023 was prioritized as the functional variant in CCR6, which may contribute to IgAN susceptibility by regulating Th17 cells, suggesting that target Th17 cells may be a promising treatment for IgAN. Intriguingly, therapies that target Th17-related cytokines, Th17 intracellular signaling pathways, and Th17-specific transcription factors are available and have been reported to be efficient in the treatment of several immune-related diseases. For example, IL-17 neutralizing monoclonal antibodies, such as secukinumab, ixekizumab, and brodalumab, have been approved for treating psoriasis, psoriatic arthritis, and ankylosing spondylitis (43–47). A phase I randomized controlled trial revealed the potential benefit of GSK2981278, a RORγt inverse agonist, which significantly inhibits the production of Th17-producing cytokines in psoriasis (48). Although there were no studies reporting the efficacy of Th17 blockade in IgAN, two recent cases of IgAN concurrent with psoriasis demonstrated a decrease in proteinuria after the treatment of secukinumab, implying the potential beneficial effect of Th17 blockade in IgAN (49, 50). Collectively, further clinical trials are needed to verify the efficacy of targeting Th17 cells in the treatment of IgAN.
Our study has several limitations. First, as the genetic association was based on a northern Han Chinese population, multiracial replication studies should be conducted in the future. Second, we failed to identify rare variants in CCR6, and the effect sizes of the susceptibility locus accounted for <10% disease risk; therefore, additional genetic and environmental risk factors should be investigated. Third, the regulation effects of rs3093023 were mainly derived from in silico predictions, more reliable experiments especially in Th17 cells are needed for validation in the future.
In conclusion, through a two-stage genetic association study, we found that the CCR6 functional variant rs3093023 is involved in IgAN susceptibility via the regulation of Th17 cells. Our study also provides genetic evidence that Th17 cells are involved in IgAN and may be useful for understanding the pathogenesis and treatment of the disease.
Data Availability Statement
The original contributions presented in the study are included in the article/Supplementary Material, further inquiries can be directed to the corresponding author/s.
Ethics Statement
The studies involving human participants were reviewed and approved by Medical Ethics Committee of Peking University First Hospital (No. [2019]320). The patients/participants provided their written informed consent to participate in this study.
Author Contributions
YZ, XZ, and HZ conceived the study, generated the original hypothesis, and designed the experiments. YZ performed the experiments. YZ and XL conducted data analysis and drafted the manuscript. YZ, XL, XZ, LL, SS, PH, JL, and HZ recruited and prepared the samples. YZ, and HZ provided financial support. All authors reviewed/edited the manuscript and approved the final version.
Funding
This work was supported by National Natural Science Foundation of China [81800636 and 82070733], and CAMS Innovation Fund for Medical Sciences [2019-I2M-5-046].
Conflict of Interest
The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.
Acknowledgments
The authors would like to thank all patients and healthy controls who participated in this study.
Supplementary Material
The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fimmu.2021.600598/full#supplementary-material
Supplementary Figure 1. Genotyping of CCR6DNP.
Supplementary Figure 2. Flow cytometry gating strategy for CCR6+ cells.
Supplementary Figure 3. Genetic association of CCR6 common variants and IgAN. IgAN, immunoglobulin A nephropathy.
Supplementary Figure 4. Genetic association of CCR6 rare variants and IgAN.
Supplementary Table 1. Details of genotyping, quality control, and imputation for each cohort.
Supplementary Table 2. Associated regulatory variants after false discovery rate adjustment in cohort 1 and cohort 2.
Supplementary Table 3. Association of IgAN susceptibility and CCR6 rare variants in cohort 3. IgAN, immunoglobulin A nephropathy.
Supplementary Table 4. Functional annotations of 68 associated regulatory variants.
Supplementary Table 5. LD analysis of rs3093023 and other CCR6 SNPs with eQTL effects. eQLT, expression quantitative trait loci.
Supplementary Table 6. Association analysis of rs3093023 and CCR6DNP.
Supplementary Table 7. Tissues/cells enrichment for regulatory effects of rs3093023 and its proxy SNPs. SNPs, single-nucleotide polymorphisms.
Abbreviations
CCL20, C-C motif chemokine ligand 20; CCR4, C-C chemokine receptor 4; CCR6, C-C chemokine receptor 6; CCR6DNP, CCR6 dinucleotide polymorphism; CXCR3, C-X-C motif chemokine receptor 3; DBP, diastolic blood pressure; eGFR, estimated glomerular filtration rate; EMSA, electrophoretic mobility shift assay; eQTL, expression quantitative trait loci; GWAS, genome-wide association study; IgAN, IgA nephropathy; IL-17A, interleukin-17A; INFG, interferon-γ; LD, linkage disequilibrium; OR, odds ratio; PBMCs, peripheral blood mononuclear cells; PBS, phosphate-buffered saline; Pou2f1, POU class 2 homeobox 1; SBP, systolic blood pressure; Scr, serum creatinine; Sox_10, SRY-box transcription factor 10; Tfh, T follicular helper; Th, T helper.
References
1. Kiryluk K, Li Y, Sanna-Cherchi S, Rohanizadegan M, Suzuki H, Eitner F, et al. Geographic differences in genetic susceptibility to IgA nephropathy: GWAS replication study and geospatial risk analysis. PLoS Genet. (2012) 6:e1002765. doi: 10.1371/journal.pgen.1002765
3. Feehally J, Farrall M, Boland A, Gale DP, Gut I, Heath S, et al. HLA has strongest association with IgA nephropathy in genome-wide analysis. J Am Soc Nephrol. (2010) 10:1791–7. doi: 10.1681/ASN.2010010076
4. Gharavi AG, Kiryluk K, Choi M, Li Y, Hou P, Xie J, et al. Genome-wide association study identifies susceptibility loci for IgA nephropathy. Nat Genet. (2011) 4:321–7. doi: 10.1038/ng.787
5. Yu XQ, Li M, Zhang H, Low HQ, Wei X, Wang JQ, et al. A genome-wide association study in Han Chinese identifies multiple susceptibility loci for IgA nephropathy. Nat Genet. (2011) 2:178–82. doi: 10.1038/ng.1047
6. Kiryluk K, Li Y, Scolari F, Sanna-Cherchi S, Choi M, Verbitsky M, et al. Discovery of new risk loci for IgA nephropathy implicates genes involved in immunity against intestinal pathogens. Nat Genet. (2014) 11:1187–96. doi: 10.1038/ng.3118
7. Li M, Foo JN, Wang JQ, Low HQ, Tang XQ, Toh KY, et al. Identification of new susceptibility loci for IgA nephropathy in Han Chinese. Nat Commun. (2015) 6:7270. doi: 10.1038/ncomms8270
8. Kochi Y, Okada Y, Suzuki A, Ikari K, Terao C, Takahashi A, et al. A regulatory variant in CCR6 is associated with rheumatoid arthritis susceptibility. Nat Genet. (2010) 6:515–9. doi: 10.1038/ng.583
9. Koumakis E, Bouaziz M, Dieude P, Ruiz B, Riemekasten G, Airo P, et al. A regulatory variant in CCR6 is associated with susceptibility to antitopoisomerase-positive systemic sclerosis. Arthritis Rheum. (2013) 12:3202–8. doi: 10.1002/art.38136
10. Akhtar M, Jamal T, Jamal H, Din JU, Jamal M, Arif M, et al. Identification of most damaging nsSNPs in human CCR6 gene: in silico analyses. Int J Immunogenet. (2019) 6:459–71. doi: 10.1111/iji.12449
11. Zhou XJ, Mu R, Li C, Nath SK, Zhang YM, Qi YY, et al. Association of variants in CCR6 with susceptibility to lupus nephritis in Chinese. Arthritis Rheumatol. (2015) 11:3091–3. doi: 10.1002/art.39268
12. Schutyser E, Struyf S Van Damme J. The CC chemokine CCL20 and its receptor CCR6. Cytok Growth Factor Rev. (2003) 5:409–26. doi: 10.1016/S1359-6101(03)00049-2
13. Koga T, Otomo K, Mizui M, Yoshida N, Umeda M, Ichinose K, et al. Calcium/calmodulin-dependent kinase IV facilitates the recruitment of interleukin-17-producing cells to target organs through the CCR6/CCL20 Axis in Th17 cell-driven inflammatory diseases. Arthritis Rheumatol. (2016) 8:1981–8. doi: 10.1002/art.39665
14. Hirota K, Yoshitomi H, Hashimoto M, Maeda S, Teradaira S, Sugimoto N, et al. Preferential recruitment of CCR6-expressing Th17 cells to inflamed joints via CCL20 in rheumatoid arthritis and its animal model. J Exp Med. (2007) 12:2803–12. doi: 10.1084/jem.20071397
15. Furue K, Ito T, Tsuji G, Nakahara T, Furue M. The CCL20 and CCR6 axis in psoriasis. Scand J Immunol. (2020) 3:e12846. doi: 10.1111/sji.12846
16. Brkic Z, Corneth OBJ, van Helden-Meeuwsen CG, Dolhain RJEM, Maria NI, Paulissen SMJ, et al. T-helper 17 cell cytokines and interferon type I: partners in crime in systemic lupus erythematosus? Arthritis Res Ther. (2014) 3:410. doi: 10.1186/ar4578
17. Krebs CF, Paust HJ, Krohn S, Koyro T, Brix SR, Riedel JH, et al. Autoimmune renal disease is exacerbated by S1P-receptor-1-dependent intestinal Th17 cell migration to the kidney. Immunity. (2016) 5:1078–92. doi: 10.1016/j.immuni.2016.10.020
18. Ruszkowski J, Lisowska KA, Pindel M, Heleniak Z, Debska-Slizieowski J, Witkowski JM. T cells in IgA nephropathy: role in pathogenesis, clinical significance and potential therapeutic target. Clin Exp Nephrol. (2019) 3:291–303. doi: 10.1007/s10157-018-1665-0
19. Lin FJ, Jiang GR, Shan JP, Zhu C, Zou J, Wu XR. Imbalance of regulatory T cells to Th17 cells in IgA nephropathy. Scand J Clin Lab Invest. (2012) 3:221–9. doi: 10.3109/00365513.2011.652158
20. Meng T, Li X, Ao X, Zhong Y, Tang R, Peng W, et al. Hemolytic Streptococcus may exacerbate kidney damage in IgA nephropathy through CCL20 response to the effect of Th17 cells. PLoS ONE. (2014) 9:e108723. doi: 10.1371/journal.pone.0108723
21. Lu G, Zhang X, Shen L, Qiao Q, Li Y, Sun J, et al. CCL20 secreted from IgA1-stimulated human mesangial cells recruits inflammatory Th17 cells in IgA nephropathy. PLoS ONE. (2017) 5:e0178352. doi: 10.1371/journal.pone.0178352
22. Levey AS, Stevens LA, Schmid CH, Zhang YL, Castro AF III, Feldman HI, et al. A new equation to estimate glomerular filtration rate. Ann Intern Med. (2009) 9:604–12. doi: 10.7326/0003-4819-150-9-200905050-00006
23. Working Group of the International Ig ANN, the Renal Pathology S, Cattran DC, Coppo R, Cook HT, Feehally J, et al. The Oxford classification of IgA nephropathy: rationale, clinicopathological correlations, and classification. Kidney Int. (2009) 5:534–45. doi: 10.1038/ki.2009.243
24. Boyle AP, Hong EL, Hariharan M, Cheng Y, Schaub MA, Kasowski M, et al. Annotation of functional variation in personal genomes using RegulomeDB. Genome Res. (2012) 9:1790–7. doi: 10.1101/gr.137323.112
25. Westra HJ, Peters MJ, Esko T, Yaghootkar H, Schurmann C, Kettunen J, et al. Systematic identification of trans eQTLs as putative drivers of known disease associations. Nat Genet. (2013) 10:1238–43. doi: 10.1038/ng.2756
26. Mathelier A, Zhao X, Zhang AW, Parcy F, Worsley-Hunt R, Arenillas DJ, et al. JASPAR 2014: an extensively expanded and updated open-access database of transcription factor binding profiles. Nucleic Acids Res. (2014) 42:D142–7. doi: 10.1093/nar/gkt997
27. Ward LD Kellis M. HaploReg: a resource for exploring chromatin states, conservation, and regulatory motif alterations within sets of genetically linked variants. Nucleic Acids Res. (2012) 40:D930–4. doi: 10.1093/nar/gkr917
28. Delaneau O, Marchini J, Genomes Project C. Integrating sequence and array data to create an improved 1000 Genomes Project haplotype reference panel. Nat Commun. (2014) 5:3934. doi: 10.1038/ncomms4934
29. Hsueh H-m, Chen JJ, Kodell RL. Comparison of methods for estimating the number of true null hypotheses in multiplicity testing. J Biopharm Stat. (2003) 4:675–89. doi: 10.1081/BIP-120024202
30. Wyatt RJ, Julian BA. IgA nephropathy. N Engl J Med. (2013) 25:2402–14. doi: 10.1056/NEJMra1206793
31. Zhang F, Lupski JR. Non-coding genetic variants in human disease. Hum Mol Genet. (2015) 1:R102–10. doi: 10.1093/hmg/ddv259
32. Chintalacharuvu SR, Yamashita M, Bagheri N, Blanchard TG, Nedrud JG, Lamm ME, et al. T cell cytokine polarity as a determinant of immunoglobulin A (IgA) glycosylation and the severity of experimental IgA nephropathy. Clin Exp Immunol. (2008) 3:456–62. doi: 10.1111/j.1365-2249.2008.03703.x
33. Lin J-R, Wen J, Zhang H, Wang L, Gou F-F, Yang M, et al. Interleukin-17 promotes the production of underglycosylated IgA1 in DAKIKI cells. Ren Fail. (2018) 1:60–7. doi: 10.1080/0886022X.2017.1419972
34. Acosta-Rodriguez EV, Rivino L, Geginat J, Jarrossay D, Gattorno M, Lanzavecchia A, et al. Surface phenotype and antigenic specificity of human interleukin 17-producing T helper memory cells. Nat Immunol. (2007) 6:639–46. doi: 10.1038/ni1467
35. Peng Z, Tian J, Cui X, Xian W, Sun H, Li E, et al. Increased number of Th22 cells and correlation with Th17 cells in peripheral blood of patients with IgA nephropathy. Hum Immunol. (2013) 12:1586–91. doi: 10.1016/j.humimm.2013.08.001
36. Honczarenko M, Glodek AM, Swierkowski M, Na I-K, Silberstein LE. Developmental stage-specific shift in responsiveness to chemokines during human B-cell development. Exp Hematol. (2006) 8:1093–100. doi: 10.1016/j.exphem.2006.05.013
37. Lee AYS, Körner H. The CCR6-CCL20 axis in humoral immunity and T-B cell immunobiology. Immunobiology. (2019) 3:449–54. doi: 10.1016/j.imbio.2019.01.005
38. Woltman AM, de Haij S, Boonstra JG, Gobin SJ, Daha MR, van Kooten C. Interleukin-17 and CD40-ligand synergistically enhance cytokine and chemokine production by renal epithelial cells. J Am Soc Nephrol. (2000) 11:2044–55.
39. Van Kooten C, Boonstra JG, Paape ME, Fossiez F, Banchereau J, Lebecque S, et al. Interleukin-17 activates human renal epithelial cells in vitro and is expressed during renal allograft rejection. J Am Soc Nephrol. (1998) 8:1526–34.
40. Krohn S, Nies JF, Kapffer S, Schmidt T, Riedel JH, Kaffke A, et al. IL-17C/IL-17 receptor E signaling in CD4(+) T cells promotes TH17 cell-driven glomerular inflammation. J Am Soc Nephrol. (2018) 4:1210–22. doi: 10.1681/ASN.2017090949
41. Wang F, Yin J, Lin Y, Zhang F, Liu X, Zhang G, et al. IL-17C has a pathogenic role in kidney ischemia/reperfusion injury. Kidney Int. (2020) 97:1219–29. doi: 10.1016/j.kint.2020.01.015
42. Lv J, Shi S, Xu D, Zhang H, Troyanov S, Cattran DC, et al. Evaluation of the Oxford classification of IgA nephropathy: a systematic review and meta-analysis. Am J Kidney Dis. (2013) 5:891–9. doi: 10.1053/j.ajkd.2013.04.021
43. Hueber W, Patel DD, Dryja T, Wright AM, Koroleva I, Bruin G, et al. Effects of AIN457, a fully human antibody to interleukin-17A, on psoriasis, rheumatoid arthritis, and uveitis. Sci Transl Med. (2010) 2:52ra72. doi: 10.1126/scitranslmed.3001107
44. Papp KA, Langley RG, Sigurgeirsson B, Abe M, Baker DR, Konno P, et al. Efficacy and safety of secukinumab in the treatment of moderate-to-severe plaque psoriasis: a randomized, double-blind, placebo-controlled phase II dose-ranging study. Br J Dermatol. (2013) 168:412–21. doi: 10.1111/bjd.12110
45. Langley RG, Elewski BE, Lebwohl M, Reich K, Griffiths CE, Papp K, et al. Secukinumab in plaque psoriasis–results of two phase 3 trials. N Engl J Med. (2014) 371:326–38. doi: 10.1056/NEJMoa1314258
46. Blauvelt A, Prinz JC, Gottlieb AB, Kingo K, Sofen H, Ruer-Mulard M, et al. Secukinumab administration by pre-filled syringe: efficacy, safety and usability results from a randomized controlled trial in psoriasis (FEATURE). Br J Dermatol. (2015) 172:484–93. doi: 10.1111/bjd.13348
47. Farahnik B, Beroukhim K, Abrouk M, Nakamura M, Zhu TH, Singh R, et al. Brodalumab for the treatment of psoriasis: a review of phase III trials. Dermatol Ther. (2016) 6:111–24. doi: 10.1007/s13555-016-0121-x
48. Kang EG, Wu S, Gupta A, von Mackensen YL, Siemetzki H, Freudenberg JM, et al. A phase I randomized controlled trial to evaluate safety and clinical effect of topically applied GSK2981278 ointment in a psoriasis plaque test. Br J Dermatol. (2018) 178:1427–9. doi: 10.1111/bjd.16131
49. Segawa Y, Ishida R, Kanehisa F, Nakai K, Morimoto M, Seno M, et al. IgA nephropathy in a patient receiving infliximab for generalized pustular psoriasis. BMC Nephrol. (2020) 21:366. doi: 10.1186/s12882-020-02015-0
Keywords: IgA nephropathy, CCR6, Th17 cells, genetic association, functional annotation
Citation: Zhang Y-m, Liu X-z, Zhou X-j, Liu L-j, Shi S-f, Hou P, Lv J-c and Zhang H (2021) A Functional Variant rs3093023 in CCR6 Is Associated With IgA Nephropathy by Regulating Th17 Cells in a North Han Chinese Population. Front. Immunol. 12:600598. doi: 10.3389/fimmu.2021.600598
Received: 30 August 2020; Accepted: 01 February 2021;
Published: 25 February 2021.
Edited by:
Kazumitsu Sugiura, Fujita Health University, JapanReviewed by:
Andras Perl, Upstate Medical University, United StatesSun Jung Kim, Northwell Health, United States
Copyright © 2021 Zhang, Liu, Zhou, Liu, Shi, Hou, Lv and Zhang. This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner(s) are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.
*Correspondence: Hong Zhang, hongzh@bjmu.edu.cn
†These authors have contributed equally to this work