- Research article
- Open Access
Mitochondrial GWAS and association of nuclear – mitochondrial epistasis with BMI in T1DM patients
BMC Medical Genomics volume 13, Article number: 97 (2020)
BMI is a strong indicator of complications from type I diabetes, especially under intensive treatment.
We have genotyped 435 type 1 diabetics using Illumina Infinium Omni Express Exome-8 v1.4 arrays and performed mitoGWAS on BMI. We identified additive interactions between mitochondrial and nuclear variants in genes associated with mitochondrial functioning MitoCarta2.0 and confirmed and refined the results on external cohorts: the Framingham Heart Study (FHS) and GTEx data. Linear mixed model analysis was performed using the GENESIS package in R/Bioconductor.
We find a borderline significant association between the mitochondrial variant rs28357980, localized to MT-ND2, and BMI (β = − 0.69, p = 0.056). This BMI association was confirmed on 1889 patients from FHS cohort (β = − 0.312, p = 0.047). Next, we searched for additive interactions between mitochondrial and nuclear variants. MT-ND2 variants interacted with variants in the genes SIRT3, ATP5B, CYCS, TFB2M and POLRMT. TFB2M is a mitochondrial transcription factor and together with TFAM creates a transcription promoter complex for the mitochondrial polymerase POLRMT. We have found an interaction between rs3021088 in MT-ND2 and rs6701836 in TFB2M leading to BMI decrease (inter_pval = 0.0241), while interaction of rs3021088 in MT-ND2 and rs41542013 in POLRMT led to BMI increase (inter_pval = 0.0004). The influence of these interactions on BMI was confirmed in external cohorts.
Here, we have shown that variants in the mitochondrial genome as well as additive interactions between mitochondrial and nuclear SNPs influence BMI in T1DM and general cohorts.
Mitochondria are organelles whose main role is energy production. They are the only organelles that contain their own genome. The mitochondrial genome is a double stranded 16.5 kb long molecule which resembles that of an alpha-proteobacterium . The two strands differ in nucleotide (G + T) composition – the guanosine-rich strand is named the heavy (H) strand; the other, cystosine-rich strand is called the L-strand (light strand). Like a free-living eubacterium, mitochondrial genome contains no introns and minimal intergenic regions, although it contains an approximately 1000 nucleotide non-coding control region (displacement loop or D-loop) where the origin of replication of the H-strand as well as promoters of transcription, both for H- and L-strands are localized. mtDNA codes for 37 genes. Among them 13 code for polypeptides, while the remainder - 2 rRNAs (12S and 16S) and 22 tRNAs – are necessary for mitochondrial protein synthesis. All 13 mRNAs code for subunits of oxidative phosphorylation (OXPHOS) complexes. The rest of the peptides needed to build the electron transport chain (ETC), as well as to maintain mitochondrial functioning, are nuclear-encoded . It has been shown that mitochondrial function correlates with cells’ metabolic state and can influence obesity .
Obesity is a severe epidemic world-wide . Current trends suggest that by the year 2030 more than 50% of Americans will be obese , while morbid obesity will affect even 10% of the UK population . Type 1 diabetics have seen an even greater increase in obesity incidence than is observed in the general population, , with obesity in type 1 diabetics becoming two times more prevalent over the last 30 years [8, 9]. Obesity in type 1 diabetics has been associated with severe complications, especially in patients undergoing intensive therapy . Such a burst in obesity prevalence can, in part, be attributed to lifestyle and to higher doses of insulin , however obesity is also highly heritable. Twin studies suggest that 40–70% of variability of BMI, the most popular measure to assess obesity, can be attributed to genetic variation . Even though large scale genome wide association studies (GWAS) including hundreds of thousands of individuals and millions of autosomal single nucleotide variants have been performed, they led to discovery of only around 100 genetic variants associated with BMI [13,14,15]. Thus, a substantial part of genetic variation that influences BMI remains to be discovered. As with most complex traits, most of the associated variants are non-coding (making explanation of their role in obesity even more difficult) .
Ultimately, obesity results from an imbalance between energy intake and its expenditure. Since interaction and communication between nuclear and mitochondrial genomes is indispensable for normal cell function [17, 18], it seems reasonable to look for interactions between SNPs in the nuclear genome and SNPs in the mitochondrial genome associated with obesity [19,20,21]. Here, we have performed mitochondrial GWAS as well as a study of genetic interactions between mitochondrial and nuclear variants which are localized to genes known to have an influence on mitochondrial functioning, associating variants in both genomes with BMI.
Patients were recruited either in Department of Metabolic Diseases University Hospital in Krakow or in Division of Reproduction Department of Obstetrics, Gynecology and Gynecological Oncology, Poznan University of Medical Sciences. All patients enrolled to the study were young women with type 1 diabetes (T1D) and on insulin treatment, who were pregnant or were trying to conceive. Whole blood samples were drawn and stored at − 80 °C. This study was approved by the Bioethical Committees of the Jagiellonian University and Poznan University of Medical Sciences and performed according to the Helsinki Declaration. Written informed consent was collected from all patients.
DNA was extracted from whole blood with the use of automated nucleic acid extraction system Maxwell (Promega). Five hundred twenty-seven samples were genotyped on Illumina Infinium Omni Express Exome-8 v1.4 arrays according to manufacturer’s instructions.
Quality control (QC)
Genotypes were called by GenomeStudio software Genotyping module (version 2.0, Illumina Inc.) according to manufacturers’ instructions (Technical Note: Genotyping, Infinium® Genotyping Data Analysis, www.illumina.com). Briefly, we removed all samples with 10%GenCall score < 0.4, Call Rate < 0.95 and discordant sex information. Then, all SNPs on chrX, chrM, chr0 were removed and hard cut off metrics was applied. Next, we manually curated 201 variants on chrM also in Genome Studio, recalculated statistics and excluded samples with Call Rate < 0.99 (180 remained in the analysis). Remaining nuclear genotypes for 940,911 SNPs and 513 samples were exported as custom report using plink (version 1.09) .
Second step of quality check was performed using custom rscript in Rstudio (version 1.1.383)  with R (version 3.4.2)  according to plink manufacturers’ instructions and KING . Briefly, we removed duplicates with lower Call Rate and samples with cryptic kinship, equivalent to third-degree relatives and higher (kinship coefficient > 0.0442). In addition, we also checked for heterozygosity rate of autosomal SNPs and removed variants with deviation from the Hardy-Weinberg equilibrium (HWE) with P < 1e10–5. The last part of QC was based on EIGENSOFT [26, 27] in order to determine population structure in our samples based on HapMap3 dataset . Samples that failed to qualify as European population (CEU) were removed. Finally, 476 individuals with 940,374 nuclear variants passed QC and were included in the analysis, however, due to lack of metadata only 435 patients were included in the final analysis. Apart from nuclear variants, 180 mitochondrial variants were included in the analysis.
The QC-filtered genotype data were checked against reference panel of the Haplotype Reference Consortium (HRC r1.1 2016, ftp://ngs.sanger.ac.uk/production/hrc/HRC.r1-1/HRC.r1-1.GRCh37.wgs.mac5.sites.tab.gz) with  “HRC/1KG Imputation Preparation and Checking Tool” (v.4.2.9, https://www.well.ox.ac.uk/~wrayner/tools) to exclude strand coding issues during the imputation step. Imputation of QC-filtered genotypes was performed on Michigan Imputation Server (using Minimac3)  with the HRC r1.1 reference panel. Phasing was performed with ShapeIT (version 2.r790) . The Minimac3 output dosage files were then converted to hard calls in PLINK. For further analysis imputed MitoExome variants were chosen based on MitoCarta2.0 gene localization . The phenotypic data were available for 435 patients which were included into the final analysis.
Framingham heart study cohort
The genotype and clinical phenotype files were downloaded from dbGaP (project #5358; accession number phs000007). Only the genotypes corresponding to the SHARe substudy were used (obtained via OMNI5M genotyping array, phg000256). In the GLMs, the following covariates were used: age (phv00177930), gender (phv00177929) and blood glucose level (phv00007558). Our access was approved by the Ohio State University IRB (Protocol #2013H0096) (see Availability of data and materials).
Tissue specific RNA sequencing data was acquired from the Genotype and Tissue Expression Project (GTEx -v7) via the website gtexportal.org. The phenotye data was accessed via dbGaP Project #5358 (dbGaP accession number phs000424).
RNA sequencing had been generated using poly-adenylated priming with reads aligned to HG19, Gencodev12. For further details see Lonsdale et al.  and the GTEx website (http://www.gtexportal.org/home/documentationPage). We considered subcutaneous adipose and thyroid gland. Testing for eQTLs was performed via the GTEx online search tool . The GLM approach was utilized to test the association between BMI (phv00169070) (adjusted for age (phv00169063) and sex (phv00169062)) and expression of selected mitochondrial and/or nuclear mRNAs (see Availability of data and materials).
Generalized linear model (GLM) analysis
The linear mixed model analysis was performed using the GENESIS package in R/Bioconductor . In brief, first the genetic relatedness matrix was estimated via the PC-AiR and PC-Relate methods, subsequently a linear model was built and Wald’s test was used to assign significance. The interactions between variants was modeled via the option ‘ivars’ in function ‘assocTestMM’. All p-values reported are non-adjusted for multiple testing.
Filtering variants likely to be involved in epistatic interactions
For the analysis of epistatic interactions between the mitochondrial and nuclear (whole genome) variants a filtering step was performed to increase the likelihood of detection and limit the number of hypotheses tested. The details of our approach are presented in the Supplementary Methods. Our concept is closely related to the notion of variance QTLs [36, 37] i.e. variants which are associated with the variability of a trait. In brief, we compare the relative dispersion (conditional on the given nuclear SNP) of the residuals in two GLMs: (1) BMI ~ covariates+mito_variant, and (2) BMI ~ covariates+mito_variant+nuc_variant. The basic assumption is that if the residuals of the second GLM have higher relative (conditional on the nuclear SNP) dispersion that the residuals from the first model, then the nuclear SNP is more likely to interact with a variable (not necessarily the mitochondrial genotype) in the context of BMI. We test this assumption (together with our filtering approach) by means of simulation and random subsampling (from real-life data) and will present complete results in a forthcoming paper (results, data and code available upon request). Below we show that this approach allows to select variants likely to interact epistaticaly with the MT-ND2 rs28357980 variant (Table 5). At the same time, a more extensive analysis based on the FHS genotype data proves that our method allows to select subsets of variants which are enriched in epistatic interactions in the context of BMI.
Single SNP analysis – mitoGWAS
The cohort consisted of 435 young women with mean age of 28.5 years. Most of them were of normal weight and have never been pregnant before being included to the study. Mean disease duration was 12 years, while median insulin daily dose equaled 40 IU. Anthropometric data are presented in Table 1.
Forty-two mitochondrial variants were included in the analysis after filtering according to alternative allele frequency (Table S1a). Borderline significant nominal association with lower BMI was shown for 4 mitochondrial variants. Three of them were localized to non-coding part of mitochondrial genome – one to mitochondrial ribosomal gene MT-RNR2 and two to MT-tRNAs (Arg, Thr). The last one – rs28357980 was localized to MT-ND2 gene (MitoA4917G), which is part of the complex I of electron transport chain. The variant leads to amino acid change from Asparagine to Aspartic acid in position 150 of MT-ND2 protein. Results of mitoGWAS on BMI in T1DM cohort are presented in Table 2.
An analogous analysis of 87 mitochondrial variants (Table S1b) on 1889 subjects gathered in Framingham Heart Study (validation cohort) was performed. For the analysis, only the first measurement of BMI was used. The cohort consisted of 1037 women and 852 men. Their mean age was 34.6, half of them were normoglycemic. Median BMI for this cohort was 24.28 (Table 3).
Nominal associations with BMI were found for 7 variants. All nominally associated variants led to BMI decrease. Four of them were coding, the rest was localized to rRNA and tRNAs.
Among coding variants, we have found the same variant rs28357980 localized to MT-ND2 gene for which an association with BMI in T1DM cohort was found. Results of mitoGWAS on BMI in FHS cohort are presented in Table 4.
MT-ND2 gene - nuclear genes interactions
Each of the analysis we performed (mitoGWAS in both cohorts) directed us towards MT-ND2 gene. Since there must be a tight communication between the mitochondrial and the nuclear genome in order to achieve effective energy management [38,39,40,41], we looked for the interactions (i.e. epistasis between genetic variants in GLMs) between the two genomes. The results are presented in Table 5.
We also performed such an analysis on FHS cohort and found 53 nuclear variants which interact with rs28357980 (MT-ND2), several of which served as eQTLs. When these genes were subjected to GoTerm analysis we found an enrichment of terms associated with mitochondria (Mitochondrion padj = 9.56e-05, Mitochondrial part padj = 2.62e-04, Mitochondrial matrix padj = 2.52e-03).
MT-ND2 gene interactions within MitoCarta
The above results (i.e. the enrichment of genes involved in the functioning of mitochondria among the ones which have variants involved in epistatic interaction with rs28357980) guided us to further explore the potential for epistasis between the MT-ND2 variant and nuclear SNPs located in the MitoCarta genes, which are known to be associated with mitochondrial pathways (Table S2a, b). Most significant interactions for MT-ND2 were associated with TCA cycle and respiratory electron transport, metabolism of amino acids and mitochondrial biogenesis (Table S3). The last category can potentially have the most crucial influence on mitochondrial functioning and can have profound significance for energy metabolism. Reactome lists 17 genes (POLRMT, TFB1M, TFAM, MTERF, PERM1, TFB2M, ATP5B, SSBP1, SIRT3, POLG2, PPARGC1A, CYCS, NRF1, ALAS1, PEO1, GABPA, ESRRA) which constitute this pathway, five of which are present in the results of our analysis. The list of MT-ND2 variants interactions within this pathway is presented in Table 6.
We looked closer into these interactions and found that variant rs3021088 of MT-ND2 gene interacted with variants in TFB2M (rs6701836) and POLMRT (rs41542013) genes, both of which are necessary for mitochondrial transcription. TFB2M is a mitochondrial transcription factor which, together with TFAM, creates transcription promoter complex and enables transcription by mitochondrial polymerase POLRMT (Figure S1).
When looking at the interaction between TFB2M nuclear variant – rs6701836 and MT-ND2 mitochondrial variant – rs3021088 we found that the combination of the two led to BMI decrease [(nuc_eff(single model) = − 0.0660, nuc_pval(single) = 0.857, joint_eff(inter) = 0.0874, nuc_pval(inter) = 0.0773, inter_eff = − 2.225, inter_pval = 0.0241]. Neither nuclear nor mitochondrial variant on its own had such an effect [p(nuc) – p = 0.8577, p(mito) – p = 0.116].
Interaction between POLRMT variant rs41542013 and MT-ND2 mitochondrial variant rs3021088 led to BMI increase [nuc_eff(single model) = 0.296, nuc_pval(single) = 0.546, nuc_eff(inter) = − 0.085, joint_pval(inter) = 0.0016, inter_eff = 4.015, inter_pval = 0.0004]. None of these variants on its own was not associated with BMI.
Thus, we have found an example of negative interaction between MT-ND2 and TFB2M and positive interaction between MT-ND2 and POLRMT variants.
Next, we assessed the functional significance of these interactions. Mitochondrial rs3021088 is a missense variant which leads to Alanine to Threonine substitution in position 331 of MT-ND2 protein.
eQTL data (GTEx) for nuclear rs6701836 showed that it influenced TFB2M expression in thyroid gland. POLRMT variant rs41542013 in subcutaneous adipose tissue leads to lower POLRMT expression (Fig. 1). In GTEx data, mRNA levels of MT-ND2 and TFB2M correlated with higher BMI [p(nuc) = 0.0269, eff_nuc = 2.465e-01, p(mito) = 0.0244, eff_mito = 2.243e-04] in liver tissue, however, their interaction led to decrease of BMI [p = 0.0308, inter_eff = − 1.009e-05]. In GTEx data, mRNA levels of MT-ND2 and POLRMT on its own correlated with lower BMI [p(nuc) = 0.0492, eff_nuc = − 2.386e-01, p(mito) = 0.0688, eff_mito = − 2.233e-04] in subcutaneous adipose tissue, however, their interaction led to increase of BMI [p = 0.0235, inter_eff = 1.023e-05]. Taken together, the interactions on the mRNA level are in line with what was discovered in the T1DM cohort.
We also used FHS as a validation cohort to confirm the epistatic interactions between SIRT3, ATP5B, CYCS, TFB2M or POLRMT and MT-ND2.
We were able to confirm the interactions and their influence on BMI for SIRT3 and ATP5B genes. For SIRT3 we have found an interaction between mitochondrial variant localized to 5460 position (MAF 10.1%) and two nuclear variants – rs11602954 (p = 0.03) and rs11606393 (p = 0.028). These variants are in LD with each other and with previously described variants in SIRT3 (in Table 6). For ATP5B gene we have found that mitochondrial variant localized to position 4769 interacted with seven variants (rs2255074, rs1107479, rs7973157, rs2950390, rs2958154, rs2290893, rs2035081) which were in LD with rs2950393.
We looked closer into interactions of MT-ND2 and TFB2M or POLRMT. For TFB2M gene, we have found two interactions. Variant rs10924779 interacted with MT-ND2 variant localized to 5460 position (MAF = 10.1%). Their interaction leads to BMI increase (p = 0.037), while each of these variants on their own did not have an influence on BMI. Variant rs4654291 interacts with the same mitochondrial variant and leads to BMI increase (p = 0.045), while each of these on its own does not affect BMI. For POLRMT gene we have identified an intronic variant rs10411491 (MAF = 3%) which interacted with mitochondrial variant positioned at 4647 (MAF = 4.8%). Their interaction led to BMI increase, although due to low intronic variant MAF it did not reach statistical significance.
Thus, we can conclude that interactions between TFB2M or POLRMT and MT-ND2 gene might affect BMI.
Whole genome MT-ND2 interactions
Next, we investigated whole genome MT-ND2 interactions. All variants which entered statistically significant interactions with MT-ND2 gene in FHS cohort and T1DM cohort were listed (Table S4a, b). These variants were compared to results gathered in GWAS catalog. Among 366 nuclear variants which interacted with mitochondrial variant rs28357980, 5 were found in GWAS Catalog and 2 of them met genome wide significance threshold. Out of 1157 variants which interacted with mitochondrial variant rs3021088, 35 were listed in GWAS Catalog and 13 were significant at the level of p < e10–09 (Table S4c). According to GWAS Catalog these variants were associated, among others, with blood metabolite levels, blood protein levels, LDL and total cholesterol, bone mineral density, age-related hearing impairment.
Functioning of mitochondria have a profound effect on whole organism. The variants within mitochondrial as well as nuclear genomes might have an influence on the levels of ATP produced within mitochondria. ATP depletion can induce endoplasmic reticulum (ER) stress and lead to reactive oxygen species (ROS) generation, which later on will result in mitochondrial DNA damage and create a vicious cycle of mitochondrial inefficacy. Here we have looked for mitochondrial variants and their interaction with nuclear variants which may be associated with BMI in the general population as well as in T1DM patients.
Our mitoGWAS analyses have shown that both protein coding and non-coding variants are associated with BMI. We have however looked closer into the former, as they might directly influence electron transport chain (ETC) and cells energy production. The analysis of T1DM cohort led to discovery of a variant in MT-ND2 gene. The predictive algorithms do not suggest that rs28357980 (G4917A) might be damaging (SIFT - 0.09, PolyPhen – 0.06), however the variant has been shown to be associated with multiple sclerosis  and colorectal cancer . The G4917A variant is one of the defining variants of haplogroup T, which in 2013 was shown to be a risk factor for morbid obesity . Moreover, the variant is non-synonymous and affects the protein by substitution of highly conserved Aspartate at amino acid 150 to an Asparagine. What is more interesting, our analysis of the validation cohort (FHS) confirmed the association of the same variant with BMI.
rs28357980 (G4917A) is localized to MT-ND2 gene which is part of complex I of ETC. It is the first site of oxidative phosphorylation. Complex I is built of 46 proteins, 7 of which are mitochondrially coded (ND1, ND2, ND3, ND4, ND4L, ND5, ND6) and form very hydrophobic subunits within mitochondrial membrane. Disruptive mutations in ND subunits are commonly found as somatic mutations in tumors, but are not found as germline mutations associated with human diseases, due to their lethality . Previous studies have already associated variants in these genes with obesity. A paper from 2014 has shown, apart from associations with two other positions, an association between three variants in complex I genes – MT-ND1, MT-ND2 and MT-ND4L . Moreover, the variant C5178A in MT-ND2 gene was shown to lead to lower incidence of autoimmune diabetes. The A allele was protective against both autoimmune and alloxan-induced free radical–mediated diabetes in mice, possibly by suppressing ROS production at the β-cell level [47,48,49]. Apart from metabolic disorders, mutations in genes of complex I of ETC were shown to be associated with childhood acute lymphoblastic leukemia or were shown to be a poor prognostic factor in oral cancer .
Since several GWAS data show that a substantial part of genetic variability of obesity is still unknown one might suspect that it is hidden in more complex associations, meaning genetic interactions. It is known that mitochondrial functioning is a result of anterograde and retrograde signaling between mitochondrial and nuclear genomes. Most of the studies performed by now, did not analyze point mitochondrial variants, but conplastic animals in which nucleus from one organism was fused with cytoplasts of the other. Such experiments have shown that introduction of exogenous mitochondria (in which mitochondrial genome differed from the original by one or few nucleotides) influenced organismal phenotype. For example, conplastic rats where shown to have impaired glucose tolerance, while mice were more resistant to experimental autoimmune encephalomyelitis , had disrupted activity of the components of TCA cycle  or altered mitochondrial and cellular adaptation during aging . Moreover, the mutation in MT-ND2 gene (C4738A) in mouse fibroblasts led to significantly higher mitochondrial complex I activity, enhanced ATP production, reduced ROS production with similar MT-ND2 protein expression levels . A lot of studies have shown that even a point mutation in mitochondrial genome which was introduced onto another nuclear background led to severe mitochondrial dysfunction. A point mutation in ATP8 gene (7778 G/T) in C57BL/6 N-mtFVB/N mice led to lower insulin secretion in isolated islets after glucose stimulation when compared with C57BL/6 N-mtAKR/J mice. It also had reduced mitochondrial function in brain, spleen and liver  as well as showed a 3-fold higher generation of mitochondrial ROS production compared to C57BL/6 N-mtAKR/J mice .
However, since mitochondrial genome mutates faster than nuclear (for example due to ROS proximity) the incompatibility between the two genomes might occur during organismal lifetime. What is more, these interactions can also be influenced by the environment.
Thus, we have looked into epistatic interactions of MT-ND2 variants. We therefore checked the potential for interactions of MT-ND2 variant located in position 4917 in T1DM cohort. Our whole genome analysis has shown that it possibly interacts with 12 variants, most of which might influence energy metabolism or processes that have been shown to be associated with obesity. One variant was located in HS1BP3 gene. HS1BP3 is known to be localized to mitochondria, be involved in autophagy by inhibiting phospholipase D and its overexpression leads to apoptosis [57,58,59]. We also found four variants, in LD with each other, located either within or next to NKAIN3 gene, which turned out to be eQTLs for GGH gene. Gamma glutamyl-hydrolase is an enzyme involved in folate metabolism, while lower folates levels were associated with reduced insulin sensitivity and obesity [60,61,62]. Moreover, we identified 5 intergenic variants located near SLC25A6P2 gene, which is a mitochondrial transporter family [63,64,65]. The analysis of rs28357980 (MT-ND2) interactions in FHS also pointed us toward analysis of nuclear mitochondrial genes.
The Reactome analysis done on MitoCarta genes has shown an enrichment of interactions that are associated with mitochondrial biogenesis. Our data show that POLRMT and TFB2M variants interact with variants in MT-ND2 and affect BMI. POLRMT and TFB2M are genes that act in concert to perform mitochondrial transcription and replication thus variation within their sequence, and simultaneously in the mitochondrial DNA sequence, can have an influence on mitochondrial and gene copy number as well as influence efficacy of the two processes [66,67,68]. Since it is very difficult to assess the significance of such interactions we also looked into all interactions of MT-ND2 variants and checked whether nuclear variants which interacted with mitochondrial MT-ND2 gene were listed as significant in any GWAS study performed until today, as we believed this would strengthen our findings. Our analysis has confirmed that some of the nuclear variants were significantly associated with traits which are part of obesity phenotype, e.g. cholesterol level, blood metabolites level or with diseases which are known to be influenced by mitochondrial deficiencies e.g. hearing loss.
In conclusion, here we find that rs28357980 localized to MT-ND2 gene of mitochondrial genome is associated with BMI both in T1DM and in general cohort. What is more, we show that genetic epistasis might influence obesity phenotype by interaction of variants in MT-ND2 gene with nuclear variants in genes responsible for mitochondrial replication and transcription.
Availability of data and materials
Sequence data has been deposited at the European Genome-phenome Archive (EGA), which is hosted by the EBI and the CRG, under accession number EGAS00001004408. Further information about EGA can be found on https://ega-archive.org “The European Genome-phenome Archive of human data consented for biomedical research” (http://www.nature.com/ng/journal/v47/n7/full/ng.3312.html).
The FHS datasets were accessed via dbGaP (project #5358; accession number phs000007). The phg000256 OMNI5M SHARe dataset (access pends appropriate approval) was used for analysis. To access phg000256 dataset please use the following link in the dbGAP advanced search: (https://www.ncbi.nlm.nih.gov/gap/advanced_search/?OBJ=genotype&TERM=phg000256&COND=%7B%22study_name_accession_combo%22:%5B%22NHLBI%20Framingham%20SNP%20Health%20Association%20Resource%20(SHARe)%20%20(phs000342.v19.p12)225D%7D), which allows to look for OBJ = genotype TERM = phg000256 in the study phs000342.v19.p12. The following variables were used for analysis: age (phv00177930) and sex (phv00177929) - both belonging to pht003099, as well as BMI (calculated from height: phv00007570 and weight: phv00007569 - belonging to pht000030).
Tissue specific RNA sequencing data was acquired from the Genotype and Tissue Expression Project (GTEx -v7) via the website gtexportal.org. The phenotye data was accessed via dbGaP Project #5358 (dbGaP accession number phs000424). The variables used for glm analysis using the GTEx project data were: BMI (phv00169070), age (phv00169063) and sex (phv00169062), all in pht002742.
The QC-filtered genotype data were checked against reference panel of the Haplotype Reference Consortium (HRC r1.1 2016, ftp://ngs.sanger.ac.uk/production/hrc/HRC.r1-1/HRC.r1-1.GRCh37.wgs.mac5.sites.tab.gz) with “HRC/1KG Imputation Preparation and Checking Tool” (v.4.2.9,https://www.well.ox.ac.uk/~wrayner/tools) to exclude strand coding issues during the imputation step. Imputation of QC-filtered genotypes was performed on Michigan Imputation Server (using Minimac3) with the HRC r1.1 reference panel.
All data used in this study are deposited in public repositories. The data are publicly available and access requires bioethics committee approval. The researchers who want to gain access to these data are encouraged to apply to the respective Institutional Review Boards. The sensitive data, including those which allow identification of individuals are available to anyone upon approval based solely on bioethics ground.
Body mass index
- T1DM :
Type 1 diabetes
Single nucleotide polymorphism
Electron transport chain
Genome wide association study
Minor allele frequency
Generalized linear model
Expression quantitative trait loci
Tricarboxylic acid cycle
Mishmar D, Zhidkov I. Evolution and disease converge in the mitochondrion. Biochimica et Biophysica Acta Bioenergetics. 2010;1797:1099–104.
Wallace DC, Chalkia D. Mitochondrial DNA genetics and the heteroplasmy conundrum in evolution and disease. Cold Spring Harbor Perspect Med. 2013;5(11):a021220. https://doi.org/10.1101/cshperspect.a021220.
Picard M, Turnbull DM. Linking the metabolic state and mitochondrial dna in chronic disease, health, and aging. Diabetes. 2013;62:672–8.
Arroyo K, Herron DM. The epidemiology of obesity. In: Bariatric Endoscopy; 2013.
Wang Y, Beydoun MA, Liang L, Caballero B, Kumanyika SK. Will all Americans become overweight or obese? Estimating the progression and cost of the US obesity epidemic. Obesity. 2008;16(10):2323–30.
Keaver L, Webber L. Future trends in morbid obesity in England, Scotland, and Wales: a modelling projection study. Lancet. 2016. https://doi.org/10.1016/S0140-6736(16)32299-1.
Liu LL, Lawrence JM, Davis C, Liese AD, Pettitt DJ, Pihoker C, et al. Prevalence of overweight and obesity in youth with diabetes in USA: the SEARCH for diabetes in youth study. Pediatr Diabetes. 2010;11(1):4–11.
Libman IM, Pietropaolo M, Arslanian SA, LaPorte RE, Becker DJ. Changing prevalence of overweight children and adolescents at onset of insulin-treated diabetes. Diabetes Care. 2003;26(10):2871–5.
Minges KE, Whittemore R, Grey M. Overweight and obesity in youth with type 1 diabetes. Annu Rev Nurs Res. 2013;31(1):47–69 Available from: http://openurl.ingenta.com/content/xref?genre=article&issn=0739-6686&volume=31&issue=1&spage=47.
Purnell JQ, Braffett BH, Zinman B, Gubitosi-Klug RA, Sivitz W, Bantle JP, et al. Impact of excessive weight gain on cardiovascular outcomes in type 1 diabetes: results from the diabetes control and complications trial/epidemiology of diabetes interventions and complications (DCCT/EDIC) study. Diabetes Care. 2017;40(12):1756–1762. https://doi.org/10.2337/dc16-2523.
The Diabetes Control And Complications Trial Research Group. Influence of intensive diabetes treatment on body weight and composition of adults with type 1 diabetes in the diabetes control and complications trial. Diabetes Care. 2001;24(10):1711–21.
Zaitlen N, Kraft P, Patterson N, Pasaniuc B, Bhatia G, Pollack S, et al. Using Extended Genealogy to Estimate Components of Heritability for 23 Quantitative and Dichotomous Traits. PLoS Genet. 2013;9(5):e1003520. https://doi.org/10.1371/journal.pgen.1003520. Epub 2013 May 30.
Locke A, Kahali B, Berndt S, Justice A, Pers T. Genetic studies of body mass index yield new insights for obesity biology. Nature. 2015;518(7538):197–206.
Speliotes EK, Willer CJ, Berndt SI, Monda KL, Thorleifsson G, Jackson AU, et al. Association analyses of 249,796 individuals reveal eighteen new loci associated with body mass index. Nat Genet. 2011;42(11):937–48.
Wen W, Cho YS, Zheng W, Dorajoo R, Kato N, Qi L, et al. Meta-analysis identifies common variants associated with body mass index in east Asians. Nat Genet. 2012;44(3):307–11.
Cotsapas C, Speliotes EK, Hatoum IJ, Greenawalt DM, Dobrin R, Lum PY, et al. Common body mass index-associated variants confer risk of extreme obesity. Hum Mol Genet. 2009;18(18):3502–7. https://doi.org/10.1093/hmg/ddp292. Epub 2009 Jun 24.
Nunnari J, Suomalainen A. Mitochondria: In sickness and in health. Cell. 2012;148:1145–59.
Wallace DC, Fan WW. The pathophysiology of mitochondrial disease as modeled in the mouse. Genes Dev. 2009;23:1714–36.
Latorre-Pellicer A, Moreno-Loshuertos R, Lechuga-Vieco AV, Sánchez-Cabo F, Torroja C, Acín-Pérez R, et al. Mitochondrial and nuclear DNA matching shapes metabolism and healthy ageing. Nature. 2016;535(7613):561–5.
Quirós PM, Mottis A, Auwerx J. Mitonuclear communication in homeostasis and stress. Nat Rev Mol Cell Biol. 2016;17:213–26.
Doynova MD, Berretta A, Jones MB, Jasoni CL, Vickers MH, O’Sullivan JM. Interactions between mitochondrial and nuclear DNA in mammalian cells are non-random. Mitochondrion. 2016;30:187–96 Available from: http://linkinghub.elsevier.com/retrieve/pii/S1567724916301325.
Purcell S, Neale B, Todd-Brown K, Thomas L, Ferreira MAR, Bender D, et al. PLINK: a tool set for whole-genome association and population-based linkage analyses. Am J Hum Genet. 2007;81(3):559–75 Available from: http://linkinghub.elsevier.com/retrieve/pii/S0002929707613524.
Citing RStudio – RStudio Support [Internet]. [cited 2018 Sep 10]. Available from: https://support.rstudio.com/hc/en-us/articles/206212048-Citing-RStudio.
R: The R Project for Statistical Computing [Internet]. [cited 2018 Sep 10]. Available from: https://www.r-project.org/.
Manichaikul A, Mychaleckyj JC, Rich SS, Daly K, Sale M, Chen WM. Robust relationship inference in genome-wide association studies. Bioinformatics. 2010;26(22):2867–73.
Patterson N, Price AL, Reich D. Population structure and eigenanalysis. PLoS Genet. 2006;2(12):2074–93.
Price AL, Patterson NJ, Plenge RM, Weinblatt ME, Shadick NA, Reich D. Principal components analysis corrects for stratification in genome-wide association studies. Nat Genet. 2006;38(8):904–9.
International T, Consortium H. The international HapMap project. Nature. 2003;426(6968):789–96.
McCarthy S, Das S, Kretzschmar W, Delaneau O, Wood AR, Teumer A, et al. A reference panel of 64,976 haplotypes for genotype imputation. Nat Genet. 2016;48(10):1279–83.
Das S, Forer L, Schönherr S, Sidore C, Locke AE, Kwong A, et al. Next-generation genotype imputation service and methods. Nat Genet. 2016;48(10):1284–7.
O’Connell J, Gurdasani D, Delaneau O, Pirastu N, Ulivi S, Cocca M, et al. A General Approach for Haplotype Phasing across the Full Spectrum of Relatedness. PLoS Genet. 2014;10(4):e1004234. https://doi.org/10.1371/journal.pgen.1004234. eCollection 2014 Apr.
Calvo SE, Clauser KR, Mootha VK. MitoCarta2.0: an updated inventory of mammalian mitochondrial proteins. Nucleic Acids Res. 2016;44(D1):D1251–7.
Lonsdale J, Thomas J, Salvatore M, Phillips R, Lo E, Shad S, et al. The genotype-tissue expression (GTEx) project. Nat Genet. 2013;45(6):580–5. https://doi.org/10.1038/ng.2653.
Aguet F, Ardlie KG, Cummings BB, Gelfand ET, Getz G, Hadley K, et al. Genetic effects on gene expression across human tissues. Nature. 2017.
Conomos MP, Gogarten SM, Brown L, Chen H, Rice K, Sofer T, et al. Type Package Title GENetic EStimation and Inference in Structured samples (GENESIS): Statistical methods for analyzing genetic data from samples with population structure and/or relatedness [Internet]. 2018 [cited 2018 Sep 13]. Available from: https://git.bioconductor.org/packages/GENESIS.
Rönnegård L, Valdar W. Recent developments in statistical methods for detecting genetic loci affecting phenotypic variability. BMC Genetics. 2012;13:63. https://doi.org/10.1186/1471-2156-13-63.
Wang G, Yang E, Brinkmeyer-Langford CL, Cai JJ. Additive, epistatic, and environmental effects through the lens of expression variability QTL in a twin cohort. Genetics. 2014;196(2):413–25. https://doi.org/10.1534/genetics.113.157503. Epub 2013 Dec 2.
Ryan MT, Hoogenraad NJ. Mitochondrial-Nuclear Communications. Annu Rev Biochem. 2007;76:701–22. https://doi.org/10.1146/annurev.biochem.76.052305.091720.
Eisenberg-Bord M, Schuldiner M. Ground control to major TOM: mitochondria–nucleus communication. FEBS J. 2017;284(2):196–210. https://doi.org/10.1111/febs.13778. Epub 2016 Jul 4.
Finley LWS, Haigis MC. The coordination of nuclear and mitochondrial communication during aging and calorie restriction. Age Res Rev. 2009;8(3):173–88. https://doi.org/10.1016/j.arr.2009.03.003. Epub 2009 Mar 27.
Moghadam AA, Ebrahimie E, Taghavi SM, Niazi A, Babgohari MZ, Deihimi T, et al. How the nucleus and mitochondria communicate in energy production during stress: nuclear MtATP6, an early-stress responsive gene, regulates the mitochondrial F1F0-ATP synthase complex. Mol Biotechnol. 2013;54(3):756–69. https://doi.org/10.1007/s12033-012-9624-6.
Andalib S, Talebi M, Sakhinia E, Farhoudi M, Sadeghi-Bazargani H, Gjedde A. Mitochondrial DNA T4216C and A4917G variations in multiple sclerosis. J Neurol Sci. 2015;356(1–2):55–60.
Li T, Yuan G, Zhang L, Ye L, Li S, Fan Y, et al. ApoG2 inhibits the antiapoptotic protein, Mcl-1, and induces mitochondria-dependent apoptosis in human colorectal cancer cells. Mol Med Rep. 2015;12(5):6976–84.
Nardelli C, Labruna G, Liguori R, Mazzaccara C, Ferrigno M, Capobianco V, et al. Haplogroup T is an obesity risk factor: mitochondrial DNA haplotyping in a morbid obese population from southern Italy. Biomed Res Int. 2013;2013:631082. https://doi.org/10.1155/2013/631082. Epub 2013 Jul 2.
Gonzalez-Halphen D, Ghelli A, Iommarini L, Carelli V, Esposti MD. Mitochondrial complex I and cell death: A semi-automatic shotgun model. Cell Death Dis. 2011;2(10):e222. https://doi.org/10.1038/cddis.2011.107.
Flaquer A, Baumbach C, Kriebel J, Meitinger T, Peters A, Waldenberger M, et al. Mitochondrial genetic variants identified to be associated with BMI in adults. PLoS One. 2014;9(8):e105116. https://doi.org/10.1371/journal.pone.0105116. eCollection 2014.
Gusdon AM, Votyakova TV, Reynolds IJ, Mathews CE. Nuclear and mitochondrial interaction involving mt-Nd2 leads to increased mitochondrial reactive oxygen species production. J Biol Chem. 2007;282(8):5171–9 Available from: http://www.ncbi.nlm.nih.gov/pubmed/17189252.
Gusdon AM, Votyakova TV, Mathews CE. Mt-Nd2a suppresses reactive oxygen species production by mitochondrial complexes I and III. J Biol Chem. 2008;283(16):10690–7.
Chen J, Gusdon AM, Piganelli J, Leiter EH, Mathews CE. Mt-Nd2a modifies resistance against autoimmune type 1 diabetes in NOD mice at the level of the pancreatic β-cell. Diabetes. 2011;60(1):355–9.
Uzawa K, Kasamatsu A, Baba T, Kimura Y, Nakashima D, Higo M, et al. Quantitative detection of circulating tumor-derived mitochondrial NADH subunit variants as a potential prognostic biomarker for oral cancer. Int J Oncol. 2015;47(3):1077–83.
Yu X, Gimsa U, Wester-Rosenlöf L, Kanitz E, Otten W, Kunz M, et al. Dissecting the effects of mtDNA variations on complex traits using mouse conplastic strains. Genome Res. 2009;19(1):159–65.
Scheffler K, Krohn M, Dunkelmann T, Stenzel J, Miroux B, Ibrahim S, et al. Mitochondrial DNA polymorphisms specifically modify cerebral β-amyloid proteostasis. Acta Neuropathol. 2012;124(2):199–208.
Niemann J, Johne C, Schröder S, Koch F, Ibrahim SM, Schultz J, et al. An mtDNA mutation accelerates liver aging by interfering with the ROS response and mitochondrial life cycle. Free Radic Biol Med. 2017;102:174–87.
Schauer M, Kottek T, Schönherr M, Bhattacharya A, Ibrahim SM, Hirose M, et al. A mutation in the NADH-dehydrogenase subunit 2 suppresses fibroblast aging. Oncotarget. 2015;6(11):8552–66 Available from: http://www.pubmedcentral.nih.gov/articlerender.fcgi?artid=4496166&tool=pmcentrez&rendertype=abstract.
Schröder T, Kucharczyk D, Bär F, Pagel R, Derer S, Jendrek ST, et al. Mitochondrial gene polymorphisms alter hepatic cellular energy metabolism and aggravate diet-induced non-alcoholic steatohepatitis. Mol Metab. 2016;5(4):283–95.
Weiss H, Wester-Rosenloef L, Koch C, Koch F, Baltrusch S, Tiedge M, et al. The mitochondrial Atp8 mutation induces mitochondrial ROS generation, secretory dysfunction, and β-cell mass adaptation in conplastic B6-mtFVBmice. Endocrinology. 2012;153(10):4666–76.
Yin Z, Klionsky DJ. HS1BP3 provides a novel mechanism of negative autophagy regulation through membrane lipids. Autophagy. 2017;13(5):779–780. https://doi.org/10.1080/15548627.2017.1305534.
Holland P, Knævelsrud H, Søreng K, Mathai BJ, Lystad AH, Pankiv S, et al. HS1BP3 negatively regulates autophagy by modulation of phosphatidic acid levels. Nat Commun. 2016;7:13889. https://doi.org/10.1038/ncomms13889.
Shi T, Xie J, Xiong Y, Deng W, Guo J, Wang F, et al. Human HS1BP3 induces cell apoptosis and activates AP-1. BMB Rep. 2011;44(6):381–6. https://doi.org/10.5483/BMBRep.2011.44.6.381.
Lind MV, Lauritzen L, Kristensen M, Ross AB, Eriksen JN. Effect of folate supplementation on insulin sensitivity and type 2 diabetes: a meta-analysis of randomized controlled trials. Am J Clin Nutr. 2019;109(4):1233. https://doi.org/10.1093/ajcn/nqz021.
Casanueva E, Drijanski A, Fernández-Gaxiola AC, Meza C, Pfeffer F. Folate deficiency is associated with obesity and anemia in Mexican urban women. Nutr Res. 2000;20(10):1389–1394. https://doi.org/10.1016/S0271-5317(00)80020-2.
Zhao JV, Schooling CM, Zhao JX. The effects of folate supplementation on glucose metabolism and risk of type 2 diabetes: a systematic review and meta-analysis of randomized controlled trials. Ann Epidemiol. 2018;28(4):249–257.e1. https://doi.org/10.1016/j.annepidem.2018.02.001.
Yang L, He Y, Kong Q, Zhang W, Xi D, Mao H, et al. Isolation, nucleotide identification and tissue expression of three novel ovine genes-SLC25A4, SLC25A5 and SLC25A6. Mol Biol Rep. 2010;37(6):2743–8. https://doi.org/10.1007/s11033-009-9812-z. Epub 2009 Sep 11.
Palmieri F. Mitochondrial transporters of the SLC25 family and associated diseases: a review. J Inherited Metab Dis. 2014;37(4):565–75. https://doi.org/10.1007/s10545-014-9708-5.
Palmieri F. The mitochondrial transporter family SLC25: identification, properties and physiopathology. Mol Asp Med. 2013;34(2-3):465–84. https://doi.org/10.1016/j.mam.2012.05.005.
Tan BG, Wellesley FC, Savery NJ, Szczelkun MD. Length heterogeneity at conserved sequence block 2 in human mitochondrial DNA acts as a rheostat for RNA polymerase POLRMT activity. Nucleic Acids Res. 2016;44(16):7817–29.
Alonso-Montes C, Castro MG, Reguero JR, Perrot A, Özcelik C, Geier C, et al. Mitochondrial transcription factors TFA, TFB1 and TFB2: a search for DNA variants/haplotypes and the risk of cardiac hypertrophy. Dis Markers. 2008;25(3):131–9.
Sánchez-Ferrero E, Coto E, Blázquez M, Ribacoba R, Guisasola LM, Salvador C, et al. Mutational screening of the mitochondrial transcription factors B1 and B2 (TFB1M and TFB2M) in Parkinson’s disease. Park Relat Disord. 2009;15(6):468–70.
The Genotype-Tissue Expression (GTEx) Project was supported by the Common Fund of the Office of the Director of the National Institutes of Health. Additional funds were provided by the NCI, NHGRI, NHLBI, NIDA, NIMH, and NINDS. Donors were enrolled at Biospecimen Source Sites funded by NCI\SAIC-Frederick, Inc. (SAIC-F) subcontracts to the National Disease Research Interchange (10XS170), Roswell Park Cancer Institute (10XS171), and Science Care, Inc. (X10S172). The Laboratory, Data Analysis, and Coordinating Center (LDACC) was funded through a contract (HHSN268201000029C) to The Broad Institute, Inc. Biorepository operations were funded through an SAIC-F subcontract to Van Andel Institute (10ST1035). Additional data repository and project management were provided by SAIC-F (HHSN261200800001E). The Brain Bank was supported by a supplements to University of Miami grants DA006227 & DA033684 and to contract N01MH000028. Statistical Methods development grants were made to the University of Geneva (MH090941 & MH101814), the University of Chicago (MH090951, MH090937, MH101820, MH101825), the University of North Carolina - Chapel Hill (MH090936 & MH101819), Harvard University (MH090948), Stanford University (MH101782), Washington University St Louis (MH101810), and the University of Pennsylvania (MH101822).
“The Framingham Heart Study is conducted and supported by the National Heart, Lung, and Blood Institute (NHLBI) in collaboration with Boston University (Contract No. N01-HC-25195 and HHSN268201500001I). This manuscript was not prepared in collaboration with investigators of the Framingham Heart Study and does not necessarily reflect the opinions or views of the Framingham Heart Study, Boston University, or NHLBI. Funding for SHARe Affymetrix genotyping was provided by NHLBI Contract N02-HL64278. SHARe Illumina genotyping was provided under an agreement between Illumina and Boston University. Funding for Affymetrix genotyping of the FHS Omni cohorts was provided by Intramural NHLBI funds from Andrew D. Johnson and Christopher J. O’Donnell.”
The bioinformatic analysis was performed using Prometheus (AGH, Krakow, Poland), Michigan Imputaton Server (Michigan, MI, USA) and Ohio Supercomputer Center (Columbus, OH, USA).
The study was funded by the National Science Centre in Poland through the Sonata Grant “Search for genetic variants influencing gestational weight gain in type 1 diabetes patients by genome wide association method” to ALS (Nr 2013/11/D/NZ5/03219).
Ethics approval and consent to participate
This study was approved by the Bioethical Committees of the Jagiellonian University and Poznan University of Medical Sciences and performed according to the Helsinki Declaration. Written informed consent was collected from all patients.
Consent for publication
The authors declare they have no competing interests.
Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.
List of mitochondrial variants included in the analysis in T1DM cohort. b. List of mitochondrial variants included in the analysis in FHS cohort.
All interactions of MT-ND2 gene with MitoCarta genes in T1DM cohort. b. Interactions of MT-ND2 gene significant after all corrections of MT-ND2 variant rs28357980 in FHS (MAF > 0.05, FDR < 0.05).
List of most significant associations of MT-ND2 gene.
All interactions of MT-ND2 gene genome wide in T1DM cohort. b. All interactions of MT-ND2 gene genome wide in FHS cohort. c. Significant associations of MT-ND2 variants in GWAS Catalog.
MitoGWAS on BMI type 1 diabetes patients and additive interactions between mitochondrial and nuclear variants in T1DM patients and FHS cohort.
About this article
Cite this article
Ludwig-Słomczyńska, A.H., Seweryn, M.T., Kapusta, P. et al. Mitochondrial GWAS and association of nuclear – mitochondrial epistasis with BMI in T1DM patients. BMC Med Genomics 13, 97 (2020). https://doi.org/10.1186/s12920-020-00752-7
- Mitochondrial-nuclear interactions