Introduction {#Sec1}
============

Type 2 diabetes mellitus (T2DM) has become one of the greatest challenges for public health worldwide in the 21st century^[@CR1]^. Diabetic complications, including cardiovascular disease, stroke, kidney disease and even cancer, contribute to about 4 million deaths per year globally^[@CR2]^. According to International Diabetes Federation, the number of diabetic patients could possibly increase from 366 million in 2011 to 552 million by 2030^[@CR3]^. In Xinjiang Uyghur Autonomous Region, the prevalence of T2DM is significantly higher in Uyghur people than in other ethnic minorities of close phylogenetic relationship^[@CR4]^. At the same body mass index (BMI) level, Uyghurs have greater risks for diabetes compared to other ethnic minorities in this area^[@CR5]^. However, the factors attributed to this high prevalence of T2DM in Uyghur people remain unexplored.

T2DM is a progressive multifactorial disease, characterized by insulin resistance and reduction in insulin secretion^[@CR6]^. Insulin resistance has developed many years before the diagnosis of T2DM^[@CR7]^. Mitochondrial DNA (mtDNA) plays an important role in insulin resistance and T2DM. Mitochondrial dysfunctions, such as reduced capacity to mitochondrial oxidative phosphorylation, submaximal ADP-stimulated oxidative phosphorylation and mitochondria plasticity, have been implicated in the pathogenesis of insulin resistance and diabetes^[@CR7]--[@CR10]^. The effects of mitochondrial haplogroups and SNPs on diabetes have also been explored. In Japanese and Korean populations, the haplogroup N9a was indicated to possibly confer resistance against T2DM^[@CR11]^. In contrast, in Chinese population, N9a is significantly associated with an increased risk of diabetic nephropathy^[@CR12]^. D4 haplogroup had borderline resistance to T2DM in Chinese Han individuals, whereas B haplogroup associated with an increased risk of T2DM^[@CR13]^. In Finnish, m.16189 T \> C (D-loop) and m.3010 G \> A (16S rRNA) was found to be more frequent in patients with diabetes than in controls^[@CR14]^. In Chinese Han population, m.3243 A \> G (16S rRNA) was associated with high prevalence of T2DM^[@CR15]^.

T2DM is a multi-factorial common disease affected by both genetic and environmental factors. For example, racial/ethnic factor can significantly modify the prevalence of T2DM. Studies performed in different populations often exhibit diverse even conflicting results due to different health behaviors, genetic factors and individual-level sociodemographics^[@CR16]^. To investigate the association between mtDNA and T2DM in Uyghur people, we conducted a hospital-based case-control study on 210 Uyghur individuals including 88 T2DM patients and 122 controls.

Results {#Sec2}
=======

Mitochondrial genome analysis {#Sec3}
-----------------------------

According to mitochondrial phylogenetic tree from Phylo-tree^[@CR17]^, we classified the samples into 12 haplogroups (A, B, C, D4, F, G, H, HV, J, T, U and W), haplogroups (D5, I, X, Z, O, W, K, M7, M8, M9, N30) with a limited number (smaller than 6, 3%) were clustered to group Others in M and group Others in N. Haplogroups A, B, C, D4, F and G were common in Asians but rare in Europeans, while haplogroups H, HV, J, T, U and W were very rare in Asians but common in Europeans. This result was consistent with human mitochondrial haplogroup distribution formed by the Human migration ([www.mitomap.org](http://www.mitomap.org)), suggesting that Uyghur people was an admixture of populations from Europe and Asia. Therefore, subject clustering was required before investigating potential association between mtDNA and T2DM. PCA analysis was performed and all 210 subjects were classified into 4 clusters (Fig. [1](#Fig1){ref-type="fig"}). The distribution of haplogroups in each cluster were shown in Fig. [2](#Fig2){ref-type="fig"}.Figure 1Subjects distribution according to the PCA analysis. We employed Smartpca from EIGENSOFT package and K-means method from Weka to classify all 210 subjects into 4 clusters. Different colors indicated different clusters. Triangle indicates case subjects. Crossover indicates control subjects. Principal component 1 (PC1) and principal component 2 (PC2) were the first two axes of PCA and represented 7.626% and 5.276% variance. Figure 2Distribution of haplogroups in 4 different clusters. Haplogroups with a sample size smaller than 6 (3%) are allocated into groups 'OtherI' and 'OtherII'. Different colors indicated different clusters.

Using multivariate logistic regression analysis, we found that people harboring haplogroup H (odds ratio \[OR\] = 1.40; 95% confidence interval \[CI\]: 1.20--1.64; *P* = 0.0005138) and D4 (odds ratio = 1.47; 95% CI: 1.22--1.77; *P* = 0.001064) was more frequent in patients with diabetes than in controls (Table [1](#Tab1){ref-type="table"}).Table 1Multivariate logistic regression analysis of mitochondrial haplogroups.HaplogroupT2DM (n = 88)Control (n = 122)Total (n = 210)Adjusted P-value\*OR (95%CI)A16710.97(0.72--1.30)B791611.14(0.93--1.39)C5141910.99(0.83--1.20)D4136190.0010641.47(1.22--1.77)F57120.281.32(1.04--1.67)G87150.3781.27(1.03--1.57)H1912310.00051381.40(1.20--1.64)HV\#0101010.86(0.66--1.10)J35810.98(0.74--1.29)T571211.16(0.92--1.46)U68140.52411.26(0.92--1.57)W07710.84(0.63--1.14)Others in M5152011.10(0.92--1.33)Others in N118190.0111581.39(1.14--1.68)Symbol \*indicates that adjusted *P*-values were adjusted by Bonferroni correction. Haplogroups (D5, I, X, Z, O, W, K, M9) with a sample size smaller than 6 (3%) were allocated into group 'Others in M' and group 'Others in N' according to the phylogenetic tree. T2DM: type II diabetes mellitus. HV\#: HV haplogroup not included H haplogroup.

In addition, the markers related to haplogroups D4 and H, including variants m.3010 G \> A (D4), m.8414 C \> T (D4), m.14668 C \> T (D4), m.2706 G \> A (H) and m.7028 T \> C (H) were more frequent in cases than in controls (Table [2](#Tab2){ref-type="table"}), by Fisher's exact test. Of these variants, m.8414 C \> T from the haplogroup D4 was non-synonymous and located in mitochondrially encoded ATP synthase protein 8 (*MT-ATP8*). The m.2706 G \> A and m.3010 G \> A were located in mitochondrially encoded 16S RNA.Table 2Analysis of SNPs related to haplogroup D4 and H in 4 clusters※.SNPhaplogroupRegionAA changeOR95% CIAdjusted *P*-value\*2706A^1^H16S RNA6.491.89--22.330.0041723010A^2^D416S RNA5.941.52--23.180.036247028C^1^H*MT-CO1*8.712.40--31.550.00088328414T^2^D4*MT-ATP8*L17F5.941.52--23.180.0362414668T^2^D4*MT-ND6*5.941.52--23.180.03624※m.8414 T were nonsynonymous. Other 2 SNPs (m.2706 A, m.3010 A) from haplogroup H and D4 might change the structure of mitochondrial 16S RNA. AA: amino acid. Symbol \*indicated that adjusted *P*-values were corrected by Bonferroni correction. Superscript in the first column indicated a specific cluster that a specific SNP belonged to.

Variant damage analysis {#Sec4}
-----------------------

We then extracted unique non-synonymous variants in cases (Table [3](#Tab3){ref-type="table"}) to explore their potential functions. Four popular predictors were applied to assess the impact of AA substitution caused by these variants on the properties of mitochondrial proteins. The predictors were PolyPhen-2^[@CR18]^, CADD^[@CR19]^, PANTHER^[@CR20]^ and CAROL^[@CR21]^.Table 3Unique SNPs from subjects with T2DM.SNPFrequencyRegionAA change7673 G3*MT-CO2*30,I \> V4924 A3*MT-ND2*152,S \> N9966 A3*MT-CO3*254,V \> I12362 G3*MT-ND5*9,T \> S13780 G3*MT-ND5*482,I \> V8961 G2*MT-ATP8*146,T \> A4012 G2*MT-ND1*236,T \> A13637 G2*MT-ND5*434,Q \> R10084 C2*MT-ND3*9,I \> T13204 A2*MT-ND5*290,V \> I5301 G2*MT-ND2*278,I \> V7859 A2*MT-CO2*92,D \> NAA: amino acid.

Interestingly, m.8414 C \> T was highly suspicious. Three of four predictors suggested that m.8414 C \> T was likely to influence the stability of the protein and may contribute to disease development (Table [4](#Tab4){ref-type="table"}).Table 4Protein properties analysis results.SNPPrediction ResultsPolyphen-2CADDCAROLPANTHER8414,C \> TDDDX4924,G \> A^†^XXXD7673,A \> G^†^XXXX9966,G \> A^†^XXXX12362,C \> G^†^XXXD13780,A \> G^†^XXXD8961,A \> G\*DXDX4012,A \> G\*XXXX5301,A \> G\*XXXX7859,G \> A\*XXXX10084,T \> C\*XXXX13204,G \> A\*XXXX13637,A \> G\*XXXXD represented damaging, deleterious or disease in Polyphen-2, CADD, CAROL and PANTHER. X represented benign and neutral. ^†^Indicated unique SNP existing in case group with a count of 3. \*Indicated unique SNP existing in case group with a count of 2.

Protein 3D structure analysis {#Sec5}
-----------------------------

We performed 3D structure analysis using CABS-fold^[@CR22]^ to assess the effect of variant p.Leu17Phe (m.8414 T \> G) on the structure of *MT-ATP8*. The 3D structure of *MT-ATP8* was displayed in Supplementary Figure 1.

In term of molecular modeling, potential energy of a system of atoms is calculated as force field, in which angles, bonds, torsions and impropers are taken into consideration. In current study, GROMOS96 force field was used to estimate the potential energy of *MT-ATP8* by Swiss-PDBviewer^[@CR23]^ (version 4.10). When *MT-ATP8* 17^th^ amino acid (AA) was mutated from Leu to Phe by Swiss-PDBviewer, a decrease in H-bond number indicated a more unstable structure (Fig. [3](#Fig3){ref-type="fig"}). Using Swiss-PDBviewer, the force field around 17^th^ AA was estimated to be increased from −40.646 to −12.161 (Fig. [4](#Fig4){ref-type="fig"}), which may lead to unstable protein structure around this area.Figure 33D structure of *MT-ATP8* protein. (**a**) *MT-ATP8* with p.17 L, red sticks are 17^th^ amino acid Leu; (**b**) *MT-ATP8* with p.17 F, red sticks are 17^th^ amino acid Phe. Green dotted line: H-bonds. Figure 4Force field curve around the 17^th^ AA between 17 F and 17 L. Grey areas indicate higher energy caused by *MT-ATP8* p.17 F. Black or line indicates the force field of *MT-ATP8* with p.17 L or p.17 F, respectively.

Protein transmembrane property analysis {#Sec6}
---------------------------------------

The subcellular location of ATP synthase is on the mitochondrial inner membrane, and the ATP synthase protein 8 is a single-pass membrane protein, according to annotations from Uniprot. So it is necessary to investigate whether the mutation p.L17F (m.8414 C \> T) in the *MT-ATP8* has an impact on the protein's transmembrane area.

TMpred was used to assess protein sequence with a mutation p.L17F in the *MT-ATP8*. The result indicated that the transmembrane area of *MT-ATP8* was changed from 8^th^-24^th^ AAs to 11^th^-29^th^ AAs, which might affect the assembly of ATP synthase (Fig. [5](#Fig5){ref-type="fig"}).Figure 5Possible transmembrane helices change in *MT-ATP8*. Red or blue vertical line indicates the transmembrane region of *MT-ATP8* with m.8414 T or m.8414 C, respectively. Only scores above 500 were considered reliable by TMpred.

RNA structure analysis {#Sec7}
----------------------

The rRNA plays an important role in human mitochondrial ribosome assembly and RNA modification, which are strongly linked with energy metabolism. RNA-fold server was used to predict the structural change caused by variant m.2706 G \> A (Fig. [6](#Fig6){ref-type="fig"}). Oligowalk from RNAstructure^[@CR24]^ (version 5.3) was employed to compare the target binding ability between 16S rRNA with and without m.2706 G \> A. Notably, the structural change caused by variant m.2706 G \> A from haplogroup H might weaken the binding between mitochondrial 16S RNA and its nucleic acid targets (Fig. [7](#Fig7){ref-type="fig"}). Therefore, structural change of mitochondrial 16S rRNA caused by variant m.2706 G \> A might affect RNA modifications and ribosome assembly, in Uyghur population.Figure 6Secondary structure of 16S RNA. (**A**) Secondary structure of 16S RNA with m.2706 G. (**B**) Secondary structure with m.2706 G \> A. Figure 7Oligo-Target binding curve. Green or red curve indicates the score of 16S rRNA with m.2706 G or m.2706 A, respectively. A more negative value indicates a tighter binding. Dashed line shows the position corresponding to m.2706 G \> A.

Discussion {#Sec8}
==========

In this case-control study, we have identified an association between mitochondrial DNA and T2DM in Chinese Uyghur population. Mitochondrial DNA sequencing was conducted in 88 T2DM patients and 122 controls. Haplogroups H and D4 were found to be more frequent in T2DM patients than in controls.

A variety of diseases, such as Alzheimer's disease, cardiovascular disease, Leber's hereditary optic neuropathy (LHON) and T2DM have been associated with mitochondrial haplogroups^[@CR25]--[@CR28]^. It was reported that racial factors played an important role in the prevalence of T2DM^[@CR16]^. Researches performed in different populations presented different even conflicting results. For Han Chinese, haplogroup N9a increased the risk of T2DM complications^[@CR12]^, while in Japanese population, N9a was associated with resistance against T2DM^[@CR11]^ or against metabolic syndrome in Japanese women^[@CR29]^. In Chinese Taiwan population, haplogroup D4 was related to a borderline reduced risk for T2DM^[@CR13]^. In contrast, haplogroup D4 was associated with an increased risk of T2DM in Korean population^[@CR11]^ and Chinese Uyghur population in present research.

In our study, haplogroup D4 was more frequent in Chinese Uyghur T2DM patients, which was consistent with the findings in Finnish population, showing that variant m.3010 G \> A was more frequent in cases than in controls. The variant m.3010 G \> A marked haplogroup D4 in Chinese Uyghur subjects and haplogroup H1 in Finnish population. This indicated that m.3010 G \> A variant was likely to be involved in the pathogenesis of T2DM. However, its biologic functions in the pathogenesis of T2DM require further investigation.

Haplogroup H was found to be associated with enhanced mitochondrial activities, including higher mitochondrial respiratory capacity (intrinsic activity)^[@CR30]^, mitochondrial oxidative damage^[@CR31]^ and intracellular ATP concentration^[@CR32]^, while its sub-haplogroup H1 was associated with T2DM in Finnish population. In present study, haplogroup H showed a significant association with T2DM. Variants m.2706 G \> A and m.7028 T \> C were related to increased incidence of T2DM in Uyghur people harboring haplogroup H. The SNP m.2706 G \> A was located in mitochondrially encoded 16S RNA, while m.7028 T \> C was synonymous.

The structure of mitochondrially encoded rRNA is critical for human mitochondrial ribosome assembly^[@CR33]^. RNA modifications also affect energy metabolism^[@CR34]^. Structural changes of mitochondrial 16S rRNA might regulate RNA modifications, protein synthesis as well as ribosome assembly. In present study, a structural change in 16S rRNA caused by m.2706 G \> A from haplogroup H was observed, and thus the binding ability of mutated 16S rRNA might be altered. This structural change might affect 16S rRNA modifications and mitochondrial encoded protein synthesis, which are important for energy metabolism.

In this study, we found that haplogroup D4 and H were more frequent in Chinese Uyghur T2DM patients. Variants m.2706 G \> A and m.8414 C \> T might be risk factors in the development of T2DM in Chinese Uyghur population.

Materials and Methods {#Sec9}
=====================

Samples and mitochondrial genome sequencing {#Sec10}
-------------------------------------------

The healthy controls and cases, residents of Heitan, Keshi, and Korla districts in Xinjiang Uyghur Autonomous Region, were recruited from May 2010 to November 2011 at the Fifth Affiliated Hospital of Xinjiang Medical University. A written informed consent has been obtained from each subject. This study has been approved by Ethics Committee at the Fifth Affiliated Hospital of Xinjiang Medical University. All procedures were conducted according to the regulations and guidelines approved by the ethics committee at the Fifth Affiliated Hospital of Xinjiang Medical University.

Demographic and clinical information was extracted from electronic medical record, including age, gender, ethnicity, medical history, smoking and alcohol drinking. Smoking was defined as a current smoker or a former smoker (who quits) with tobacco/cigarettes consumption of more than 10 packs per year. Alcohol drinking was defined as alcohol consumption of at least 50 g per week for at least 6 months. Body weight and height were measured with light clothing and in bare feet. Blood samples were collected after overnight fasting. Glucose, Total Cholesterol, LDL, HDL and Triglyceride were measured with standard biochemical assays. T2DM was diagnosed according to WHO 1999 Diagnostic Criteria. Hypertension was defined as systolic pressure ≥140 mmHg and/or diastolic pressure ≥90 mmHg (from an average of 2 measurements or pre- and post-therapy). Hyperlipidemia was defined as fasting venous total cholesterol ≥200 mg/dL or LDL-C ≥ 130 mg/dL.

A total of 210 Uyghur objects, including 88 subjects with T2DM and 122 control subjects, were recruited in this study. Genomic DNA was extracted from peripheral whole blood using a DNA Blood Mini Kit 250 (QIAGEN Translational Medicine Co., Ltd., German). Fragments of mitochondrial DNA were amplified by PCR. The reaction mixture (100 μl) contained 20 ng of genomic DNA, 30 pmol of each primer, 5 mmol of each dNTP, and 2.5 U of DNA polymerase in a 10× PCR buffer. The cycling conditions were as follows: 94 °C for 5 min, followed by 35 cycles of denaturation at 94 °C for 30 sec, annealing at 65 °C for 40 sec, and extension at 72 °C for 4 min, with a final extension at 72 °C for 10 min. The fragments were purified and sequenced by Applied Biosystems 3730 DNA automatic sequencer at the Beijing Genomics Institute (Beijing, China). Primers for sequencing and PCR were shown in Supplementary Tables [1](#MOESM1){ref-type="media"} and [2](#MOESM1){ref-type="media"}, respectively. The sequence data reported in this paper have been deposited in the Genome Sequence Archive of Beijing Institute of Genomics, Chinese Academy of Sciences. The information about age, gender and BMI of samples were shown in Supplementary Table [3](#MOESM1){ref-type="media"}.

Mitochondrial haplogroups classification {#Sec11}
----------------------------------------

The resultant sequence data are aligned against with the Cambridge reference sequence (GenBank No., NC_012920). We determine the haplogroups of mitochondrial DNA according to the mtDNA tree^[@CR17]^ and Mitomaster^[@CR35]^.

Statistical analysis {#Sec12}
--------------------

We applied smartpca program from the EIGENSOFT version 3.0 to analyze the PCA and used K-means method from software Weka to cluster the PC1, PC2 and PC3 from PCA analysis.

We performed multivariate logistic regression analyses to demonstrate the association between the haplogroups and T2DM by glm function from R package (version 3.3.2), considering T2DM as a dependent variable, whereas information on haplogroups, the first two PCs from PCA, BMI, age and sex of each subject as independent variables. In multivariate logistic regression analyses, haplogroups with sample size smaller than 6 (3%) are allocated into group 'Others in M' and group 'Others in N' according to the phylogenetic tree. All *P*-values are corrected by Bonferroni method.

We used Fisher's Exact Test to compare the frequencies of SNPs between cases and controls within each cluster. We calculated the *P*-values by two-tailed Fisher's Exact Tests from R version 3.3.2, along with odds ratio (OR) and 95% confidence interval (CI). All *P* values were adjusted by Bonferroni correction. For Bonferroni correction, *P* values are adjusted according to the amount of haplogroups and SNPs in different clusters.

Variant damage analysis {#Sec13}
-----------------------

We used a series of predicting tools, including PolyPhen-2^[@CR18]^, CADD^[@CR19]^, PANTHER^[@CR20]^ and CAROL^[@CR21]^, to analyze a specific protein with special mutations. These tools focus on different areas, disease related SNPs, including conservation in the phylogenetic tree and structure features in transmembrane, intramembrane, signal peptide, binding region, etc. A variant was likely to be suspicious as disease-causing if more than 2 tools suggest deleterious.

Protein 3D structure analysis {#Sec14}
-----------------------------

We used CABS-fold server^[@CR22]^ to predict the 3D structure of MT-ATP8 protein, due to lack of available X-ray 3D structure of human MT-ATP8 molecule. Afterwards, we used SWISS-Pdbviewer (version 4.10) to analyze mutations and to estimate force field.

RNA structure analysis {#Sec15}
----------------------

The RNA-fold sever (<http://rna.tbi.univie.ac.at/cgi-bin/RNAfold.cgi>) was employed to predict the secondary structure of mitochondrially encoded 16S RNA. We adopted the minimum free energy prediction model. Oligowalk algorithm from RNAstructure program^[@CR24]^ was employed to assess the target binding ability.

Data availability {#Sec16}
-----------------

The sequence data reported in this paper have been deposited in the Genome Sequence Archive of Beijing Institute of Genomics, Chinese Academy of Sciences, with an accession number CRA000365.

Electronic supplementary material
=================================

 {#Sec17}

Supplementary information

Wenxi Jiang and Ronghui Li contributed equally to this work.

**Electronic supplementary material**

**Supplementary information** accompanies this paper at 10.1038/s41598-017-17086-7.

**Publisher\'s note:** Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

This work was supported by grants from the Natural Science Foundation of Xinjiang Autonomous Region (2010211A33, 2012211A033) and the National Natural Science Foundation of China (81260057).

J.Y., W.J. and M.G. coordinated and managed the project. P.W., T.W., J.L. conducted the experimental work. R.L. and Y.Z. conducted data analysis. R.L., Y.Z., W.J. and M.G. participated in manuscript writing and editing. All authors have read and approved the final version of the manuscript. We declare that this work, as original research conducted in China, in whole or in part has never been published previously, and is not under consideration for publication elsewhere.

Competing Interests {#FPar1}
===================

The authors declare that they have no competing interests.
