Genome-wide association study of kidney function decline in individuals of European descent.

Genome-wide association studies (GWASs) have identified multiple loci associated with cross-sectional eGFR, but a systematic genetic analysis of kidney function decline over time is missing. Here we conducted a GWAS meta-analysis among 63,558 participants of European descent, initially from 16 cohorts with serial kidney function measurements within the CKDGen Consortium, followed by independent replication among additional participants from 13 cohorts. In stage 1 GWAS meta-analysis, single-nucleotide polymorphisms (SNPs) at MEOX2, GALNT11, IL1RAP, NPPA, HPCAL1, and CDH23 showed the strongest associations for at least one trait, in addition to the known UMOD locus, which showed genome-wide significance with an annual change in eGFR. In stage 2 meta-analysis, the significant association at UMOD was replicated. Associations at GALNT11 with Rapid Decline (annual eGFR decline of 3 ml/min per 1.73 m(2) or more), and CDH23 with eGFR change among those with CKD showed significant suggestive evidence of replication. Combined stage 1 and 2 meta-analyses showed significance for UMOD, GALNT11, and CDH23. Morpholino knockdowns of galnt11 and cdh23 in zebrafish embryos each had signs of severe edema 72 h after gentamicin treatment compared with controls, but no gross morphological renal abnormalities before gentamicin administration. Thus, our results suggest a role in the deterioration of kidney function for the loci GALNT11 and CDH23, and show that the UMOD locus is significantly associated with kidney function decline.

Introduction

Chronic kidney disease (CKD) is an important public health problem affecting up to 10% of adults world-wide [1-3]. Faster rates of decline in estimated glomerular filtration rate (eGFR), and entry into CKD stages of increasing severity are associated with an increased risk of cardiovascular and all-cause mortality [4-9]. Thus, recently issued guidelines on the evaluation and management of patients with CKD have highlighted the importance of evaluating longitudinal measures of renal function in addition to determining eGFR and urinary albumin excretion at discrete time points [3].

Traditional risk factors for CKD include diabetes and hypertension, but these do not fully account for CKD risk [10]. There is evidence for considerable clustering of CKD within families [11] and the heritability of eGFR has been estimated at up to 36–75% in population-based studies [12]. Using genome-wide association studies (GWAS), multiple loci have been identified in association with eGFR and CKD in both European [13-16] and non-European populations [17,18] using data from one time point. However, multiple lines of evidence suggest that there may be unique genetic contributions to renal function decline above and beyond baseline renal function. First, there is substantial variability in the rate of eGFR decline in studies of healthy persons as well as among those with CKD [3,4,19,20]. Second, we have previously shown that some genetic loci associated with cross-sectional eGFR are also associated with incident CKD even after accounting for baseline eGFR [21]. Finally, genetic background has been shown to affect CKD progression in animal models [22,23].

Taken together, these data suggest that unique loci may exist for renal function decline in addition to those identified for a one-time measure of eGFR. Thus, we conducted a genome-wide association study (GWAS) meta-analysis among participants from 16 cohorts with serial kidney function measurements within the CKDGen Consortium, followed by independent replication among additional participants from 13 cohorts.

Results

Study participants

Changes in renal function over time were derived from 45,530 individuals who participated in stage 1 meta-analysis of study-specific GWAS, and an additional 18,028 independent individuals who participated in stage 2 meta-analysis (Table 1). Details on study design and genotyping are provided in Supplementary Tables 1 and 2 respectively.

Table 1

Stage 1 and Stage 2 cohort characteristics.

	n	Age atbaseline,years	% women,n	%HTN atbaseline,n	%DM atbaseline,n	%CKD atbaseline,n	eGFRbaseline,ml/min/1.73m2	eGFR follow-up,ml/min/1.73m2	Durationbetweenbaseline andfollow-up(Years)mean, SD
Stage 1 discovery cohorts
AGES	3219	54.2 (8.98)	58.0, 1867	21.6, 694	3.2, 103	3.2, 104	89.6 (19.3)	73.0 (20.0)	22.2 (6.7)
Amish	458	46.9 (14.3)	54.4, 249	8.7, 40	4.2, 19	NA	95.9 (24.3)	89.2 (19.4)	5.2 (2.6)
ARIC	9049	54.5 (5.7)	52.9, 4793	26.9, 2427	8.6, 780	2.9, 260	89.7 (17.0)	83.5 (17.1)A80.6 (17.1)B	8.0 (2.2)
ASPS	848	65.8 (8.1)	56.8, 482	72.5, 615	9.8, 83	12.7, 108	80.2 (20.2)	74.6 (15.1)	5.5 (2.0)
CHS	2820	71.9 (5.0)	61.3, 1729	34.5, 966	11.0, 310	7.9, 224	77.3 (20.8)	75.4 (20.2)	5.9 (1.8)
CoLaus	1918	53.9 (10.8)	54.9, 1053	35.4, 679	6.31, 121	5.0, 95	91.5 (20.4)	84.9 (18.2)	5.58 (0.29)
FHS (Offspring and Cohort)	2523	58.1 (8.6)	55.6, 1405	36.8, 927	8.1, 206	8.3, 210	88.7 (25.5)	79.5 (19.1)	11.1 (3.6)
GENOA	1041	54.7 (10.3)	55.6, 579	71.7, 746	8.8, 92	4.3, 45	92.1 (20.7)	89.8 (22.8)	4.0 (1.1)
HABC	888	73.4 (2.8)	48.5, 431	33.1, 294	8.5, 75	21.4, 190	71.8 (13.2)	72.9 (21.2)	9.0 (0.3)
JUPITER	8780	66.1 (7.8)	32.2, 2826	63.8, 5602	0.6, 54	11.5, 1008	80.1 (18.1)	78.2 (17.7)	2.6 (1.0)
KORA3	1641	52.5 (10.1)	49.5, 813	38.3, 629	4.3, 71	3.2, 53	91.3 (18.0)	83.9 (21.0)	10 (0)
KORA4	1807	53.8 (8.9)	51.3, 927	33.7, 606	3.7, 66	3.9, 70	89.5 (17.5)	85.1 (20.2)	7.1 (0.4)
MESA	2324	63.2 (10.1)	51.7, 1201	37.3, 866	5.5, 135	13.3, 310	74.2 (13.9)	70.7 (15.1)	4.8 (0.5)
The Rotterdam Study (RS-I)	2422	66.5 (7.0)	60.2, 1459	50.8, 1230	7.5, 182	7.7, 186	79.9 (15.5)	73.7 (15.8)	6.5 (0.4)
SHIP	3203	49.2 (15.3)	51.8, 1659	24.3, 778	7.0, 225	3.7, 119	92.4 (19.8)	90.6 (23.6)	5.3 (0.7)
Three Cities (3C)	2589	73.0 (4.5)	61.9, 1602	76.7, 1986	8.6, 223	18.9, 489	73.1 (16.1)	71.0 (16.8)	3.8 (0.3)
Stage 2 replication cohorts
ADVANCE	2034	67.0 (6.6)	31.9, 649	55.2, 1123	100, 2034	16.0, 325	84.1 (28.1)	84.8 (34.5)	4.9 (0.9)
BMES	1304	62.9 (7.7)	60.1, 784	67.1, 875	5.4, 71	23.9, 312	82.6 (31.7)	75.5 (34.9)	10.4 (0.6)
COLAUS	2238	53.1 (10.4)	53.9, 1207	24.0, 538	4.1, 91	3.5, 79	90.5 (19.5)	88.7 (18.7)	5.5 (0.3)
HYPERGENES	651	53.4 (7.5)	45.3, 295	13.9, 91	0	0.61, 4	107.4 (23.5)	103.4 (35.1)	5.6 ( 3.2)
KORA3	1494	51.6 (13.3)	52.5, 785	29.4, 437	5.1, 76	2.6, 39	98.0 (20.1)	92.4 (21.3)	9.6 (0.6)
KORA4	1200	49.2 (15.4)	52.4, 629	13.3, 159	4, 48	5.8, 70	97.4 (21.7)	92.6 (22.4)	7.2 (0.5)
NESDA	1270	43.3 (12.3)	67.2, 854	32.4, 411	5.4, 69	0.9, 11	97.8 (20.4)	95.5 (18.5)	2.0 (0.3)
popgen	577	60.2 (9.4)	42.1, 243	50.4, 288	5.2, 30	6.1, 35	84.5 (17.2)	79.9 (18.4)	4.8 (0.8)
PREVEND (4 year follow-up)	791	53.0 (13.3)	50.8, 402	40.8, 323	4.2, 33	6.0, 47	89.2 (19.5)	106.4 (34.2)	4.3 (0.6)
PREVEND (9 year follow-up)	2169	48.0 (11.1)	48.0, 1040	28.2, 612	3.1, 66	1.7, 37	93.7 (17.7)	86.9 (18.6)	9.4 (0.84)
RS-II	1243	61.8 (5.2)	54.6, 679	21.4, 266	7.8, 97	8.6, 186	81.2 (17.0)	73.7 (15.8)	10.6 (0.4)
SAPHIR	1374	51.6 (6.0)	39.0. 536	54.9, 754	2.6, 36	0.9, 13	91.5 (15.8)	88.0 (16.0)	4.6 (0.7)
YFS	1683	31.9 (4.9)	56.0, 943	8.3, 139	1.1, 18	0.2, 3	105.4 (16.4)	100.4 (16.0)	6.0 (0)

Unless indicated otherwise, values are given as mean (SD)

A: eGFR at visit 2; B: eGFR at visit 4; DM: Diabetes; HTN: hypertension; CKD: Chronic kidney disease (eGFR<60ml/min/1.73m2); eGFR: estimated glomerular filtration rate; SD: standard deviation

At the baseline examination, the prevalence of CKD, defined as eGFR<60 ml/min/1.73m2, ranged from 3.2% to 21.4% in stage 1 cohorts and from 0.2% to 23.9% in stage 2 replication cohorts. As expected, cohorts with lower mean age at baseline tended to have a lower baseline prevalence of CKD. Four kidney function decline traits were derived from serial eGFR values in each study participant to model mechanisms underlying different rates of kidney function change over time: 1) annual decline of eGFR (eGFRchange, in ml/min/1.73m2 decline per year; a positive value represents a decline in eGFR, whilst a negative value represents a rise in eGFR over time), 2) incident CKD to select individuals with a decline in kidney function to the clinical outcome CKD stage 3 or higher (CKDi, cases defined as those free of CKD at baseline but eGFR<60ml/min/1.73m2 during follow-up), 3) incident CKD with additionally at least a 25% eGFR decline from baseline to select individuals reaching CKD stage 3 after a sizeable decline in kidney function (CKDi25,) [24], and 4) rapid eGFR-decline to select individuals with the highest risk of adverse outcomes (Rapid Decline, cases defined as those with annual eGFR-decline≥3ml/min/1.73m2) [5]. Most cohorts showed a decline in kidney function over time (Table 1). The distribution of all four traits in stage 1 and stage 2 cohorts can be found in Supplementary Table 3.

Heritability of eGFR change

The heritability of eGFR change in the Framingham Heart Study was estimated as 38%, after adjusting for age, sex, and baseline eGFR.

Stage 1 meta-analysis of GWAS of measures of kidney function change over time

Stage 1 GWAS meta-analysis was performed in all samples for all four traits. Two secondary association analyses were performed to account for potentially different rates of kidney function decline in those with and without CKD: 1) eGFRchange stratified by baseline CKD status and 2) Rapid Decline in only those without baseline CKD; too few individuals with CKD fulfilled Rapid Decline criteria to perform this analysis. Supplementary Figure 1 shows the Manhattan and QQ-plots of the stage 1 meta-analysis of each trait. The genomic control factor ranged from 1.007 – 1.05, suggesting negligible evidence for population stratification.

In GWAS meta-analysis of stage 1 cohorts, the minor T allele of rs12917707 at the UMOD locus, previously identified by GWAS to be associated with higher eGFR in cross-sectional analysis [14], was associated with an increase in eGFR over time at a genome-wide significant level (p=2.6×10−14, Table 2), and showed at least nominally significant, direction consistent association with all other analyzed phenotypes (Supplementary Table 4). In addition, SNPs at the novel CDH23, GALNTL5/GALNT11, MEOX2, IL1RAP/OSTN, C2orf48/HPCAL1 and NPPB/NPPA loci were associated with at least one of the analyzed traits at a significance level of p<10−6 (Table 2). Thus, a total of 7 SNPs were moved forward to stage 2 meta-analysis. These SNPs mostly showed high imputation quality in each cohort or were genotyped de-novo (Supplementary Table 5), and low between-study heterogeneity (I2<25%).

Table 2

Genetic association results of SNPs identified in stage 1 meta-analysis

								discovery stage 1		Replicationstage 2		Stage 1 and stage 2combined
SNPID	trait	Chr	Position(b36)	Locus	codedallele	non-codedallele	AFcodedallele	beta	PVal2GC	beta	one-sidedpval1GC	beta	two-sidedpval1GC	totalsamplesize
rs12917707	eGFRchange overall	16	20275191	UMOD, PDILT	t	g	0.18	−0.15	2.6×10−14	−0.12	4.7×10−5	−0.14	1.8×10−17	59373
rs11803049	eGFRchange CKD	1	11851482	NPPB, NPPA, KIAA2013, CLCN6	a	g	0.07	−0.57	3.6×10−7	0.02	0.43*	−0.27	6.2×10−4	4116
rs875860	eGFRchange CKD	10	72979535	CDH23	t	c	0.12	−0.49	6.2×10−7	−0.15	0.047*	−0.31	4.6×10−6	4116
rs11764932	Rapid Decline overall	7	15699643	MEOX2	a	g	0.36	0.12	6.8×10−8	0.04	0.14	0.09	3.6×10−7	61078
rs1019173	Rapid Decline overall	7	151341480	GALNTL5, GALNT11, MLL3, CCT8L1	a	g	0.63	−0.12	3.0×10−7	−0.06	0.04	−0.10	2.1×10−7	61077
rs9814367	Rapid Decline noCKD	3	192075180	IL1RAP, OSTN	t	c	0.92	−0.20	4.1×10−7	0.02	0.39	−0.13	7.3×10−5	56687
rs759341	CKDi25	2	10297660	C2orf48, HPCAL1, RRM2	a	g	0.31	0.18	1.5×10−6	0.06	0.27$	0.16	2.7×10−6	41122

“Locus” is based on build 36, hg18. The gene closest to the SNP is listed first and is in boldface if the SNP is located within the gene.

studies included: ADVANCE, BMES, COLAUS, RS-II

studies included: ADVANCE, BMES, RS-II

Stage 2 meta-analysis

Of the seven loci moved forward for stage 2 meta-analysis, only rs12917707 at UMOD was significantly associated with the stage 1 trait after correcting for multiple testing (p=4.7*×10−5). Two further SNPs showed suggestive significance (one-sided p<0.05) with their respective stage 1 trait: rs875860 in CDH23 with eGFRchange in those with CKD at baseline, and rs1019173 at GALNTL5/GALNT11 with Rapid Decline (Table 2). There was no significant heterogeneity between studies for these two SNPs (rs875860: I2=9.7%, p=0.34; rs1019173: I2=32.4%, p=0.12) or for the other SNPs analyzed in stage 2 meta-analysis (I2 <30.0%).

The SNP rs1019173 is located in an intron in the GALNTL5 gene, and lies in a linkage disequilibrium (LD) block spanning the genes GALNT11, MLL3, CCT8L, and part of the GALNTL5 gene (Figure 1a). The SNP in CDH23, rs875860, is an intronic SNP in an LD block whose boundaries lie within the coding region of the CDH23 gene (Figure 1b).

Figure 1

Regional Association Plots of the novel loci identified by GWAS of kidney function decline traits. Negative log10 pvalues are plotted versus genomic position (build 36, hg18). The lead SNP in each region is labeled. Other SNPs in each region are color-coded based on their LD to the lead SNP. Light blue lines indicate recombination rate (cM/Mb). (A) GALNTL5/GALNT11 locus. (B) CDH23 locus.

In the combined meta-analysis of these three SNPs from both stage 1 and stage 2 cohorts, there was no evidence of between-study heterogeneity in the combined metaanalysis (I2<25%). Only the SNP at UMOD showed genome-wide significant association (rs12917707, p=1.2×10−16) in the combined stage 1 and stage 2 analysis, whereas there was suggestive evidence of significance for the two novel loci identified in stage 1 (rs875860 in CDH23: p=1.5×10−6 for the association with eGFRchange in those with CKD; rs1019173 at GALNTL5/GALNT11: OR=0.91 for the A allele, p=2.2×10−7 for the association with Rapid Decline).

Functional validation of novel loci in zebrafish

To investigate the role of the two suggestive novel loci in vertebrate kidney development and function and to bolster confidence in the nominally significant statistical associations in the replication studies, we knocked down the corresponding genes in the zebrafish using antisense morpholino (MO) technology. We focused on the CDH23 region and the block containing GALNTL5, GALNT11, MLL3 and CCT8L1. For the latter region, we focused on GALNT11 and MLL3, because there are no zebrafish GALNTL5 and CCT8L1 orthologs. Further, we investigated the effect of MO knockdown of umod. Following MO injection at the 1-cell stage, we performed in situ hybridization for the established renal markers pax2a (global kidney) and nephrin (podocytes) at 48 hours post-fertilization (hpf). Compared to control embryos, cdh23, galnt11, mll3a, mll3b and umod morphants did not display significant defects in glomerular or tubule gene expression (Figure 2A, n>25 embryos per MO injection).

Figure 2

Cdh23 and galnt11 knockdowns exacerbate nephrotoxic injury in zebrafish embryos

(A) Whole mount in situ hybridization for the global kidney marker pax2a (arrowhead denotes the glomerulus, bracket denotes the tubule) and the podocyte marker nephrin demonstrates that morpholino (MO) knockdowns of cdh23, galnt11, mll3a, and mll3b do not result in changes in kidney gene expression compared to control embryos at 48 hours post-fertilization (hpf). Similar results were obtained for MO knockdowns of umod (images not shown). (B) Morpholino knockdown of cdh23 and galnt11 causes embryos to develop edema at a higher frequency than control embryos following gentamicin challenge. Data are presented as number of observed abnormalities per total number of embryos scored at 5, 24, 48, and 72 hours post-gentamicin injection (hpi), normalized to control experiments. *p < 0.05, **p < 0.005 by Fisher’s exact test. (C) Graphical representation of edema prevalence in embryos injected with gentamicin in (B). (D) Control embryos develop minor (cardiac) edema whereas cdh23 and galnt11 MO-injected embryos develop severe (cardiac, intestinal, ocular) edema 72 hours after gentamicin injection.

It is possible that morphant embryos develop a kidney function decline phenotype only after exposure to a nephrotoxin, despite observing no differences in renal marker expression at 48 hpf. Accordingly, after MO injection, we injected embryos with gentamicin at 48 hpf and observed edema prevalence and severity over the next three days. In control embryos, gentamicin injection predictably resulted in a majority of embryos developing minor (cardiac) edema by 24 hours post-injection (hpi) (Figure 2B–D). In comparison, cdh23 and galnt11 morphants developed significantly more severe (cardiac, intestinal, and ocular) and more frequent edema (Figure 2B–D). Specifically, whereas 10% of control embryos developed severe edema by 72 hpi, 43% of cdh23 morphants (p=0.009) and 55% of galnt11 morphants (p=0.001) developed severe edema at this time point. Additionally, a significant proportion of cdh23 (33%, p=0.035) and galnt11 morphant embryos (46%, p=0.005) injected with gentamicin developed edema earlier compared to controls at 5 hpi. In contrast, knockdown of mll3 or umod affected neither kidney development nor susceptibility to gentamicin (Figure 2B–C). Taken together, these data demonstrate that knockdown of cdh23 and galnt11 results in altered renal function after a nephrotoxic insult.

Interrogation of novel loci in eSNP databases and the CRIC Study

We interrogated eSNP data bases for evidence of SNPs at the CDH23 and GALNTL5/GALNT11 loci to evaluate an effect on gene expression [25] but found no relevant associations. Similarly, annotation information provided by ANNOVAR [26] did not yield genetic variants of potential functional interest within 500kb of and in linkage disequilibrium (r2 > 0.8 based on HapMap release 22) with the index SNPs.

In Caucasian participants of the Chronic Renal Insufficiency Cohort (CRIC) study, a prospective study of patients with CKD at baseline [27], neither SNPs in GALNTL5/GALNT11 or CDH23 were associated with eGFRchange (n=1476) or time to a composite renal event that consisted of incident end stage renal disease or halving of eGFR (n=1585, with a total of n=178 events; results not shown).

Discussion

Key findings

Our key findings are fourfold. First, we estimate the heritability of eGFR decline as being 38% in the general population of European descent, providing a rationale to search for genetic variants associated with kidney function decline. Second, we extend evidence of a known locus (UMOD) previously associated with incident CKD and ESRD [21,28] by showing genome-wide significant association with kidney function change. Third, we have identified two novel genetic loci (CDH23 and GALNTL5/GALNT11) with suggestive association with kidney function decline phenotypes. Finally, we show that knock-down of the two novel loci in zebrafish renders the nephron susceptible to a nephrotoxic insult.

Our findings in the context of the literature

We extend the current literature by performing the first large-scale GWAS of renal function decline traits in the general population. Previous studies analyzing progression of renal disease in African Americans [29−32], individuals of European descent [21], healthy nurses [33], and patients with diabetes [34,35], hypertension [31], IgA nephropathy [36,37] and ESRD [21] focused only on candidate genes.

The SNP in UMOD has previously been identified in a GWAS of eGFR measured at one time point [14], and was significantly associated with incident CKD and ESRD in a candidate gene study [21] and with salt-sensitive hypertension and kidney damage in rodents and humans [38]. Our data extend this knowledge base by providing strong evidence that genetic variation at the UMOD locus affects different definitions of kidney function decline.

For Rapid Decline, the associated region on chromosome 7 contains the genes GALNTL5, GALNT11, MLL3, and CCT8L1, with our zebrafish data suggesting GALNTL5 and GALNT11 as the genes of interest. GALNTL5 encodes the putative polypeptide N-acetylgalactosaminyltransferase-like protein 5, which by similarity has a presumed role in O-linked oligosaccharide biosynthesis. Polypeptide N-acetylgalactosaminyltransferase 11, encoded by GALNT11, is a glycosyl transferase that catalyzes the initial reaction in O-linked oligosaccharide biosynthesis. Studies in Xenopus support a role of the gene product in left-right patterning by modulating Notch1 signaling and thus establishing the crucial balance between motile and immotile cilia, and it is also expressed in the developing kidney of zebrafish [39,40]. Our data suggest that galnt11 is not essential for kidney development, but protects against susceptibility from nephrotoxins.

The region of chromosome 7 also contains a locus (rs7805747 in PRKAG2) that was previously identified in a GWAS meta-analysis of cross-sectional eGFR [15]. However, this SNP is independent of rs1019173 (r2=0.002, D’=0.061 in the 1000 Genomes Pilot Version 1, hg18); therefore, the novel locus identified in the present study is unlikely tagging the PRKAG2 locus. Moreover, conditional analysis using genotypes from both SNPs from individual level data from the ARIC study showed that the association between rs1019173 and Rapid Decline is unchanged when controlling for rs7805747 (data not shown).

The other locus identified from this study is an intronic SNP in CDH23 that is nominally associated with eGFR change in those with CKD at baseline. CDH23 encodes cadherin 23, a glycoprotein of the cadherin family. Cadherin 23 and protocadherin 15, encoded by PCDH15, form the tip-links spanning the stereocilia of the inner ear’s hair cells. These tip-links are key contributors to the mechanosensory transduction in hair cells required for hearing [41]. Rare mutations of CDH23 cause progressive, nonsyndromic deafness (DFNB12, MIM # 601386) [42−44] or Usher Syndrome 1D, characterized by profound deafness, vestibular dysfunction and retinitis pigmentosa (MIM # 601067). The transmembrane protein cadherin 23 is expressed in many tissues, including the kidney [44,45], where it is found predominantly in the tubulointerstitium [46]. While a kidney phenotype has not been reported for patients with DFNB12 or Usher syndrome, our zebrafish data provide evidence that cadherin 23 plays a role in protecting from susceptibility to nephrotoxins, while not being essential for nephrogenesis.

Implications

Our GWAS findings point towards two novel gene loci, CDH23 and GALNTL5/GALNT11, and one previously identified locus (UMOD) as being associated with kidney function decline. The zebrafish experiments support a role of the two newly identified loci in increasing renal susceptibility to nephrotoxic insults and may indicate that a perturbation model could serve as a model of longitudinal kidney function decline. In previous work, we have shown that knockdown of two genes identified by GWAS of cross-sectional eGFR, mpped2 and casp9, resulted in abnormal kidney development, with susceptibility to gentamicin only in casp9 knockdown [16]. Taken together, our current and previous data highlight the differential role of genes in affecting kidney development, function and susceptibility to damage.

Strengths and Limitations

Strengths of this study include the large sample size of renal function decline traits, follow-up in independent samples, analysis of several definitions of kidney function decline and validation in zebrafish. Some limitations warrant mention. Even though we addressed inter-assay differences of serum creatinine measurement by calibrating creatinine to representative NHANES standards, several other factors causing imprecision in defining kidney function decline phenotypes may have reduced our statistical power to identify genome-wide significant associations: 1) despite our use of different renal function decline definitions all featured in current guideline statements [3], there is no standard definition of renal function decline, 2) kidney function trajectories are less well-defined with two vs. several serum creatinine measurements given that renal function change may not be linear over time [3] and there may be day-to-day alterations in GFR, 3) GFR estimation equations are known to be imprecise especially at a GFR>60 ml/min/1.73m2, 4) we observed heterogeneity in design between studies including a wide range of length of follow-up. We cannot rule out that low statistical power also accounts for the negative finding in the CRIC study. Further, our findings, obtained mainly in general population cohorts, provide novel insights into mechanisms of kidney function decline, but may not be generalizable to cohorts enriched for CKD. This limitation deserves particular attention due to the unexpected observation that in most cohorts, the subgroup with baseline CKD (defined as eGFR<60 ml/min/1.73m2) showed a mean increase in eGFR over time irrespective of length of follow-up interval. This may indicate that in the CKD subgroup of these cohorts, a baseline eGFR<60 ml/min/1.73m2 may not represent progressive CKD with active disease but rather stable disease or imprecise GFR estimation. This highlights that more work with expanded datasets and functional models are necessary to further elucidate the genetics of CKD initiation and progression in population-based studies. Finally, the role of genes contributing to aging and chronic disease in humans may not be entirely modeled by transient morpholino knockdown and observation of a developmental phenotype: while zebrafish allows high throughput modeling of the effects of gene knockdown in a vertebrate organism, the developmental role of specific genes may well be different from homeostatic organ maintenance in the adult. Specifically, umod may not play a relevant role in zebrafish renal development or toxin susceptibility.

Conclusion

In a large GWAS of kidney function decline phenotypes in individuals of European descent, we showed that a SNP in UMOD is associated with kidney function decline phenotypes, and that there is suggestive statistical evidence for two novel loci (GALNTL5/GALNT11 and CDH23). Zebrafish experiments at the two novel loci suggest roles in the deterioration of kidney function after acute injury. Given the complexity of the kidney function decline phenotype, further interrogation of these regions is warranted.

Materials and Methods

Ethics Statement

In all studies, all participants gave informed consent. All studies were approved by their responsible Research Ethics Committees.

Phenotype definition

Serum creatinine was measured at a minimum of 2 time points spaced several years apart (2.0 – 22.2 years, median 5.6 years). In almost all studies, there were only two serum creatinine measurements in total. To be consistent across studies, we used each individual’s two creatinine measurements with the longest follow-up in between for phenotype creation in all cohorts (see below). Baseline and follow-up serum creatinine were calibrated to the US nationally representative National Health and Nutrition Examination Study (NHANES) data in all discovery and replication studies to account for between-laboratory variation [47]. In order to be consistent with our prior work, GFR based on serum creatinine (eGFRcrea) was estimated using the four-variable MDRD Study Equation. eGFRcrea values <15 ml/min/1.73m2 were set to 15, and those >200 were set to 200 ml/min/1.73m2.

Several phenotypes were used to model different mechanisms involved in change of renal function over time, using each individual’s two serum creatinine measurements with the longest follow-up. The continuous phenotype eGFRchange, modeling annual change in kidney function, was calculated by subtracting the eGFR at follow-up from the eGFR at baseline, and then dividing by the number of years of follow-up for each participant. Thus, a positive value of eGFRchange corresponds to a decline in kidney function over time, whereas a negative value of eGFRchange corresponds to an increase in kidney function over time. Three dichotomous phenotypes were calculated to model kidney function decline phenotypes with different clinical implications [5,24]: For Rapid Decline, cases were defined as individuals with a rapid decline in kidney function >= 3ml/min/1.73 m2 per year, and controls as those with a kidney function decline < 3ml/min/1.73 m2 per year [6]. For incident CKD (CKDi), cases were defined as participants with eGFR at baseline >= 60 ml/min/1.73m2 declining to an eGFR at follow-up < 60ml/min/1.73 m2; a more stringent definition of incident CKD (CKDi25) is restricted to incident CKD cases with a decline of eGFR >= 25% at follow-up. For both CKDi and CKDi25, controls were defined as those with an eGFR >= 60ml/min/1.73 m2 at baseline and follow-up.

Heritability of eGFR in the Framingham Heart Study

Heritability of eGFRchange was calculated with family data of the Framingham Heart Study using the variance components analysis implemented in SOLAR [48]. eGFRchange was calculated by taking follow-up eGFR (obtained between 2005−2008) and subtracting baseline eGFR (obtained in 1995−1998), divided by the number of years of follow-up. Residuals were created after adjusting for age, sex, baseline eGFR, and principal components as necessary. With residuals as response variable, a variance components model with an additive genetic and a random environmental variance components was fitted, where the correlation among relatives attributable to the genetic component is assumed proportional to the kinship coefficient matrix. Heritability is calculated as the ratio of the estimated genetic variance to the total phenotypic variance.

Definition of strata

Kidney function decline is known to differ depending on level of baseline eGFR. Thus, eGFRchange was analyzed (A) in the overall sample [eGFRchange overall], (B) in those with eGFR >= 60 ml/min/1.73m2 at baseline [eGFRchange noCKD], and (C) in those with eGFR < 60 ml/min/1.73m2 at baseline [eGFRchange withCKD]. Rapid Decline was analyzed in the overall sample [Rapid Decline overall] and in those with eGFR >= 60 ml/min/1.73m2 at baseline [Rapid Decline noCKD]. CKDi and CKDi25 were analyzed in the overall sample only.

Stage 1 genome-wide association analyses

All participating studies used a uniform analysis plan and each trait was created using standard programming commands that were provided to collaborating studies. The continuous trait (eGFRchange) was analyzed by linear regression, the dichotomous traits by logistic regression (Rapid Decline, CKDi, CKDi25). Models included the allelic dosage at each marker from imputed study data consisting of 2.5 million HapMap-II SNPs [49] on average, based on imputations with different programs and reference panels. Details of genotyping and imputation in each study are shown in Supplementary Table 2. We used the additive genetic model, adjusted for age and sex, baseline eGFR and, where applicable, for study site and principal components.

Stage 1 meta-analysis

For our stage 1 analysis, we used aggregated statistics of 16 population-based GWA studies of individuals of European ancestry for each of the longitudinal traits: eGFRchange overall, eGFRchange noCKD, eGFRchange with CKD, Rapid Decline overall, Rapid Decline noCKD, CKDi and CKDi25. All 16 stage 1 studies contributed data to every trait, except for the AMISH study, which provided data to eGFRchange overall and eGFRchange no CKD only due to low number of CKD cases at baseline and follow-up.

All input files underwent quality control using the GWAtoolbox package in R (www.eurac.edu/GWAtoolbox.html) [50], before including them into meta-analysis. Study data was meta-analyzed assuming fixed effects and using inverse-variance weighting. Thus the pooled effect βpooled is estimated as , where β and SE are the effect and standard error of the SNP on the outcome in the ith study. The meta-analyses were performed by METAL. We performed genomic control correction if the inflation factor λ in the study files was greater than 1 (1st GC correction) or if it was greater than 1 in the meta-analysis result (2nd GC correction) [51].

Next, we created a list of independent SNPs (pairwise r2<0.2, HAPMAP II release 22) that had a genomic control corrected p-value <10−6 and minor allele frequency > 5% in stage 1 meta-analysis and were present in at least 85% of the contributing studies.

The stage 2 meta-analysis of SNPs identified in stage 1 was performed on the same phenotypes and using the same analysis plan as the stage 1 analysis, and was based on in silico genetic data or on de novo genotyped variants. Details on each stage 2 study’s genotyping and imputation platforms are shown in Supplementary Table 2. In addition, we also performed a combined inverse-variance weighted fixed-effects stage 1 and stage 2 meta-analysis using individual study files as input. Studies with less than 50 cases for a dichotomous trait or with an overall sample size of less than 50 for a continuous trait were excluded from the meta-analyses of the corresponding trait. SNPs with a stage 2 meta-analysis one-sided p-value <0.05 and effect direction consistency with the stage 1 meta-analysis effect direction were defined as showing nominally significant evidence of replication. The I2 statistic was computed to assess heterogeneity between studies.

Zebrafish functional experiments

Zebrafish were maintained according to established IACUC protocols. Zebrafish were injected at the 1-cell stage with 2 nl of 400 uM morpholinos (MO; GeneTools, Philomath, OR) designed to block the ATG start site or an exon-intron splice site of the target gene (Supplementary Table 6). Embryos were fixed in 4% PFA at the appropriate stages for in situ hybridization using well established protocols (http://zfin.org/ZFIN/Methods/ThisseProtocol.html). Renal gene expression was visualized using established markers for pax2a (global kidney) and nephrin (podocytes) [52,53]. The number of embryos displaying abnormal renal gene expression was compared to uninjected control embryos, and statistical significance was determined by Fisher’s exact test. For the gentamicin nephrotoxin experiment, embryos were injected with MO at the 1-cell stage and then injected with 5 nl of 10 mg/ml gentamicin prepared from one stock solution in the cardiac sinus venosus at 48 hpf after being anesthetized in a 1:20 dilution of 4 mg/ml Tricaine in embryo water. Live embryo development and edema prevalence was documented over the next three days.