Abstract
Our study aimed to investigate if genetic variants around 16p13.3’s HBA1 locus, associated with erythrocyte indices and HbA1c levels, predict α-thalassemia-related erythrocyte indices, cardiometabolic traits, and diabetes risk in Taiwanese individuals. We analyzed Taiwan Biobank data, including whole-genome sequencing from 1,493 participants and genotyping arrays from 129,542 individuals. First, we performed regional association analysis using whole-genome sequencing data to identify genetic variants significantly associated with erythrocyte indices, confirming their linkage disequilibrium with the α0 thalassemia --SEA deletion mutation, a common cause of α-thalassemia in Southeast Asian populations. Deletion mutation sequencing further validated these variants’ association with α-thalassemia. Subsequently, we analyzed genotyping array data, revealing associations between specific genetic variants and cardiometabolic traits, including lipid profiles, HbA1c levels, bilirubin levels, and diabetes risk. Using Mendelian randomization, we established causal relationships between α-thalassemia-related erythrocyte indices and cardiometabolic traits, elucidating their role in diabetes susceptibility. Our findings highlight genetic variants around the α-globin genes as surrogate markers for common α-thalassemia mutations in Taiwan, emphasizing the causal links between α-thalassemia-related erythrocyte indices, cardiometabolic traits, and heightened diabetes risk.
Introduction
Microcytic anemia is the most common type of anemia in both childhood and adulthood (Richardson, 2007; DeLoughery, 2014; Urrechaga et al, 2015; Zanetti et al, 2021). Microcytic anemia is characterized by the presence of small-sized erythrocytes because of the decreased availability of hemoglobin components, such as globins, iron, and heme; this, in turn, causes a reduction in the hemoglobin content in RBC precursors, subsequently resulting in delayed erythroid differentiation (DeLoughery, 2014; Cappellini et al, 2020). Microcytic anemia may be caused by genetic or nongenetic factors, and three subgroups of inherited microcytic anemia have been classified based on the causative defects: hemoglobinopathies (thalassemia) caused by defects in globin chain synthesis, defects in heme synthesis, and defects in iron availability or acquisition by erythroid precursors (DeLoughery, 2014; Cappellini et al, 2020). Although most of the genetic causes of microcytic anemia are rare in terms of occurrence and each mutation contributes to a small proportion of patients with microcytic anemia, some mutations may be dominant in distinct ethnic populations, such as α0 thalassemia deletion --SEA, --Fil, and --THAI in the Southeast Asian region (Chen et al, 2002; Iolascon et al, 2009; Harteveld & Higgs, 2010; Chao et al, 2014; Lee et al, 2015; Wang et al, 2017; Farashi & Harteveld, 2018; Mettananda & Higgs, 2018; Cappellini et al, 2020).
Thalassemia is among the most common human monogenic disorders (Weatherall & Clegg, 1996). α-thalassemia, which is caused by defective α-globin gene mutations on chromosome 16p13.3, is a hemoglobinopathy characterized by deficits in α-globin chain synthesis and microcytic hypochromic erythrocytes (Harteveld & Higgs, 2010; Higgs et al, 2012; Piel & Weatherall, 2014; Mettananda & Higgs, 2018; Taher et al, 2018). Recently, Mendelian randomization (MR) studies have reported that genetically predicted erythrocyte indices are causally related to cognitive function deficits, Alzheimer’s disease, and venous thromboembolism (Winchester et al, 2018; Luo et al, 2020). Multiple erythrocyte indices were identified to be associated with HbA1c levels and diabetes mellitus (DM) development (Wang et al, 2021). Genome-wide association studies (GWASs) have revealed that chromosome 16p13.3 as a gene locus is associated with erythrocyte traits and hemoglobin A1c levels (Kichaev et al, 2019; Vuckovic et al, 2020; Sinnott-Armstrong et al, 2021). However, causal relationships between α-thalassemia-related erythrocyte indices and cardiometabolic traits remain to be completely elucidated. The Taiwan Biobank (TWB) is a large-scale population-based cohort study including individuals aged between 30 and 70 yr who have no history of cancer (Fan et al, 2008; Chen et al, 2016). More than 120,000 participants in the TWB study underwent whole-genome genotyping array analysis, and a subgroup of them also underwent whole-genome sequencing (WGS) analysis. Our study was primarily designed to uncover and elucidate the causal relationships between α-thalassemia-related erythrocyte indices, cardiometabolic traits, and DM. We initiated by conducting a regional association analysis within the 16p13.3 region and performed deletion mutation sequencing on WGS data to pinpoint genetic variants serving as proxies for α-thalassemia. Subsequently, our investigation expanded to encompass genotyping array data from an extensive cohort of 129,542 participants. This enabled us to comprehensively explore associations between specific genetic variants and a spectrum of cardiometabolic traits, along with assessing the risk of DM. In the final phase, we harnessed MR to estimate the causal effect of α-thalassemia-related erythrocyte indices on various cardiometabolic traits and DM, thereby enhancing our understanding of these intricate interactions.
Results
Baseline characteristics of TWB participants
Fig 1 depicts the flowchart of participant enrollment with whole-genome genotyping array and WGS data. Table 1 summarizes the demographic data, clinical and biochemical data, and lipid profiles of TWB participants with WGS and GWAS array data. The average ages of participants selected for WGS and GWAS array analysis were 50 and 51 yr of age, respectively. The numbers of men and women were equal for participants with WGS data, whereas the number of women was higher for participants with array data.
Regional association study for RBC parameters
We analyzed the association between genetic variants at positions between 0.06 and 0.68 Mb on chromosome 16p13.3 and the RBC parameters in TWB participants with WGS data. Our findings revealed that three single-nucleotide variations (SNVs), namely NPRL3 rs191086839, LUC7L rs372755452, and PGAP6 rs375498857, were specific to Asians (with minor allele frequencies ranging from 0.0177 to 0.0218 among TWB and East Asian populations versus all <0.0001 among other ethnic populations) (Table S1) and significantly associated with four RBC traits, namely the erythrocyte count, Hb level, mean corpuscular volume (MCV), and mean corpuscular hemoglobin (MCH), with the lowest P-values for each variant being 8.68 × 10−93, 2.30 × 10−17, 1.98 × 10−101, and 1.76 × 10−122, respectively (Fig 2). MCV and MCH were thus selected for in-depth examination with cardiometabolic traits based on their most significant association with 16p13.3 and their impact on erythrocyte count and Hb level. The LD map revealed strong LD among these three variants (Fig 3) in our study population.
Association of the three variants on chromosome 16p13.3 with RBC parameters
In total, the data of 1,493 TWB volunteers were included in the study on the association of genotypes and phenotypes with RBC parameters (Table 2 and Fig 4). By performing general linear regression by using an additive model, we determined that after adjustment for age, sex, body mass index (BMI), and smoking status, participants harboring the minor allele of the three variants (i.e., the C allele for NPRL3 rs191086839, a single-nucleotide deletion for LUC7L rs372755452, and the A allele for PGAP6 rs375498857) exhibited a genome-wide significant association (P < 5 × 10−8) with a higher erythrocyte count and lower Hb, MCV, MCH, and mean corpuscular hemoglobin concentration values and a subthreshold significant association (P between 5 × 10−8 and 1 × 10−4) with lower HCT values (Table 2). Furthermore, participants harboring the minor allele of each variant had significantly higher risks of the microcytic hypochromic trait, anemia, and microcytic anemia (Table 2 and Fig 4A–I). Stepwise linear regression analysis for MCV and MCH revealed that the association with these erythrocyte traits markedly diminished for the PGAP6 rs375498857 genotype, but not for the other two variants (Table 3). These results suggest that the association between PGAP6 rs375498857 and the erythrocyte traits is because of the strong LD with NPRL3 rs191086839 and LUC7L rs372755452. Through a logistic regression analysis of all the three variants, we determined a significant association of NPRL3 rs191086839 and LUC7L rs372755452 with the risk of the microcytic hypochromic trait (odds ratio [OR] = 40.87 and 20.2; 95% confidence interval [CI] = 7.26–229.99 and 2.13–191.3, respectively; Table S2) and of NPRL3 rs191086839 with the risk of microcytic anemia (ORs = 132.74; 95% CIs = 11.80–1,493.49). These results revealed that the NPRL3 rs191086839 variant was independently associated with all the erythrocyte traits analyzed.
Searching for deletion mutations in the α-globin gene region in participants with WGS data: BreakDancer and PCR analysis with direct DNA sequencing
α-thalassemia is most likely caused by large deletion mutations in the α-globin gene (Fig 5A and B) (Richardson, 2007; DeLoughery, 2014; Urrechaga et al, 2015; Zanetti et al, 2021). Thus, we first used BreakDancer v1.3.6 to screen for common and large deletions in the α-globin gene cluster region in Taiwanese participants. When a deletion was present, the section of DNA that was absent in the participant’s genome was identified through a comparison with the reference genome. When pairs from a section of DNA spanning the deletion are aligned to the genome, the inferred insert size will be larger than expected. BreakDancer uses these pairs of reads to detect deletions. Thus, the sequencing depth affects the resolution of deletion positions. The average depth of the WGS in TWB participants was ∼30X. Thus, we did not anticipate to obtain numerous pairs of reads, which means that the detected location of the deletion may be approximate. In this study, we selected only deletions involving the HBA1 and HBA2 gene regions for analysis. A total of 128 participants were enrolled in BreakDancer analysis, including 20 with normocytic erythrocytes and 108 with microcytic hypochromic erythrocytes. 62 participants with microcytic hypochromic erythrocytes and one participant with normocytic hypochromic erythrocytes had a large deletion involving the α-globin gene region: 58 individuals had an ∼20-kb deletion, 4 had an ∼30-kb deletion, and 1 had an ∼33-kb deletion (Fig 5C–F). Because the α0 thalassemia --SEA 20-kb deletion in the α-globin gene is the most common cause of α-thalassemia in the Taiwanese population (Chen et al, 2002; Chao et al, 2014; Lee et al, 2015; Wang et al, 2017), we conducted tests to detect the presence of this deletion in 1,474 participants who had WGS data and genomic DNA available for analysis (Fig 5G). By performing PCR, we detected the α0 thalassemia --SEA 20-kb deletion in the α-globin gene in 60 participants, including 58 participants screened using BreakDancer (57 participants had an MCV of <80 fl and an MCH of <25 pg/RBC and one had an MCV of >80 fl and an MCH of <25 pg/RBC) and two participants who lacked BreakDancer data (one participant had an MCV of >80 fl and an MCH of <25 pg/RBC and another participant had an MCV of <80 fl and an MCH of >25 pg/RBC; Fig 5L). PCR findings confirmed that the 30-kb deletion detected by BreakDancer was because of the α0 thalassemia --FIL deletion and the 33-kb deletion was because of the α0 thalassemia --THAI deletion (Fig 5H). The PCR results were confirmed through direct DNA sequencing (Fig 5I–K), with the exception of the 30-kb deletion mutation, in which multiple single-nucleotide repeated sequences were detected in PCR products that caused difficulty in performing DNA sequencing.
Association between three significant SNVs on chromosome 16p13.3 and the α0 thalassemia --SEA deletion mutation detected through PCR
The α0 thalassemia --SEA deletion is the most common cause of α-thalassemia in Taiwan. We analyzed the association among the three significant SNVs and the α0 thalassemia --SEA deletion mutation detected through PCR in 1,474 TWB participants with WGS data (Fig 5M–O). For NPRL3 rs191086839 risk allele carrier, the sensitivity, specificity, positive predictive value (PPV), and negative predictive value (NPV) were 85.00%, 99.93%, 98.08%, and 99.37%, respectively. For LUC7L rs372755452 risk allele carrier, the sensitivity, specificity, PPV, and NPV were 90.00%, 99.79%, 94.74%, and 99.58%, respectively. For PGAP6 rs375498857 risk allele carrier, the sensitivity, specificity, PPV, and NPV were 85.00%, 99.65%, 91.07%, and 99.37%, respectively (Table S3).
Regional association analysis of cardiometabolic traits with chromosome 16p13.3
Among the three significant variants located on chromosome 16p13.3, only PGAP6 rs375498857 was used in genotype–phenotype analysis performed by employing GWAS imputation data from 115,926 TWB participants for all study phenotypes (Table S4). In brief, the imputed PGAP6 rs375498857 genotypes completely matched with genotypes derived from the WGS data of 1,493 TWB participants. Furthermore, the imputation data from the GWAS CHB-1 Array (TWBv1.0) and the GWAS CHB-2 Array (TWBv2.0) were consistent, with 1,206 participants undergoing analysis with both the arrays. The other two variants did not meet the aforementioned criteria for imputation and thus were not included in further analysis performed using the imputation data. In addition to RBC traits, genotype–phenotype analysis revealed a genome-wide significant association of the rs375498857 genotype with total high-density lipoprotein (HDL)– and low-density lipoprotein (LDL)–cholesterol levels and HbA1c and total bilirubin levels (P = 4.57 × 10−36, 1.92 × 10−12, 5.91 × 10−24, 5.13 × 10−40, and 2.25 × 10−14, respectively; Table S4). Furthermore, the rs375498857 genotype was significantly associated with the risk of DM (P = 2.00 × 10−6). By performing a regional association study, we determined that rs375498857 is the only variant associated with the aforementioned phenotypes and disease (Figs S1, S2, S3, S4, S5, S6, and S7).
Correlation between RBC parameters and cardiometabolic traits
Various metabolic traits were significantly associated with MCV and MCH, and most of the study traits had a P-value threshold of 10−4. Consistent associations were observed between MCV and MCH and nearly all studied traits, with MCV generally having a stronger effect. Elevated MCV and MCH were associated with increased uric acid levels and lipid profiles, including total, HDL and LDL cholesterol and triglyceride levels. Elevated MCV reduced the risk of metabolic syndrome and most of the metabolic syndrome-related components, such as systolic blood pressure, hypertension, fasting plasma glucose level, HbA1c level, DM status, eGFR, albuminuria, and microalbuminuria. All liver function-related test results showed the same direction with changes in MCV and MCH (Table 4).
MR study with two-stage least square (2SLS) for the association between α-thalassemia-related erythrocyte indices and cardiometabolic traits and DM using rs375498857 variant as an IV
For MR analyses, we selected rs375498857-related cardiometabolic traits (with a P-value threshold of 10−4). In the 2SLS IV analysis for the direction and causality of α-thalassemia-related erythrocyte indices and cardiometabolic traits, the association of PGAP6 rs375498857 genotypes with cardiometabolic traits remained significant even after adjustment for multiple parameters associated with MCV or MCH (Table 5): the association between the rs375498857 genotype and LDL cholesterol levels and DM subsided after adjustment for MCH. Moreover, the association between the rs375498857 genotype and total cholesterol and total bilirubin levels and HbA1c subsided after adjustment for MCV (all P > 0.05). The association between the rs375498857 genotype and HDL cholesterol levels did not totally abolish after either adjustment for MCH or MCV (P = 2.22 × 10−10 and P = 0.0008, respectively). These results suggest that the association between the rs375498857 genotype and cardiometabolic traits and DM, with the exception of HDL cholesterol levels, is dependent on MCH or MCV. Moreover, consistent findings were observed in the same 2SLS analysis, using either NPRL3 rs191086839 or LUC7L rs372755452 as the instrumental variable (Table S5).
Discussion
This study explored the association of chromosome 16p13.3 variants with erythrocyte indices in the Taiwanese population. Our findings revealed that three significant SNVs at this chromosome location around the HBA1 gene were closely associated with erythrocyte parameters, namely RBC counts, Hb, MCV, and MCH levels. TWB participants who had the minor allele of these three variants in their erythrocytes were both microcytic and hypochromic. Using BreakDancer v1.3.6 for screening, followed by PCR amplification and direct sequencing, we confirmed that the α-thalassemia deletion mutation --SEA, the most common α-thalassemia mutation in the Taiwanese population, exhibited strong LD with these SNVs, with their NPVs ranging from 0.994 to 0.996 and PPVs from 0.912 to 0.981. Thus, these SNVs can be considered as crucial surrogate genetic markers for the α-thalassemia deletion mutation --SEA. We performed an MR study by using PGAP6 rs375498857 as IV in participants with whole-genomic genotyping data. We observed causal relationships between MCV/MCH and rs375498857-related cardiometabolic traits and DM (Fig 3). This study is the first to demonstrate that α-thalassemia caused by the --SEA deletion mutation can be analyzed using SNVs as the surrogate marker, obviating the need to directly genotype the deletion mutation. This approach can be applied to a larger population using genotyping analyses, such as array data, for the mass screening of α-thalassemia in adults or newborns, making it a highly powerful and effective tool for epidemiological research on α-thalassemia in Taiwan.
Linking three significant SNVs located on chromosome 16p13.3 with a deletion mutation of α-thalassemia
Previous studies have reported that chromosome 16p13.3 as a gene locus is associated with erythrocyte traits (Kichaev et al, 2019; Vuckovic et al, 2020; Sinnott-Armstrong et al, 2021). By including TWB participants in their study, Lee et al (2022) revealed that the PGAP6 rs375498857 genotype is associated with several RBC traits. In our study, we observed that NPRL3 rs191086839, LUC7L rs372755452, and PGAP6 rs375498857 are associated with multiple erythrocyte indices, exhibit strong LD, and are specific to the Asian population. The minor allele frequencies of all the three variants were <0.0001 in European populations but were between 0.0188 and 0.0218 in Asian populations in the 1,000 Genome Project (Table S1). Our findings indicate that the three significant SNVs may be linked to the common deletion mutation of α-thalassemia in the Taiwanese population. First, α-thalassemia is characterized by a deficit in α-globin chain synthesis and most commonly caused by the deletion in the HBA1 gene, which is encoded by α-globin and is localized on chromosome 16p13.3 where the three SNVs are located. Second, thalassemia is prevalent in some regions of Asia and around the Mediterranean region but not in most European countries; these three SNVs are specific to the Asian population. Third, most TWB participants having the minor alleles of these three variants were both microcytic and hypochromic, with more than 90% of the participants having the microcytic hypochromic trait. Finally, more than 100 genetic forms of α-thalassemia have been identified in which α0-thalassemia is usually the most common clinically relevant form (Giardine et al, 2014). In Taiwan, the --SEA type of the α0-thalassemia mutation is the most common deletion mutation, accounting for 69–91% cases; this finding indicates that the prevalence of this mutation is ∼4.0–4.5% in Taiwanese individuals (Chen et al, 2002; Lee et al, 2015; Wang et al, 2017). All the three SNPs are specific to the Asian population and highly linked to the microcytic hypochromic trait, with a minor allele frequency of ∼1.77–1.94%. This translates to a heterozygous genotype frequency of ∼3.5–3.8% (Table 2), which is close to the predicted prevalence of α0-thalassemia mutation --SEA. The associations were further confirmed using both BreakDancer screening and PCR with direct sequencing, which revealed considerably high NPV (0.9937–0.9958) and PPV (0.9107–0.9808) values for the three SNVs. With the increasing health and economic burden of thalassemia in recent years because of population growth and epidemiologic transition, identifying surrogate genetic markers can be helpful for future mass screening and developing preventive medicine strategies for thalassemia.
The NPRL3 rs191086839 variant is the strongest and independent genetic surrogate marker for the α0-thalassemia deletion mutation --SEA in Taiwanese individuals
NPRL3 is a highly conserved gene located upstream of the HBA gene cluster. The α-globin locus of all mammalian species analyzed lies within a region of 135–155 kb of conserved synteny, with α-like genes arranged along the chromosome in the order 5′-ξ-α-α-3′ (Fig 3). The HBA cluster is located between NPRL3 and LUC7L genes in almost all mammals except mouse, in which the HBA cluster no longer has LUC7L downstream of the globin genes (Vernimmen, 2014; Philipsen & Hardison, 2018). The NPRL3 gene contains the enhancers of the HBA gene cluster. The erythroid-specific multispecies conserved sequences (MCSs) identified by DNase-hypersensitive sites are numbered from MCS-R1 to MCS-R4 (Hughes et al, 2005). Three of these elements (MCS-R1, MCS-R2, and MCS-R3) lie within the body, and MCS-R4 lies upstream of the promoter of the NPRL3 gene. MCS-R2 has multiple roles, and these roles may be applicable to any other enhancer: the recruitment of polymerase II and key transcription factors at the promoter, formation of a looped structure involving several remote regulatory elements, and removal of repressive complexes, such as PcG. By performing functional analysis, Miyata et al (2020) demonstrated that an erythroid-specific enhancer is located in the intron 7 of vertebrate NPRL3, which indicates the presence of a remote enhancer on nprl3 in multiglobin gene expression. In multivariate analysis, NPRL3 rs191086839 was observed to have the strongest effect on various erythrocyte indices, including MCV and MCH, and the microcytic hypochromic trait and microcytic anemia (Tables 3 and S2). Therefore, considering the robust evidence of co-inheritance, as reflected in the linkage disequilibrium analysis and the context of surrogate markers traditionally being associated with sensitivity and specificity, we propose NPRL3 rs191086839 as a strong and independent candidate genetic surrogate marker for the α0-thalassemia deletion mutation --SEA in Taiwanese individuals.
Association between α-thalassemia-related erythrocyte indices and metabolic traits
Our findings revealed that an elevated MCV was associated with lower risks of metabolic syndrome and some of its components and complications, such as DM, hypertension, microalbuminuria, and lower HDL cholesterol levels. These results are compatible with those reported previously indicating that an elevated MCV was associated with lower risks of metabolic syndrome and visceral obesity (Tanaka et al, 2020, 2021). Metabolic syndrome and its related components have been reported to be associated with various adverse cardiovascular and cancer outcomes (Lakka et al, 2002; Emery et al, 2022). In addition, a high MCV was determined to be associated with elevated total and LDL cholesterol, uric acid, and triglyceride levels in our study. Previous studies have reported that macrocytosis was associated with the severity or poor prognosis of various cardio-renal diseases (Solak et al, 2013; Ueda et al, 2013; Hsieh et al, 2017; Wang et al, 2020). Similar trends of associations were noted between MCH and most of the metabolic traits analyzed. These results suggest the diverse and bidirectional effects of MCV and MCH, which result in both favorable and unfavorable cardiometabolic outcomes. Additional studies may be necessary to elucidate the effect of the interaction between erythrocyte indices and various cardiometabolic risk factors on the outcome of cardiovascular diseases.
Causal relationships among α-thalassemia-related RBC traits, metabolic traits, and DM
Previous GWASs have revealed associations between chromosome 16p13.3 variants and HbA1c and LDL cholesterol levels (Sinnott-Armstrong et al, 2021; Hsu et al, 2022; Lee et al, 2022). Thus, we examine the association between PGAP6 rs375498857 and various metabolic and biochemical traits in 115,926 TWB participants with whole-genome genotyping array data. We observed that the rs375498857 genotype was closely related to many cardiometabolic traits, such as lipid profiles (total cholesterol, LDL cholesterol, and HDL cholesterol levels), HbA1c levels, and total bilirubin levels, and DM. By performing an MR study, we determined that the correlation between this SNV and cardiometabolic traits (low total and LDL cholesterol levels) is mainly mediated by erythrocyte traits, such as MCV and MCH. This study is the first to indicate the importance of α-thalassemia-related erythrocyte traits in causing cardiometabolic traits and increasing DM risk (OR = 1.01–1.04); however, exact mechanisms underlying these associations are not fully understood. Total and LDL cholesterol levels in β0-thalassemia carriers (thalassemia minors) were lower than those in age- and sex-matched controls (Maioli et al, 1989). The LDL-lowering effect of β0-thalassemia was proposed to be because of (1) reactive mild erythroid hyperplasia increasing LDL removal by the bone marrow and (2) the chronic activation of the monocyte–macrophage system resulting in the increased secretion of some cytokines that affect hepatic secretion and receptor-mediated removal of apolipoprotein B-containing lipoproteins (Deiana et al, 2000). The latter mechanism resembles mild chronic inflammation and has been linked to DM development through an increase in oxidative stress and insulin resistance. Moreover, our study confirms the findings of earlier reports indicating that having the α-thalassemia trait significantly increases the likelihood of developing gestational diabetes among Chinese women (Lao & Ho, 2001). Moreover, individuals with α-thalassemia in Iran had a higher risk of impaired glucose tolerance (either prediabetes or diabetes) than did the general population (Bahar et al, 2015). Therefore, health-care providers should screen individuals with α-thalassemia for these conditions and provide interventions to reduce their risk, including lifestyle modifications and medication. The causal relationship sheds light on underlying mechanisms that may be studied in the future, such as how variations in RBC production and function affect inflammation, oxidative stress, and endothelial dysfunction, all of which are linked to the emergence of cardiometabolic traits and DM. Future large transethnic prospective studies should be conducted to determine the roles of these SNPs in thalassemia screening and epidemiology and their association with cardiometabolic disorders.
Limitation
Because of the existence of ethnic heterogeneity, the results of this study should be interpreted with caution when applying them to different races. This includes investigating whether the relationship between SNVs and deletion mutations is the same among different races and whether other deletion mutations can be identified using matching SNVs.
In conclusion, our results revealed that three Asian-specific chromosome 16p13.3 variants, namely NPRL3 rs191086839, LUC7L rs372755452, and PGAP6 rs375498857, exhibit strong LD and can be used as surrogate genetic markers for α-thalassemia or the α0 thalassemia --SEA deletion, with high specificity and sensitivity. PGAP6 rs375498857 in GWAS arrays revealed a significant association with multiple metabolic traits. MR indicated that many α-thalassemia-related metabolic traits, such as total and LDL cholesterol levels, HbA1c, and DM, are causally related to erythrocyte traits, such as MCV and MCH. These results may be beneficial in the future preventive medicine for population screening and understanding the underlying genetic causes of thalassemia-related cardiometabolic traits and DM.
Materials and Methods
TWB participants
The TWB recruited adults from centers across Taiwan between 2008 and 2020. The initial cohort for the regional association study of erythrocyte indices consisted of 1,494 TWB participants who underwent WGS; one individual from this cohort was excluded because of incomplete WGS data. Subsequently, to expand the study cohort, we included 129,542 participants who underwent genotyping with Axiom Genome-Wide CHB Arrays for a regional association study. From this cohort, 13,616 participants were excluded for the following reasons: the data of 7,216 participants were excluded as quality control (QC) for the GWAS, 4,647 participants fasted for <6 h, and 1,753 participants had failed genotyping for the rs375498857 polymorphism. During QC for the GWAS, all excluded participants were related as second-degree relatives or closer, with an identity-by-descent value of >0.187. Participants with a history of hypertension, hyperlipidemia, DM, and gout were excluded when respective parameters were analyzed. Hypertension, DM, obesity, current smoking, microalbuminuria, and metabolic syndrome were defined as reported previously (Wu et al, 2022). This study was approved by the Research Ethics Committee of Taipei Tzu Chi Hospital, Buddhist Tzu Chi Medical Foundation (approval numbers: 05-X04-007 and 10-XD-056) and the Ethics and Governance Council of the TWB (approval numbers: TWBR10507-02 and TWBR10611-03). Each participant was asked to sign an approved informed consent form.
DNA extraction and genotyping
Genomic DNA from blood samples was isolated using the PerkinElmer Chemagic 360 instrument (PerkinElmer). SNV genotyping was performed either through WGS or by using custom TWB chips. SNV genotyping was performed using the Axiom Genome-Wide Array Plate System (Affymetrix).
Hematological indices
The following hematological parameters were analyzed: RBC count and white blood cell count, platelet count, and hematocrit (Hct) and hemoglobin (Hb) concentrations. Briefly, blood cell indices were calculated using a hematology analyzer, which provided complete count data. Other parameters were calculated from the same indices, as follows: MCV, which was calculated as Hct divided by the RBC count, and MCH, which was calculated as the Hb concentration divided by the RBC count. The microcytic hypochromic trait was considered when MCV was <80 fl and MCH was <25 pg. Microcytic anemia was defined as the microcytic hypochromic trait with Hct <40% for men and <36% for women or Hb < 13 g/dl for men and <12 g/dl for women.
Regional association studies performed using TWB WGS data for gene regions surrounding the HBA1 gene
To determine the association of the significant SNVs located on chromosome 16p13.3 with hematological indices, we performed a regional association analysis of the WGS data of TWB participants (Fig 1). WGS data from a subgroup of TWB participants were evaluated through an ultrafast whole-genome secondary analysis on Illumina sequencing platforms with Illumina HiSeq 2500 or Illumina Hiseq 4000 sequencers (Raczy et al, 2013). The resulting reads were aligned to the hg19 reference genome with iSAAC 01.13.10.21. iSAAC Variant Caller 2.0.17 was used to perform SNP and insertion–deletion variant discovery and genotyping (Raczy et al, 2013). A total of 1,493 participants with 18,198 SNVs located at positions between 0.06 and 0.68 Mb on chromosome 16p13.3 were recruited. An in-house protocol written in shell script was used to combine 1,493 vcf files. A union table of all detected variants was used for regional association analysis. The association between SNPs and clinical and laboratory parameters was analyzed using the GWAS method.
Detection of deletion mutation by using deletion mutation detection software BreakDancer
The bam files of WGS data were provided by TWB (Juang et al, 2021). The reference genome version used was hg19. Samtools v1.3.1 (Li et al, 2009) was used to select the region of interest from original WGS bam files. BreakDancer v1.3.6 (Chen et al, 2009) was used to detect long deletion mutations in the region of interest.
Polymerase chain reaction amplification and direct sequencing
To detect candidate α-thalassemia deletion mutations, the WT and α0 thalassemia deletion --SEA sequences were genotyped through PCR by using the primer pairs (forward) 5′-GCGATCTGGGCTCTGTGTTCT-3′, (reverse) 5′-GTTCCCTGAGCCCCGACACG-3′ (Chang et al, 1991), and (reverse) 5′-GCCTTGAACTCCTGGACTTAA-3′, respectively (Sanguansermsri et al, 1999). The α0 thalassemia deletion --THAI sequence was genotyped through PCR amplification by using the primer pairs (forward) 5′-CACGAGTAAAACATCAAGTACACTCCAGCC-3′ and (reverse) 5′-TGGATCTGCACCTCTGGGTAGGTTCTGTACC-3′. The α0 thalassemia deletion --FIL sequence was genotyped through PCR amplification by using the primer pairs (forward) 5′-AAGAGAATAAACCACCCAATTTTTAAATGGGCA-3′ and (reverse) 5′-GAGATAATAACCTTTATCTGCCACATGTAGCAA-3′ (Liu et al, 2000).
Clinical phenotypes and laboratory examinations
To determine cardiometabolic traits, we examined the following clinical phenotypes: age; BMI; waist circumference; waist–hip ratio; and systolic, mean, and diastolic blood pressure. We also collected the following biochemistry data: lipid profiles, namely total, HDL, and LDL cholesterol and triglyceride levels; glucose metabolism parameters, namely fasting plasma glucose and HbA1c levels; and liver and renal functional test–related parameters, namely serum creatinine, uric acid, aspartate aminotransferase, alanine aminotransferase, γ-glutamyl transferase, albumin, and total bilirubin levels. BMI and the estimated glomerular filtration rate were calculated as reported previously (Hsu et al, 2021). Because of the absence of data on the urine creatinine level, only the spot urine albumin level was used to evaluate albuminuria. The hematological parameters analyzed also included white blood cell and platelet counts.
Regional association analysis of cardiometabolic traits
We performed regional association studies including the Axiom Genome-Wide CHB array data of participants, as reported previously (Hsu et al, 2022; Yeh et al, 2022). A total of 115,926 participants were included in the analysis.
Statistical analysis
Categorical data are presented as frequencies and numbers. Continuous variables are expressed as medians and interquartile ranges. Lipid parameters and urine albumin levels were logarithmically transformed when analyzed using regression. A general linear regression was performed in the association study between the phenotypes and genotypes, with adjustment for potential confounders, such as age, sex, BMI, and current smoking status. Logistic regression was performed to investigate the effect of variants on microcytic traits, anemia, and other risk factors. Stepwise multivariable regression was performed to determine the independent correlates of hematological parameters. Statistical analyses were performed using IBM SPSS Statistics for Windows (Version 22; SPSS). We used PLINK (version 1.07; Shaun Purcell, https://zzz.bwh.harvard.edu/plink/, accessed on 3 January, 2022) for regional association analysis. LocusZoom (http://csg.sph.umich.edu/locuszoom/, accessed on 4 January, 2022) was employed for drawing regional plots. A P-value of < 5 × 10−8 was defined as indicating genome-wide significance. For P-values between 5 × 10−8 and 1 × 10−4, a subthreshold that suggested GWAS significant association of SNVs was considered. The subthreshold locus has been demonstrated to be effective for predicting studied phenotypes confirmed using luciferase reporter assays (Wang et al, 2016). Linkage disequilibrium (LD) between SNVs was analyzed using LDmatrix software (https://ldlink.nih.gov/?tab=ldmatrix, accessed on 4 January, 2022).
MR approaches
We performed an instrumental variable (IV) regression analysis by using 2SLS methods to examine whether rs375498857, which is a α-thalassemia-related MCV and MCH-determining genetic variant, is associated with other cardiometabolic traits and DM development through their associations with MCV or MCH. The 2SLS method is widely used for evaluating continuous and binary exposures and outcomes (Bowden et al, 2016; Burgess et al, 2017; Hsu et al, 2022). In standard MR, the first stage involved regressing MCV- and MCH-determining SNV to generate predicted MCV and MCH values. The second stage involved regressing study parameters on rs375498857 to predict the risk of cardiometabolic traits. The strength of the instruments was assessed using the F statistic that was calculated using the equation F = R2(n − 2)/(1 − R2), where R2 is the proportion of the variability in genetically determined MCV or MCH accounted for by the SNV, and n is the sample size (Palmer et al, 2012). An F statistic of >10 indicates a relatively low risk of weak instrument bias in MR analyses (Palmer et al, 2012). Statistical power for 2SLS MR was calculated using the noncentrality parameter-based approach described previously (Brion et al, 2013).
Data Availability
The data presented in this study are available on request from the corresponding author.
Ethics statement
This study was approved by the Research Ethics Committee of Taipei Tzu Chi Hospital, Buddhist Tzu Chi Medical Foundation (approval numbers: 05-X04-007 and 10-XD-056) and the Ethics and Governance Council of the TWB (approval numbers: TWBR10507-02 and TWBR10611-03). Each participant was asked to sign an approved informed consent form.
Acknowledgements
We greatly appreciate the technical support of the Core Laboratory of the Taipei Tzu Chi Hospital, Buddhist Tzu Chi Medical Foundation, expert statistical analysis assistance from Tsung-Han Hsieh, and Kimforest LTD Taiwan for bioinformatics support. This study was supported by grants from Buddhist Tzu Chi Medical Foundation (TCMF-EP 111-02), grants from the Ministry of Science and Technology (MOST 108-2314-B-303 -026 -MY3), and Taipei Tzu Chi Hospital, Buddhist Tzu Chi Medical Foundation (TCRD-TPE-112-02) to Y-L Ko.
Author Contributions
L-A Hsu: formal analysis, validation, methodology, and writing—original draft, review, and editing.
S Wu: software, validation, methodology, and writing—review and editing.
M-S Teng: methodology.
Y-L Ko: conceptualization, resources, formal analysis, supervision, funding acquisition, validation, project administration, and writing—original draft, review, and editing.
Conflict of Interest Statement
The authors declare that they have no conflict of interest.
- Received June 7, 2023.
- Revision received September 15, 2023.
- Accepted September 19, 2023.
- © 2023 Hsu et al.
This article is available under a Creative Commons License (Attribution 4.0 International, as described at https://creativecommons.org/licenses/by/4.0/).