Functional genomic annotation of genetic risk loci highlights inflammation and epithelial biology networks in CKD

Abstract

Genome-wide association studies (GWASs) have identified multiple loci associated with the risk of CKD. Almost all risk variants are localized to the noncoding region of the genome; therefore, the role of these variants in CKD development is largely unknown. We hypothesized that polymorphisms alter transcription factor binding, thereby influencing the expression of nearby genes. Here, we examined the regulation of transcripts in the vicinity of CKD-associated polymorphisms in control and diseased human kidney samples and used systems biology approaches to identify potentially causal genes for prioritization. We interrogated the expression and regulation of 226 transcripts in the vicinity of 44 single nucleotide polymorphisms using RNA sequencing and gene expression arrays from 95 microdissected control and diseased tubule samples and 51 glomerular samples. Gene expression analysis from 41 tubule samples served for external validation. 92 transcripts in the tubule compartment and 34 transcripts in glomeruli showed statistically significant correlation with eGFR. Many novel genes, including ACSM2A/2B, FAM47E, and PLXDC1, were identified. We observed that the expression of multiple genes in the vicinity of any single CKD risk allele correlated with renal function, potentially indicating that genetic variants influence multiple transcripts. Network analysis of GFR-correlating transcripts highlighted two major clusters; a positive correlation with epithelial and vascular functions and an inverse correlation with inflammatory gene cluster. In summary, our functional genomics analysis highlighted novel genes and critical pathways associated with kidney function for future analysis.

Keywords: CKD; gene expression; genetics of complex trait.

Figure 1.

Correlation between CRAT expression in glomeruli and renal function. The y axis shows the relative normalized glomerular expressions of (A) FAM47E, (B) PLXDC1, (C) VEGFA, and (D) MAGI2. The x axes show the eGFR for each sample. Each dot represents one individual miscrodissected glomerular sample. The lines represent the fitted linear correlation values. Immunohistochemistry shows the protein expression in human glomeruli ([E] FAM47E, [F] PLXDC1, [G] VEGFA, and [H] MAGI2). Scale bars, 100 μm. Reprinted from www.proteinatlas.org.

Figure 2.

Correlation between CRAT expression in tubules and renal function. Expressions of (A) SLC34A1, (B) SLC7A9, and (C) ACSM5 correlate with eGFR in tubule samples. The x axes represent eGFR (ml/min per 1.73 m²), whereas the y axes represent the normalized gene expression values of the transcript. Each dot represents transcript levels and eGFR values from a single kidney sample. The lines are the fitted correlation values. Immunohistochemistry shows tubular-specific expression of (D) SLC34A1, (E) SLC7A9, and (F) ACSM5. Scale bars, 100 μm. Reprinted from www.proteinatlas.org.

Figure 3.

UMOD and ACSM2A expressions correlate with renal function. The expressions of (A) UMOD and (F) ACSM2A correlate with eGFR in tubule samples. The x axes represent eGFR (ml/min per 1.73 m²), whereas the y axes represent the normalized gene expression values of the transcript. Each dot represents transcript levels and eGFR values from a single kidney sample. The lines are the fitted correlation value (Pcorr, P value after Benjamini-Hochberg multiple testing correction). Immunohistochemistry of the samples with low and high mRNA expression showed differences of (B–E) the UMOD and (G–J) the ACSM2A expression on protein level. Scale bars, 50 μm.

Figure 4.

Tubule-specific transcript levels correlate with renal function near the UMOD locus (rs4293393, rs12917707, and rs11864909 polymorphisms) and the disabled homolog 2 (DAB2) locus (rs11959928). The x axes represent the genomic positions of each gene on (A and C) chromosomes 16 and (D and F) 5. The y axes represent the negative logarithms of the corrected P values (significances) between the expressions of each gene and eGFR (ml/min per 1.73 m²). (A and D) Color coding represents the baseline expression of the transcripts in human kidney on the basis of the RNA sequencing data. Red, high expression in the kidney; yellow, medium expression in the kidney; green, low expression in the kidney. (B and E) On the basis of the results of the Illumina Body Map (www.ebi.ac.uk), a heat map was generated from the FPKM values of the CRATs near these SNPs. High expression values (90th percentile) are marked red, and low expression values (<10th percentile) are marked blue. Expressions with FPKM values<0.1 are marked white. *Genes without probe set identifications on the Affymetrix arrays. QRT-PCR validation confirmed the significant correlation with eGFR of the following transcripts: (C) glycoprotein 2 (GP2), UMOD, ACSM5, ACSM2A, and ACSM2B and (F) FYN binding protein (FYB) and DAB2. A shows a strong correlation between UMOD expression and eGFR, whereas the expressions of ACSM5 and -2A/2B also highly correlate with renal function. (D) At the rs11959928 locus, not only the transcript DAB2 but also, the FYB show high correlation with eGFR (PDILT, Protein disulfide isomerase-like, testis expressed; C9, Complement component 9).

Figure 5.

The expression of VEGFA correlates with renal function. The expression of VEGFA is significantly lower (*P=0.025) in samples homozygous for A alleles (A/A; n=7) at the rs881585 locus compared with samples with minor alleles (A/G; n=7 or G/G; n=7) at this locus. (A) Only control samples (eGFR>85 ml/min per 1.73 m²) were used for the analysis. (B) Microarray-based transcript levels of VEGFA correlate with renal function in tubule samples (R²=0.219, P=1.7×10⁻⁶). (C) QRT-PCR–based VEGFA transcript levels (R²=0.228, P=7.8×10⁻⁴) confirm its correlation with kidney function. (D) VEGFA protein expression (by immunohistochemistry) correlates with transcript levels. Counterstained with hematoxylin. Scale bars, 50 μm.

Figure 6.

Kidney function-correlating CRATs form tight networks. (A) CRATs showing negative correlation with eGFR (green with P corrected<0.05) clustered around TNF and TGF-β. (B) CRATs showing positive correlation with eGFR (red with P corrected<0.05) centered around VEGFA and ERBB2 [erythroblastic leukemia viral oncogene homolog 2 (EGFR2, epidermal growth factor receptor 2)] (Ingenuity Systems).

Figure 7.

Schematic representation of the experimental design. GWASs examine the relationship between genetic variants (SNP) and disease state (CKD). The eQTL examines the relationship between transcript levels and genetic variation in control samples. Here, we investigated the relationship between transcript levels around CKD risk variants and kidney function by examining the contribution of genetic and environmental factors.