First genome-wide association study in an Australian aboriginal population provides insights into genetic risk factors for body mass index and type 2 diabetes

Abstract

A body mass index (BMI) >22kg/m2 is a risk factor for type 2 diabetes (T2D) in Aboriginal Australians. To identify loci associated with BMI and T2D we undertook a genome-wide association study using 1,075,436 quality-controlled single nucleotide polymorphisms (SNPs) genotyped (Illumina 2.5M Duo Beadchip) in 402 individuals in extended pedigrees from a Western Australian Aboriginal community. Imputation using the thousand genomes (1000G) reference panel extended the analysis to 6,724,284 post quality-control autosomal SNPs. No associations achieved genome-wide significance, commonly accepted as P<5x10-8. Nevertheless, genes/pathways in common with other ethnicities were identified despite the arrival of Aboriginal people in Australia >45,000 years ago. The top hit (rs10868204 Pgenotyped = 1.50x10-6; rs11140653 Pimputed_1000G = 2.90x10-7) for BMI lies 5' of NTRK2, the type 2 neurotrophic tyrosine kinase receptor for brain-derived neurotrophic factor (BDNF) that regulates energy balance downstream of melanocortin-4 receptor (MC4R). PIK3C2G (rs12816270 Pgenotyped = 8.06x10-6; rs10841048 Pimputed_1000G = 6.28x10-7) was associated with BMI, but not with T2D as reported elsewhere. BMI also associated with CNTNAP2 (rs6960319 Pgenotyped = 4.65x10-5; rs13225016 Pimputed_1000G = 6.57x10-5), previously identified as the strongest gene-by-environment interaction for BMI in African-Americans. The top hit (rs11240074 Pgenotyped = 5.59x10-6, Pimputed_1000G = 5.73x10-6) for T2D lies 5' of BCL9 that, along with TCF7L2, promotes beta-catenin's transcriptional activity in the WNT signaling pathway. Additional hits occurred in genes affecting pancreatic (KCNJ6, KCNA1) and/or GABA (GABRR1, KCNA1) functions. Notable associations observed for genes previously identified at genome-wide significance in other populations included MC4R (Pgenotyped = 4.49x10-4) for BMI and IGF2BP2 Pimputed_1000G = 2.55x10-6) for T2D. Our results may provide novel functional leads in understanding disease pathogenesis in this Australian Aboriginal population.

Fig 1. BMI by age plotted using the R package SITAR.

(A) plot of all records for self—reporting Aboriginals (N = 1020) in the Aboriginal Health Service’s Communicare database; and (B) plot of all records for the 391 genotyped individuals contributing to association analyses for BMI and T2D. Separate lines trace multiple measurements over time per individual; females (pink) and males (blue). The heavy lines in (A) show the polynomial quintic (power of 5) curves for females (pink heavy line) and males (blue heavy line) that best fit the data. Fitting separate (by gender) curves to the data provided a significantly better fit (P = 10^-9) than fitting a single curve. The extreme outlier was not used in fitting the female polynomial curve. In (B), individuals with T2D are shown in heavy lines.

Fig 2. Manhattan plots of genome-wide association results for BMI undertaken using FASTA in GenABEL.

(A) results for genotyped data; and (B) results for 1000G imputed data. SNPs in red show the region of the top association in this discovery GWAS.

Fig 3. Regional association plots (LocusZoom [35]) of the signal for BMI association in the region SLC28A3 to NTRK2 on chromosome 9.

(A) is the plot for genotyped data; and (B) is the plot for 1000G imputed data. The −log₁₀ P-values are shown on the upper part of each plot. SNPs are colored (see key) based on their r ² with the labeled hit SNP (purple), calculated in the 146 unrelated genotyped individuals. The bottom section of each plot shows the genes marked as horizontal lines. The second Y axis is for recombination rate, as shown in blue on the plot.

Fig 4. Imputation accuracy for 402 genotyped individuals imputed against the 1000G reference panel.

Imputation accuracy is measured as average r ² across all autosomes for SNPs of different MAFs (see key).

Fig 5. SYNPLOT [38] for a section of the intergenic region (NCBI Build 37: 87,140,000bp to 87,190,000bp) ~94kb upstream of NTRK2 on chromosome 9q21.33.

Since regulatory elements are usually found with conserved regions of the genome, we interrogated the region of our top hits using this in silico analysis to look for conserved regions across multiple species. The multiple alignments of (top to bottom) human, cow and mouse sequences were performed across the complete intergenic region SLC28A3 to NTRK2 using LAGAN. The segment of the alignments shown here is the region that contains the top association hits (P<5×10^-6). The central plotted curves show the degree of conservation of sequence across all three species, on a vertical scale 0–1 (= 100%), such that the peaks represent CNS. CNS are defined here as regions with a nucleotide sequence conservation level of ≥0.7, i.e. ≥ to the least conserved exon sequence in the two genes (not shown on this plot) SLC28A3 and NTRK2 flanking the intergenic region of interest. Blue boxes indicate repetitive sequence in all 3 species. The human sequence is also annotated with positions of: (i) LINE-1 elements (yellow); (ii) the top genotyped SNPs (red vertical bars) including the top SNP rs1086204, and the SNP of interest rs1866439; (iii) the top imputed SNPs (blue vertical bars) including the top SNP rs1140653; (iv) the CpG island-like element identified using CpG Island Seacher (green); and (v) positions of peaks of mono-methylation of lysine 4 of histone H3 (H3K4Me1) as measured in NHEK or NHLF cell lines (mauve) or in H1-hESC human embryonic stem cells (orange) by the ENCODE project, as identified using the UCSC browser (see S9 Fig.).