Abstract

Main

The architecture of complex traits in humans can be fully understood only in the context of history. Present-day Mexico covers seven cultural regions, including much of Mesoamerica, with rich civilizational histories⁹. Archaeological and anthropological approaches have been used to regionalize Mexico into the north of Mexico, the north of Mesoamerica, the centre, occident and Gulf of Mexico, Oaxaca (referring here to the Oaxaca cultural region) and the Mayan region¹⁰ (Fig. ​(Fig.1a).1a). These regions are based on specific Indigenous civilizations and cultures, which flourished early in the Mayan region, Oaxaca, and the occident and the Gulf of Mexico, and later in the centre and north of Mesoamerica. Such histories have also been used to classify Mesoamerican chronology into preclassical, classical, postclassical, colonial and postcolonial periods¹¹.

Open in a separate window

Fig. 1

Mosaic ancestral patterns in the MXB and the genetic diversity within Mexico.

a, Sampling for the MXB (n = 5,812 individuals with latitude and longitude values), showing Mexico regionalized into Mesoamerican regions according to an anthropological and archaeological context. b, Unsupervised clustering using ADMIXTURE and global reference panels (n = 9,007 including MXB) from the 1000 Genomes Project, the Human Genome Diversity Project and the Population Architecture using Genomics and Epidemiology Study. c, Uniform manifold approximation and projection (UMAP) analysis of MXB (n = 5,622) coloured by Mesoamerican region. d, Archetypal analysis of MXB (n = 5,833) with reference global individuals as in b, coloured by region (top) or in grey (bottom). This approach determines each individual’s position in a ten-dimensional space that in this visualization is reduced to two dimensions. Reference individuals (bottom) are coloured using ADMIXTURE inferred clusters from b. For example, for the Americas (1000 Genomes) and Middle East, where multiple clusters are inferred, a colour combining these cluster colours is used.

In the past 500 years, Spanish colonization has left an indelible mark on this Indigenous tapestry. In a colonial and postcolonial context, genetic ancestries that trace principally to Western Europe, West Africa and East Asia can be identified in present-day Mexicans^12–16. These genetic ancestries vary in structure and timing between Mesoamerican regions and give rise to extensive fine-scale population substructure and ancestry sources across Mexico^12–16. Further, such varying genetic histories, as captured by ancestry distributions, have been shown to affect variation in complex traits such as lung force capacity¹², and a number of other complex traits and diseases¹⁷.

Nevertheless, a large gap remains in the representation of Mexicans from across Mexico in cohorts with linked genotypes and phenotypes. Such representation could enable finer-scale studies of genetic histories and a better understanding of complex trait architectures among individuals with diverse ancestries from the Americas and those living in rural areas¹⁸. Past analyses on complex traits have been limited to studying individuals from the USA and Mexico City^12,17. They have also not simultaneously modelled the influence on complex trait variation of a rich array of genetic and environmental factors as is possible with a nationwide biobank.

To bridge this gap, we launched the Mexican Biobank (MXB) project, densely genotyping 6,057 individuals from 898 localities distributed nationwide (Supplementary Figs. 1 and 2) recruited by the National Institute of Public Health (Instituto Nacional de Salud Pública) across all 32 states of Mexico. To select the samples for genomic and biochemical characterization, we enriched for those individuals that speak an Indigenous language while maximizing geographic coverage and the inclusion of rural localities (about 70% of the MXB; Supplementary Figs. 2–5). Of the participants in the MXB, 70% are female, and it comprises data for individuals born between 1910 and 1980 (Supplementary Table 1) who were genotyped at about 1.8million single nucleotide polymorphisms (SNPs) and have linked information for complex traits, sociocultural and biogeographical markers (Supplementary Table 2).

Here, we leverage rich archaeological and anthropological information to guide a regionalized analysis of Mexico, and harness the power of genome-wide local ancestry estimation and identity-by-descent (IBD) segments to decipher fine-scale genetic histories using ancestry-specific approaches to denote origins and historical population size changes^4,19. We reveal a very heterogeneous landscape of both, painting a genetically informed picture of varying demographic trajectories in Mesoamerican regions, including colonial migrations and dynamics. We further investigate the role of these evolutionary histories as captured by proxies of genetic ancestries in shaping genetic variation and complex trait patterns in Mexico today. We show that these histories result in marked geographic and ancestry-specific patterns in the distributions of runs of homozygosity (ROH) and of the genomic burden of rare deleterious variants. We carry out GWAS analyses across 22 binary and quantitative traits, and compare the prediction performance of polygenic scores computed using our GWAS or UK Biobank (UKB) GWAS data. Last, given that evolutionary histories (captured by genetic ancestries) could associate specific trait-relevant genotypes with certain genetic backgrounds, we study the impact of genetic ancestries, portions of the genome in ROH, polygenic scores and other sociocultural and biogeographical factors on creating variation in complex and medically relevant traits in Mexico.

Diverse ancestries across timescales

We begin by analysing the population structure in the MXB at different geographic resolutions and timescales (see the section entitled ‘Note on genetic ancestries’ in the Methods; Fig. ​Fig.11 and Supplementary Figs. 6–24). Given the history of Mexico, in which genetic lineages are expected to trace back to disparate geographic regions (for example, the Americas, Western Europe, West Africa and East Asia) in the past approximately 500years, we first analyse each individual in a framework that infers proportions of genetic ancestries on the basis of genetic similarity to other individuals (using ADMIXTURE²⁰) in a global reference sample. We use a similar approach to label local segments across the genomes of the study individuals. We use the term ‘ancestries from the Americas’ when referring to genetic ancestries that derive from genetic ancestors living in the Americas before European colonization; these have also been referred to as Indigenous ancestries, and in some places below we also use this term (Fig. ​(Fig.1b,1b, Supplementary Figs. 11 and 12 and Supplementary Table 3).

Higher proportional ancestries from the Americas are inferred in Mexico’s central and southern states, compared to the northern states, and ancestries from West Africa are observed in every state²¹ (Supplementary Table 3) in agreement with historical records of shipping voyages from the transatlantic slave trade^21,22 (Supplementary Fig. 14). We note the presence of a small but substantial proportion of ancestries from East Asia in almost every state (0–2.3%), the highest in the state of Guerrero (2.3%), and an even more modest proportion of ancestries from South Asia in most states as well (0–0.8%). These probably reflect migrations from Asia to Mexico dating to the Manila Galleon trade in the sixteenth and seventeenth centuries^16,23–26, and later nineteenth- and twentieth-century migrations from China and Japan, especially to the north of Mexico^27–29.

We observe the most significant genetic differentiation along a north-to-southeast cline in Mexico (measured using F_ST, which is an index quantifying the proportion of the total genetic variance contained in subpopulations (S) relative to the total genetic variance (T); Supplementary Figs. 15–18). When considering autosomes of only individuals with ≥90% proportion of ancestries from the Americas (inferred using ADMIXTURE), the Mayan region of Chiapas, Tabasco, Yucatan, Quintana Roo and Campeche show relatively larger F_ST values with the other regions (Supplementary Figs. 17 and 18). This distinction is also apparent using ADMIXTURE-inferred ancestral clusters (Supplementary Fig. 13) and dimensionality reduction techniques highlighting this population substructure within Mexico (Fig. 1c,d and Supplementary Figs. 7–20). Individuals from the Mayan region tend to cluster mostly together, but overlap with individuals from the Gulf of Mexico and central Mexico, consistent with oral histories. In the rest of the regions, subtle substructure mirroring Mesoamerican geography is visible in the MXB, probably reflecting both unique local demographic histories of Indigenous ancestries and the effects of movement and mating among the different regions. Compared to previous sampling and analyses that focused on Indigenous groups with varying degrees of isolation in Mexico¹², the MXB reveals lower average levels of F_ST and substructure, probably owing to the broader sampling (although the substructure presented by the Mayan region is more apparent in the MXB). The method of ref.⁵ further highlights the ancestral diversity reflected by the MXB samples that are represented as mixtures of multiple sources (Supplementary Fig. 22) in the presence of global references (Fig. ​(Fig.1d1d and Supplementary Figs. 21–23). Individuals from the same region (for example, the Mayan region) are modelled as mixtures of several sources, reflecting the diversity of ancestry variation within this and other Mesoamerican regions. Given this variation among ancestries from the Americas and the unique power given by the MXB to explore its impact on complex trait variation, we also obtain an axis of variation within ancestries from the Americas (Supplementary Fig. 24 and Supplementary Table 4).

Genetic histories inferred within Mexico

Contemporary Mexicans derive ancestries predominantly from diverse lineages found in the Americas, Western Europe and West Africa. These ancestral sources have different demographic histories before their arrival in present-day Mexico and probably after their arrival within different Mesoamerican regions. To reveal the history of effective population sizes (N_e) of these three ancestries in the MXB, we analyse IBD segments^4,30 stratified by local ancestry inference for each Mesoamerican region⁴ (Fig. ​(Fig.22).

Open in a separate window

Fig. 2

Effective population size (N_e) values across ancestries and geographies reveal the histories present within Mexico.

a, Mesoamerican chronology colouring different periods in Mesoamerican history using an anthropological and archaeological context. b, Ancestry-specific effective population size (N_e) changes over the past 200 generations across Mexico (n = 5,436) inferred using IBD tracts, coloured by chronology from a assuming 30 years per generation. c, Ancestry-specific effective population size (N_e) changes over time for ancestries from the Americas in different regions of Mexico (see Supplementary Figs. 25–29 for other generation intervals and ancestries). n = 1,177, 640, 952, 590, 820, 315 and 938 for the north of Mexico, north of Mesoamerica, centre of Mexico, Gulf of Mexico, occident of Mexico, Oaxaca region and the Mayan region, respectively.

We observe fine-scale variation in N_e trajectories for Indigenous lineages which we interpret in the context of the different cultural histories of Mesoamerican regions⁹ (Fig. ​(Fig.2).2). As generational time can vary, we present our analysis at two extremes of 20 and 30 years per generation³¹ (Supplementary Figs. 25 and ​and2c,2c, respectively). Chronologically speaking, archaeologists document that Mesoamerican civilizations flourished first in the Mayan region, in Oaxaca, in the occident and in the Gulf of Mexico. In these regions, we observe large N_e already in the classical period (250–900ce)³². For example, in the Gulf, where we observe high N_e since the preclassical period (2500bce–250ce), there is archaeological evidence, among a myriad of other groups, of the Olmecs in the preclassical period, the Totonacs in the classical period and the Huastecs in the postclassical period (900–1521ce)³³. In Oaxaca, we observe N_e rapidly growing in the preclassical to the classical period, in line with archaeological inferences that the Zapotecs were already starting to create sedentary settlements in the preclassical period followed by a rise in social and political structures in the classical period. The subsequent postclassical period was characterized by militarism and warfare³⁴, and our genetic evidence suggests a population decline towards the end of the postclassical period. In the Yucatan peninsula, the Maya had a prominent civilizational spread in the classical period (peak N_e observed). They started going through a slow decline only in the postclassical period due to what archaeologists have inferred as a combination of different political and ecological factors, and this trajectory is supported in the N_e trend³².

These patterns contrast with those of the centre and north of Mesoamerica, where the Aztec empire had a stronghold most recently; there we see increasing N_e in the postclassical right before the arrival of the Spaniards and into part of the colonial period, after which we start to see a population decline in N_e. The decrease in N_e after the arrival of the Spaniards is most prominent in the centre and north of Mesoamerica. In Oaxaca and the Mayan region, where Indigenous ancestries from the Americas are most prevalent today as evidenced by the ADMIXTURE analysis (Supplementary Table 3), the decrease in N_e is followed by an increase in the postcolonial period.

Concurrently, we observe that ancestries from Western Europe that entered the contemporary Mexican gene pool went through a sharp decline in effective population size during the colonial period. The extent of the founder effect varied by region, with the strongest effect seen in Oaxaca and the Mayan region (Supplementary Figs. 26 and 27). Ancestries from West Africa in Mexico revealed stronger founder effects that varied by region, with N_e ranging between 10³ and 10⁴ in the colonial period. The population size in the postcolonial period continued to grow in some regions such as the occident and north of Mexico and the Mayan region, compared to others (Supplementary Figs. 28 and 29). Consistent with previous results on self-identified Indigenous groups^13,14, our results on the MXB individuals highlight the heterogeneity of group histories across the Mesoamerican regions as well as the expansion of Indigenous lineages in the postcolonial period in several regions.

We further generated ‘admixture graphs’⁶ for individuals from the Mesoamerican regions to investigate their shared history by using an ancestry-specific approach and limiting the analysis to genomic segments with ancestries from the Americas. The admixture graph approach models the different Mesoamerican regions as populations in a progression of splits (Extended Data Fig. ​Fig.1a1a and Supplementary Fig. 30), providing information about the genetic relationships among the different regions. We can observe a clear progression of splits among the populations from north to south, with the north of Mexico splitting first, followed by the common ancestor of the north of Mesoamerica and the occident of Mexico, followed by the common ancestor of the remaining regions. Notably, the centre of Mexico and the Mayan region are related, consistent with previous suggestions based on IBD¹² and our population structure results, and both share a common ancestral source with Oaxaca and the Gulf of Mexico. These results further strengthen evidence for an Atlantic coastal corridor of gene flow between the Yucatan peninsula and central Mexico and the Gulf of Mexico previously posited in ref.¹². As demographic histories can affect patterns of genetic variation, such as distributions of ROH and of the genomic burden of deleterious variants, we next evaluate these metrics.

Open in a separate window

Extended Data Fig. 1

Genetic histories and polygenic prediction in the MXB.

A) Admixture histories of individuals in different cultural regions using an AdmixtureBayes approach. Here the admixture graph with the highest posterior probability is shown, inferred using genomic regions with ancestries from the Americas. Internal inferred ancestral node populations are colored grey. The tree is rooted using the Han as an outgroup. B) Trait variance explained and p-value threshold of best predictive polygenic score using MXB-GWAS-based or UKB-GWAS-based prediction. Polygenic scores were computed using SNPs significant at five different p-value thresholds (0.1, 0.01, 0.001, 0.00001, 10⁻⁸). A linear null model was created for each trait including age, sex and 10 principal components as covariates. A second polygenic score model was created adding the polygenic score to the null model. We computed the R² of the polygenic score by taking the difference between the R² of the polygenic score model and the R² of the null model. The maximum R² was used to the pick the p-value threshold for the best predictive polygenic score shown in the table.

Impact of genetic histories on variation

We analyse the patterns of ROH in the MXB including how they vary across geography and genetic ancestry proxies (inferred from ADMIXTURE). ROH patterns help further illuminate the demographic and mating histories of Mexicans³⁵, and are especially relevant for variation in complex traits when trait-relevant variation is affected by partially recessively acting alleles³⁶. We identify ROH (≥1Mb) in the MXB and observe that both the number of ROHs and the total length of ROH per individual increase as we move from north to southeast in the country (Supplementary Fig. 31). We confirm that this is primarily due to individuals with a higher inferred proportion of genetic ancestries from the Americas also having more ROH, particularly small ROH (smaller than those expected from recent consanguinity; for example, <8Mb), in their genomes (Fig. ​(Fig.3a,3a, Supplementary Figs. 32 and 33 and Supplementary Table 5). The appearance of many small ROHs indicates coalescences occurring at a period in the more distant past; for example, due to an ancient bottleneck or relatively small historical population size³⁷.

Open in a separate window

Fig. 3

Demographic histories affect patterns of genetic variation in Mexico.

a, Small ROH prevalence is correlated with ancestry proxies inferred from ADMIXTURE reflecting an ancient bottleneck or relatively small population size in the past (n = 5,833 individuals). b, Sum of ROH per individual as a function of birth year (n = 5,833 individuals). Solid lines show ROH overall, and dashed lines indicate ROH overlapping ancestries from the Americas (AMR). ROH are divided into small, medium and large ROH, as in a. Smoothed conditional mean lines are shown using the locally estimated scatterplot smoothing method. Error bands represent 95% confidence intervals. c, Mutation burden in different ancestries shows the effects of bottleneck events in causing loss of rare variants (n = 5,818 individuals). Rare variants are correlated with levels of ancestries from the Americas, Western Europe or West Africa for rare variants (derived allele frequency≤5%). Smoothed conditional mean lines are shown using a linear model. Error bands represent 95% confidence intervals. Spearman correlation values are shown (R and two-sided Pvalues) for all ancestries. Analysis of whole-genome sequences from 1000 Genomes MXL samples shows that the rare mutation burden result is robust to ascertainment bias of Illumina’s Multi-Ethnic Global Array (Supplementary Figs. 39 and 40). Variants were annotated using the Variant Effect Predictor tool, and nonsynonymous (deleterious) variants are a combined set of missense variants predicted to be damaging by polyphen2 along with splice, stop lost and stop gained variants.

Further, we observe that ROH found on Indigenous genomic segments are more frequent in younger individuals compared to older individuals (Spearman’s ρ=0.31, P=0.016; Fig. ​Fig.3b).3b). We corroborate that this correlation with birth year primarily derives from small ROH (ρ=0.35, P=0.006), and small ROH found on Indigenous genomic segments (ρ=0.39, P=0.002; Fig. ​Fig.3b).3b). The result is at least partly due to younger individuals having higher proportions of Indigenous ancestries compared to older individuals, especially in the rural localities (Supplementary Figs. 34 and 35), and agrees with recent observations about ancestry and ROH made in MexicanAmericans¹⁷. We also confirm that this observation is not due to sampling bias (see the sampling ascertainment note in the Supplementary Information). The observation of higher ancestries from the Americas in younger individuals in rural areas may be due to higher fertility rates in rural areas or individuals with other ancestries moving out from rural to urban areas.

We also investigate the effects of demographic histories on the frequency distribution of genetic variants. This analysis is motivated by previous theoretical and empirical work showing that undergoing a bottleneck changes the allele frequency distribution in the group that experienced the bottleneck^38–40, while leaving the overall sum of deleterious alleles per individual (‘deleterious mutational burden’) unchanged^39,41,42. In particular, rare variants are lost or increase in frequency after the bottleneck.

We evaluate this effect by computing the genome-wide sum of intergenic, synonymous and putatively deleterious (predicted-damaging missense and loss of function) alleles per individual. When considering only rare alleles (derived allele frequency≤5%), we observe that individuals with higher ancestry proportions from the Americas carry fewer rare derived alleles across variant types (strongest effect observed for intergenic variants) (Fig. ​(Fig.3c)3c) in contrast to other ancestries. We verified these observations with whole-genome sequences from a 1000 Genomes Project cohort (Mexican Ancestry in Los Angeles, California or MXL) (Supplementary Fig. 39), as well as with 50 genomes sequenced as part of the MXB project (Supplementary Fig. 40), to rule out ascertainment biases due to the array genotyping. Our result probably reflects primarily founder events during the peopling of America or subsequent genetic drift leading to loss of rare variants and/or their rise to higher frequencies.

GWAS and polygenic prediction in the MXB

To understand trait-associated locus transferability, we conduct GWAS analyses across 22 binary and quantitative traits (Supplementary Table 6). We identify genome-wide significant loci passing Bonferroni correction (P<2.27×10⁻⁹) on chromosomes 1, 9, 11 and 16 associated with lipid levels in blood (Fig. ​(Fig.4a).4a). Fine-mapping of independent signals within these loci reveals variants in or near CELSR2 (low-density lipoprotein (LDL): rs7528419), ABCA1 (high-density lipoprotein (HDL): rs9282541 and rs2065412), the LINC02702–BUD13–ZPR1–APOA1–APOA4–APOA5–APOC3–SIK3 locus (HDL: rs180326 and rs200905431; LDL: rs66505542; triglycerides: rs947989, rs66505542 and rs5104), HERPUD1–CETP (HDL: rs57502215, rs56129100, rs193695, rs56228609 and rs117427818; cholesterol: rs57502215, rs56228609 and rs118146573) and APOE (LDL: rs7412; triglycerides: rs440446), which have all previously been associated with lipid levels in European and Hispanic groups (Supplementary Table 7). Notably, we replicate the association of the ABCA1*C230 allele that has previously been associated with decreased HDL cholesterol levels (β=−0.219, s.e.=0.030, P=1.64×10⁻¹³; Fig. ​Fig.4a),4a), and is found almost exclusively in Indigenous groups from the Americas⁴³. This association was replicated in the subset containing >90% inferred Indigenous ancestries although it did not reach genome-wide significance (β=−0.210, s.e.=0.055, P=1.22×10⁻⁴). Restricting the GWAS cohort to individuals with >90% inferred ancestries from the Americas did not identify any genome-wide significant loci.

Open in a separate window

Fig. 4

Illustrative examples of GWAS and polygenic prediction in the MXB.

a, Manhattan plots showing GWAS results for HDL cholesterol (top, n = 4,484) and triglycerides (bottom, n = 4,483) in the full MXB dataset. Fine-mapped genes are labelled (Methods). To aid with visualization, 1 in 200 SNPs with P>0.01 were sampled for the Manhattan plots. b, Prediction performance is measured by the correlation between polygenic score (the sum of all alleles associated at P<0.1 weighted by their estimated effect sizes) and trait value (as measured by Pearson correlation R and its associated two-sided P value) for HDL cholesterol (top, n = 1,327) and triglycerides (bottom, n = 1,326). According to the schematic in Supplementary Fig. 41, for b, GWAS was carried out in two-thirds of the MXB, and the remaining one-third of the MXB was used to compute polygenic scores and test their ability to predict complex traits. Smoothed conditional mean lines are shown using a linear model. Error bands represent 95% confidence intervals. Scores were computed using TOPMed-imputed MXB genotypes. Traits were normalized using an inverse normal transform (INT) for both a and b. For further evaluation of prediction performance, see Extended Data Figs. ​Figs.1b1b and ​and22–10 and Supplementary Tables 8 and 9.

To assess transferability in the prediction of quantitative traits using polygenic scores, we re-perform a GWAS in only 4,000 randomly selected individuals from the MXB and construct polygenic scores in the remaining 1,778 individuals (Supplementary Fig. 41). We compute polygenic scores using both genotype data and imputed genotypes using TOPMed. To assess the impact of using different GWAS summary statistics on prediction performance, for comparison we also compute polygenic scores using pan-ancestry GWAS from the UKB in light of the varying ancestry sources in Mexico (Fig. ​(Fig.4b4b and Extended Data Fig. ​Fig.1b).1b). We observe that MXB-based prediction works better or as well as UKB-based prediction, despite much lower sample size, for glucose, creatinine, cholesterol and diastolic blood pressure (Extended Data Figs. ​Figs.1b1b and ​and22–10 and Supplementary Tables 8 and 9). Triglycerides, HDL and LDL cholesterol levels are also almost as well predicted by the MXB GWAS (Fig. 4b). These results indicate that further gains in prediction power would be achieved by increasing the sample size further. Although many factors are probably involved in differential polygenic score portability by trait, some trait architecture features are probably relevant, such as the strength of stabilizing selection that the trait is under, its mutational target size and heritability per causal site^44,45. Using estimated mutational target sizes from previous GWAS studies⁴⁵, we observe that traits with smaller inferred mutational target sizes (creatinine and triglycerides) are predicted better with SNPs discovered in MXB compared to traits inferred to have larger target size (height and body mass index (BMI))⁴⁵. UKB-based predictors are used in our complex trait modelling below, as these can be computed for all MXB individuals.

Open in a separate window

Extended Data Fig. 2

Prediction performance of MXB-GWAS-based or UKB-GWAS-based polygenic scores computed for height in the MXB.

Traits are inverse normalized. Prediction performance is measured by the correlation between polygenic score (the sum of all alleles associated at p<0.1 weighted by their estimated effect sizes) and trait value (Pearson correlation R and associated two-sided p-value). Smoothed conditional mean lines are shown using a linear model. Error bands represent 95% confidence intervals.

Open in a separate window

Extended Data Fig. 10

Prediction performance of MXB-GWAS-based or UKB-GWAS-based polygenic scores computed for diastolic blood pressure in the MXB.

Traits are inverse normalized. Prediction performance is measured by the correlation between polygenic score (the sum of all alleles associated at p<0.1 weighted by their estimated effect sizes) and trait value (Pearson correlation R and associated two-sided p-value). Smoothed conditional mean lines are shown using a linear model. Error bands represent 95% confidence intervals.

Complex trait architectures in the MXB

Last, we assess the contribution of genetic variation resulting from variable demographic and environmental histories or causal variant distributions towards affecting variation in complex traits or diseases in Mexico (Supplementary Fig. 42). We focus on several quantitative traits: height, BMI, triglycerides, cholesterol, glucose, blood pressure and others. Aiming to understand how the traits are distributed geographically and relative to single model covariates, we first visualize average trait values by units of our biogeographical and sociocultural factors to understand the dimensions of trait variation (Fig. ​(Fig.5a5a and Supplementary Figs. 43–51).

Open in a separate window

Fig. 5

An analysis of the factors influencing height and other complex trait variation.

a, Bottom: map of average height in Mexico (n = 5,770). Height was normalized using an INT. Top: box plots of height (INT) variation in each state from northwest to southeast. The box plots show the median value and the quartiles. Whiskers extend to the minimum and the maximum values. The dots represent outliers. n=5,846 biologically independent samples were used for the analysis. b, Explanatory model for height variation implicates the role of genetics and environment. The plot shows effect-size estimates and confidence intervals (1.96×s.e.m.) from a mixed-model analysis. All quantitative predictors are centred and scaled by 2 standard deviations. Asterisks show significance at false discovery rate<0.05 across traits and predictors analysed⁵⁰. n=4,625 biologically independent samples were used for the analysis. c, Height as a function of birth year in quartiles of ancestries from the Americas (n = 5,598). Error bands represent 95% confidence intervals. d, Trait profiles for BMI (left), triglycerides (middle) and glucose (right). Results of mixed-model analysis, as in b. The plot shows effect-size estimates and confidence intervals (1.96×s.e.m.) from a mixed-model analysis. n=4,607, 3,664 and 3,613 biologically independent samples were used for the analysis for BMI, triglycerides and glucose, respectively. For b and d, PS are polygenic scores computed using UKB summary statistics (SNPs significant at P < 10^–8), A(Africa/East Asia/Americas) refers to ancestry proportions from that region as inferred from ADMIXTURE, and MDS1(A(Americas)) and MDS2(A(Americas)) refers to multidimensional scaling (MDS) axes within ancestries from the Americas as inferred using a MAAS-MDS analysis (Supplementary Fig. 24). Educational (Edu.) attainment is on a scale from 0 to 8 (low to high educational attainment), and altitude is measured in metres (low to high).

Next we use a mixed model to estimate the contribution of genetic factors to trait variation jointly modelled with the environmental factors (Fig. 5b,d and Supplementary Figs. 52–61). Genetic ancestry proxies can be associated with complex traits due to genetic factors or due to non-genetic factors that covary with genetic ancestries such as differential experiences of discrimination, dietary nutrition and socioeconomic status (Supplementary Fig. 42). The genetic factors that vary with genetic ancestry proxies can be different distributions of ROH or other differential patterns of genetic variation caused by demographic and environmental histories that vary among ancestries. ROH have also previously been shown to have associations with a broad range of complex traits such as height, weight and cholesterol, pointing towards a recessive architecture of these traits^36,46. As shown above, genetic ancestry proxies in the MXB are correlated with the number and length of ROH (Fig. ​(Fig.3a).3a). We, therefore, develop a mixed model for the association of genetic factors such as ancestry proxies, ROH and polygenic scores with trait variation. We consider in our model several environmental factors to improve power and to query the role of genetic factors reflected in ancestry proxies compared to environmental factors. We include variables available in the MXB related to discrimination, socioeconomic opportunities and living environment (collectively called sociocultural and biogeographical factors), as well as unobserved random effects to model cryptic relatedness and potential unmodelled environmental factors. In this model, a significant association with ancestry proxies could reflect the association of particular causal genotypes with those ancestries or associated unmodelled environmental factors such as nutrition. Our combined model explains 66.6% of the variance for height, 30.4% for BMI, 44.3% for triglycerides, 30.9% for cholesterol and 30.91% for glucose.

As an illustrative example, height values show a clear increasing pattern from southeast to northwest in the MXB (Fig. ​(Fig.5a).5a). Even though height values in every state exhibit a large variance (Fig. ​(Fig.5a),5a), height is significantly correlated with longitude (Fig. ​(Fig.5b5b and Supplementary Figs. 43 and 45). We find that individuals with a higher proportion of Indigenous ancestries from the Americas are significantly shorter ($\beta =-0.45,P<2.2\times {10}^{-16}$) whereas individuals with a higher proportion of ancestries from West Africa are significantly taller ($\beta =0.07,P<0.005$; Fig. ​Fig.5b).5b). Further, considering ancestries at a finer resolution, we observe decreased height with a change in ancestries from the north of Mexico (for example, Huichol and Tarahumara) to those from the Mayan region (for example, Tojolabal and Maya; $\beta =-0.156,p=6.13\times {10}^{-6}$; Supplementary Fig. 55). Total length of ROH is also significantly associated with shorter height ($\beta =-0.08,P=0.01$). Simultaneously, younger individuals across the ancestry spectrum are taller than older individuals with the same ancestries (Fig. ​(Fig.5c),5c), exhibiting the impact of non-genetic factors (improving nutrition in the birth year range studied or effects of ageing) on height variation as well.

Obesity is a public health issue in Mexico⁴⁷ and has been suggested to be related to higher genetic risk associated with Indigenous ancestries⁴⁸. Contrary to this hypothesis, in the MXB as a whole, when considered univariately, Indigenous genetic ancestries and speaking an Indigenous language actually correlate with lower BMI (Supplementary Figs. 49 and 60). In our joint model with covariates, although those associations disappear, ROH in a genome (which are more prevalent in Indigenous genetic ancestries) are also associated with lower BMI. By contrast, as living in an urban environment is associated with higher BMI (Fig. ​(Fig.5d),5d), our results suggest a focus on factors related to an urban environment such as diet and sedentarism to help tackle the obesity issue in Mexico. Further segmented analysis considering only individuals in urban environments suggests the same: we observe individuals that speak an Indigenous language associating with higher BMI only in urban environments (Supplementary Fig. 61).

By contrast, some other traits show a correlation with an individual’s proportion of inferred genetic ancestries from the Americas: creatinine ($\beta =-0.13,P=0.0095$), LDL ($\beta =-0.141,P=0.013$), triglycerides ($\beta =0.16,P=0.001$) and blood glucose level ($\beta =0.19,P=0.0005$; Supplementary Fig. 54). In the MXB, the amount of an individual’s genome in ROH is associated with lower BMI ($\beta =-0.18,P=7.11\times {10}^{-5}$), triglycerides ($\beta =-0.13,P=0.004$) and blood glucose level ($\beta =-0.12$, $P=0.01$; Supplementary Fig. 56). We also find that polygenic scores computed using genome-wide significant SNPs from the UKB pan-ancestry GWAS are a significant predictor for complex trait variation for all traits analysed (Supplementary Fig. 57). Blood pressure is associated with environmental factors and polygenic scores but not other genome-wide genetic factors (Supplementary Fig. 53). Notably, among the environmental factors investigated, living in an urban environment is associated with higher height, BMI, cholesterol and creatinine levels, whereas living at high altitudes is significantly associated with higher triglyceride, glucose, cholesterol, creatinine and blood pressure levels (Supplementary Fig. 58). Higher educational attainment is associated with higher height, LDL, HDL and lower triglyceride levels, whereas speaking an Indigenous language is associated with lower creatinine and cholesterol levels (Supplementary Fig. 58).

Previous work has implicated the ABCA1*C230 allele (rs9282541) in decreasing HDL levels and shown that this allele is apparently exclusive to Indigenous genetic ancestries from the Americas (found in 29 of 36 Native American groups, but not in European, Asian or African individuals)⁴³. In the MXB, we similarly observe the ABCA1*C230 allele to be in higher frequencies in individuals with a higher proportion of ancestries from the Americas, and observe that individuals with higher ABCA1*C230 allele frequencies have lower HDL levels (Supplementary Fig. 59a). Nevertheless, overall, Indigenous ancestries are not associated with HDL levels after accounting for other covariates. In fact, genetic variants collectively on the Indigenous genetic background are associated with lower LDL levels (Supplementary Fig. 59b). These results illustrate how the interplay between cultural and diet factors and genetic factors are essential for different cholesterol outcomes. They also imply that although some functional variants may be specific to regions or genetic backgrounds, these are few (about 1,000 such variants estimated in the Americas with a frequency of 40% or higher from the Human Genome Diversity Project sampling of diverse genomes⁴⁹), and caution against using an individual’s global ancestry proportion as a predictor of the effect of a single functional variant. Our results overall support that functional variants with variable frequencies or environmental interactions are partially responsible for variation in a range of complex traits in Mexico⁴³.

Conclusion

Our work demonstrates the value of generating genotype–phenotype data on underrepresented groups to reveal lesser-known genetic histories and generate findings of biomedical relevance. It is also an illustration of the joint modelling of genetic and environmental effects to reveal the aetiology of complex traits and disease. In this project, we ensure diverse Indigenous and rural presence in our sampling strategy, consider the fluidity of ancestries from different local and global regions in our analyses, and evaluate their reflection in genetic and disease-relevant complex trait variation. By leveraging the largest nationwide genomic biobank in Mexico, we find diverse sources of ancestries in Mexico in light of its unique history, and infer demographic and admixture histories and ROH using ancestry-specific haplotype identity that reveal an elaborate fine-scale structure in the country. Observing a larger number of small ROH in younger individuals in the MXB and in genomic segments of Indigenous ancestries is relevant for parsing the genetic architecture of complex traits and diseases, especially those with a recessive component. We also show that demographic history affects the frequency distribution of genetic variants, thus changing how many rare variants individuals with different ancestries carry. We demonstrate the value of GWAS carried out on a resource such as the MXB for predicting complex traits. The MXB GWAS exhibits utility for polygenic score computation in independent Mexican cohorts, as well as for meta-analysis with other GWAS cohorts to increase prediction power further. Last, we observe a significant impact of genetic ancestries at different timescales, ROH, polygenic scores and sociocultural and biogeographic variables on various complex traits implicating the importance of both genetic and environmental factors in explaining complex trait variation and in considerations of potential public health interventions. Our results exhibit the added importance of considering genetic factors for preventive and personalized medicine above and beyond environmental factors. Our results will inform the design of future genetic and complex trait studies in Mexico and Latin America, and will hopefully motivate additional efforts to strengthen local research capacity across Latin America and benefit underserved groups globally.

Methods

Online content

Any methods, additional references, Nature Portfolio reporting summaries, source data, extended data, supplementary information, acknowledgements, peer review information; details of author contributions and competing interests; and statements of data and code availability are available at 10.1038/s41586-023-06560-0.

Acknowledgements

Author contributions

Study conception, design and funding: A.M.-E., L.G.-G., S.L.F.-V., A.J.M. and A.V.S.H. Genotyping and data quality control: C.D.Q.-C., A.Y.C., M.J.P.-M. and M.S. Phenotypic dataset generation and curation: L.P.C.-H., L.F.-R., G.D.-S., N.M.-R., S.C.-Q., E.F.-G. and A.Y.C. Biochemical measurements: M.T.T.-L., H.M.-M., C.A.A.-S. and L.F.-R. Groupings and classification: M.S., V.A.-A., J.N. and A.M.-E. Population genetics: M.S., M.J.P.-M., A.Y.C., C.B.-J., S.G.M.-M., A.R., C.R.G., G.L.W., J.N. and A.M.-E. Genotypic imputation: A.J.-K., A.Y.C., A.C., K.A., A.J.M. and A.V.S.H. GWAS: A.Y.C., A.J.M., K.A., A.C. and A.V.S.H. Complex trait modelling: M.S., G.L.W., J.N., A.M.-E. and L.G.-G. Manuscript writing: M.S., with input and edits from A.Y.C., G.D.-S., A.G.I., S.L.F.-V., M.T.T.-L., A.J.M., J.N., L.G.-G. and A.M.-E. Project coordination: A.M.-E. and L.G.-G.

Peer review

Data availability

The genotype and phenotype datasets for the 6,057 newly genotyped individuals from the MX Biobank Project are available at the European Genome-phenome Archive (EGA) through a Data Access Agreement with the Data Access Committee (EGA accession number for study: EGAS00001005797; datasets: EGAD00010002361 (Mexican_Biobank_Genotypes) and EGAD00001008354 (Mexican Biobank 50 Genomes)). Data can be accessed only for academic research and non-commercial use. GWAS summary statistics generated as part of this study are available at 10.5281/zenodo.7420254.

Competing interests

The authors declare no competing interests.

Footnotes

Contributor Information

Mashaal Sohail, Email: xm.manu.gcc@laahsam.

Alexander J. Mentzer, Email: [email protected].

Lourdes García-García, Email: xm.psni@ragicrag.

Andrés Moreno-Estrada, Email: [email protected].

References