Ascertainment bias

The average correlation coefficients of the first five principal components (PCs) between one given subset and the entire set with all markers are shown in Additional file 1: Figure S1. The correlation coefficients between the subset and the entire set sharply increased from 0.65 for a marker number of 500 to 0.97 for a marker number of 700. A second sharp increase emerged when the marker number increased from 800 to 1000, with a corresponding increase in the correlation coefficient from 0.97 to 0.99. Furthermore, the correlation coefficient did not significantly change when the marker number increased from 2000 to 43,252. The results indicated that 1,000 SNPs might be sufficient for population structure analyses.

Model-based population structure

The subpopulations of 1857 accessions based on the admixture model-based algorithm were analyzed in depth using the even distribution of 5000 SNPs. The results are depicted in Fig. 1. The delta K (∆K) peak was maximized when k = 2 (Fig. 1a), indicating that the accessions could be categorized into two groups: tropical/subtropical (TS) germplasm and temperate germplasm (Fig. 1b k = 2). A second peak of ∆K emerged at k = 4 (Fig. 1a), indicating that this panel could be further divided into four subgroups: SS, NSS, Modified Introduction in China (MICN), and TS I (Fig. 1b k = 4). Notably, MICN formed during the long history of maize breeding in China because Chinese maize breeders have developed a number of inbred lines derived from Chinese landraces and U.S. hybrids. These varieties significantly differ from U.S. inbreds [19]. A third peak of ∆K was observed at k = 7 (Fig. 1a), indicating that this panel could be comprehensively categorized into seven subpopulations, each including one of the following representative lines: B73, Huangzaosi (HZS), 207, Oh43, Mo17, Shen137, and some from TS regions (Fig. 1b, k = 7). Detailed information for each accession is listed in Additional file 2: Table S1.

Clustering analysis

A neighbor-joining tree was constructed based on the modified Euclidean distance and is shown in Fig. 2. The 1857 accessions were clustered into two major groups according to their origins: the TS and Tem-tropic subpopulation. The TS subpopulation contained 525 accessions, including 195 accessions from Mexico, 187 from the U.S., 77 from China, 17 from Sudan, 10 from Thailand, 9 from Canada, 9 from Tanzania, 6 from Nigeria, 3 from Somalia, 3 from Benin, 3 from Zambia, 3 from Chad, 2 from Spain, 1 from Ghana, 1 from Germany, 1 from Yugoslavia, and 1 from Egypt (Additional file 2: Table S1). The Tem-tropic subpopulation contained 1,332 accessions, which could be further clustered into four subpopulations, SS, NSS, Iodent (IDT) and TS, according to their origins and pedigrees. A further analysis showed that the accessions from these four subpopulations could be clustered into 13 subgroups, with the following representative lines: Reid Yellow Dent, Oh43, A634, 207, B37, B73, PHG39, Shen137, Huangzaosi (HZS), Ye8112, E28, GB and Mo17 (Fig. 2).

Principal component analysis (PCA)

The PCA results showed comprehensive patterns of subpopulation and a good agreement with both model-based population structure and clustering analyses (Fig. 3). The entire panel of 1857 accessions exhibited moderate differentiation and some overlap between the temperate and TS germplasm; representative lines from the TS and temperate region significantly differed, e.g., B73 from the temperate and Ki3 from the TS region of Thailand, but the accessions from the adjacent regions did not markedly differ. Which may be resulted by the lager introgression existing between temperate and tropical/subtropical accessions and lower power of PCA in population structure analysis by using only two PCs. The accessions from the temperate subpopulation were further categorized into the B73 subpopulation according to the results of model-based structure analysis (Fig. 3b) or the Ye8112, B37 and A634 subpopulations based on the results of modified Euclidean distance (Fig. 3c). Based on the pedigrees, most lines were from the U.S. and China (Fig. 3d and Additional file 2: Table S1). In addition, the TS population was further divided into the HZS, 207, Oh43, Mo17 and Shen137 subpopulations based on the model-based population structure, which corresponded to HZS, GB, Shen137, Mo17, and Reid Yellow Dent based on a clustering analysis (Fig. 3c). These subpopulations contained inbred lines of a TS lineage in their pedigrees or lines from CIMMYT, Mexico and other tropical regions (Fig. 3a and d). Moreover, many accessions were categorized into new groups, such as the PHG39, 207, A634, Oh43, B37 and E28 subpopulations; most accessions in these groups originated from regions between temperate and TS zones (Fig. 3) due to the introgression of TS genotypes into regions of temperate germplasms.

Summary statistics of genetic diversity

The accessions of the entire panel of 1857 accessions were moderately similar, with more than 96.22 % of the pair-kinship coefficients varying from 0.30 to 0.53 (Fig. 4a). The average linkage disequilibrium (LD) distance was 30 kilo-bases (kb), varying from 20 to 50 kb, with an r² exceeding 0.1 (Fig. 4b). Combining the results of both the model-based population structure and genomic variation analyses indicated pronounced patterns of genetic variation among different subpopulations. These patterns were fixed by artificial or natural selection and resulted in the division of subpopulations during breeding. The TS subpopulation was more genetically diverse than the temperate subpopulation, with gene diversities (GDs) of 0.364 and 0.284, respectively, and polymorphism information contents (PICs) of 0.281 and 0.231, respectively (Table 1). Similar trends were validated with a smaller proportion of SNPs in LD for TS when comparing with a larger proportion of SNPs in LD for the temperate subpopulation (Fig. 4c).

Index	Total	K = 2		K = 4				K = 7
		TS	Temperate	SS	NSS	MICN	TS I	B73	Mo17	PH207	Oh43	Shen137	HZS	TS II
Gene Diversity	0.365	0.364	0.284	0.301	0.361	0.306	0.348	0.268	0.299	0.360	0.294	0.311	0.272	0.345
PIC	0.291	0.291	0.231	0.244	0.289	0.247	0.279	0.219	0.242	0.287	0.239	0.250	0.223	0.277
Heterozygosity	0.046	0.048	0.025	0.027	0.058	0.047	0.037	0.023	0.033	0.065	0.028	0.049	0.034	0.033

Note: K is the number of subpopulations. ‘TS’ is an abbreviation of Tropical/Subtropical subpopulation. ‘SS’ is an abbreviation of Stiff Stalk Synthetic subpopulation. ‘NSS’ is an abbreviation of Non-Stiff Stalk. ‘MICN’ is an abbreviation of Modified Introduction of China

Genomic differentiation between subpopulations

The proportion of genetic variance due to subpopulations (D_EST) was measured to interpret the genomic variation between subpopulations (Table 2, Fig. 4(d), Fig. 5 and Additional file 1: Figure S2). The D_EST indicated different patterns of genomic differentiation between the subpopulations, ranging from 0 to 0.39 between TS and Temperate (average 0.08), from 0 to 0.45 between TS I and SS (average 0.09), from 0 to 0.45 between SS and NSS (average 0.07), from 0 to 0.41 between NSS and MICN (average 0.05), from 0 to 0.38 between MICN and TS I (average 0.06), from 0 to 0.30 between NSS and TS I (average 0.03), and from 0 to 0.57 between SS and MICN (average 0.08). The SS and TS I varieties were more differentiated, with 332 genomic regions having a large D_EST (twice the average level) (Fig. 5a). Furthermore, 250 genomic regions were highly differentiated between SS and MICN, 235 were highly differentiated between TS and Temperate, 92 were highly differentiated between MICN and TS I, 51 were highly differentiated between NSS and MICN, and 8 were highly differentiated between NSS and TS I, with a D_EST of more than twice the average level. Most importantly, 85 highly differentiated regions with a D_EST exceeding 0.2 were identified between the TS and the temperate subpopulations. Of these 85 regions, 68 were located within the interval of plant architecture or FT-related QTL, and two regions were closely linked to vgt1 and zcn8 (Additional file 2: Table S2 and S3). Furthermore, a number of special genomic regions were also found to be highly differentiated. In particular, subpopulation pairs and common regions were identified among different population pairs (Fig. 5b). In total, 303 genomic regions with a high D_EST of more than 0.2 were detected, and these regions were located within 197 FT- or plant architecture-related QTL. For example, the region containing the tag SNP PZE.108075114 differed more between the TS and temperate subpopulations and was associated with a D_EST of 0.32; this region was located within an FT-related QTL cluster and contained the flanking markers PHTi060 and bnlg1599 (Additional file 2: Table S3).

K = 2	Temperate	TS
Temperate	0.000	0.170
TS		0.000
K = 4	SS	NSS	MICN	TS I
SS	0.000	0.164	0.146	0.172
NSS		0.000	0.134	0.053
MICN			0.000	0.121
TS				0.000
K = 7	207	B73	HZS	Mo17	Oh43	Shen137	TS II
207	0.000	0.251	0.283	0.260	0.242	0.113	0.059
B73		0.000	0.295	0.309	0.246	0.202	0.247
HZS			0.000	0.242	0.194	0.199	0.208
Mo17				0.000	0.180	0.182	0.226
Oh43					0.000	0.156	0.194
Shen137						0.000	0.110
TS II							0.000

Genome-wide study of FT-related traits

The phenotypes of FT-related traits were significantly positively correlated between the environments (Additional file 1: Figure S3). Thus, the BLUPs for each accession across the three environments were calculated, and the phenotype-genotype associations were analyzed. To validate the putatively adaptive function of highly differentiated target regions, we used the FT-related traits DTT, DTS, and DTP to perform a GWAS with 43,252 SNPs as a case study. The results indicated that some highly differentiated genomic regions were associated with FT-related traits. For example, the SNP of PZE-108070380 was significantly associated with DTT (P = 7.05 × 10⁻¹⁴), DTP (P = 2.57 × 10⁻⁹) and DTS (P = 2.12 × 10⁻⁸) (Fig. 6). This SNP was located within the zcn8 gene, which is involved in maize migration from tropical to temperate regions [31]. The SNP PZE-108076585 was significantly associated with DTS (P = 1.43 × 10⁻¹¹) (Fig. 6). This SNP was located within the vgt1 gene, which is involved in maize adaptation [32]. Furthermore, twelve other SNPs were also strongly associated with FT-related traits (Fig. 6), and the regions surrounding these SNPs were more differentiated than the rest of the genome (Fig. 4d, Additional file 1: Figure S2, Additional file 2: Table S2).

Moderate SNPs are reliable in interpreting population structure division

Previous reports compared the effect of different marker systems and concluded that the subdivision of populations depended on the marker size and population [18, 33–35]. For instance, when 884 SNPs were used in one association panel of 154 inbred lines, more than 26.4 % of lines were allocated to the mixed group. This rate was higher than the 20.6 % rate identified by using 84 simple sequence repeat (SSR) markers [35]. Comparing the effect between 847 SNPs and 89 SSRs in one panel of 254 inbred lines yielded similar results [36], they proposed that many more SNPs would be required to study population structure. Here, we compared the average correlation coefficients of division for subpopulations between one given subset with different marker sizes and the entire set with all markers; we used SNPs varying from 500 to 43,252 in a panel of identical samples. The results showed that 700 SNPs are sufficient to reliably divide subpopulations in this panel, with an average correlation coefficient of the first five PCs of 0.97 between the subsets and the entire set of SNPs. The average correlation coefficient could be increased to 0.99 by increasing the number of SNPs to 1000 (Additional file 1: Figure S1). Yu et al. [37] reported moderate genetic diversity with a PIC of 0.24 for a sample size of 274. We herein report a similar PIC of 0.29 for a sample size of 1,857. Yu et al. [37] demonstrated that more than 1000 SNPs are necessary needed to robustly estimate the genetic differentiation of a naturally diverse panel, and this requirement exceeded the 700–1000 SNPs found to be necessary herein. Thus, a larger sample size is expected to significantly improve the detection power of subdivisions in the populations. These results were consistent with those reported by Morin et al. [33], who compared the subpopulation differentiation for sample sizes ranging from 10 to 100. The results reported herein suggested that a moderate SNP marker number (700–1000) is sufficient to divide population structures in this panel.

Comprehensive patterns of population structure in maize inbreds worldwide

The analysis of population structure is an important step in dissecting the genetic basis of complex traits via association analyses [38]. Such an analysis can result in false positive errors [34]. In the last several decades, a number of studies have evaluated the population structure of specific germplasms using limited sample sizes and sources. These studies independently provided specific information about the subpopulation differentiation of approximately 600 Minnesota maize germplasms [22], 172 Dent germplasms from Hohenheim [39], 400 maize lines from North America [23], 367 elite lines from China [19] and 527 lines representing TS and temperate backgrounds [40]. Here, we integrated maize germplasms from America, Africa, Europe and Asia, including 1857 accessions from more than 16 countries worldwide, to investigate subpopulation differentiation. The outputs of STRUCTURE V2.3.3 identified seven subpopulations: including B73, HZS, 207, Oh43, Mo17, Shen137, and TS II (Fig. 1). These results provided much more information about maize subpopulation differentiation than previous studies. In fact, the B73 (SS), Mo17 (NSS), Oh43, and 207 (IDT) subpopulations were identified using SSR markers and an Illumina MaizeSNP50 Beadchip [22]. HZS (TSPT), Shen137 (PA derived from Pioneer hybrid 78599), and TS I subpopulations were also identified in previous reports [18, 19, 41]. In addition, the findings this study was also consistent with known pedigrees. For example, LH61 shared 87.5 % of its nuclear genetic material with Mo17 [42] and clustered into the Mo17 subpopulation with an ancestry membership of 0.91 (Additional file 2: Table S1). These results were consistent with those reported by Lorenz et al. [42]. Furthermore, the clustering analysis identified many more clusters, including Mo17, GB, E28, Ye8112, HZS, Shen137, PHG39, B73, B37, 207, A634, Oh43, and Reid Yellow Dent (Fig. 2). The identification of these clusters indicated that our clustering analysis increased the resolution of the categorization of accessions into subpopulations compared with the model-based method, which commonly identifies six subpopulations, Mo17, B73, HZS, Oh43, 207, and Shen137. For instance, PB80 and A632 shared 75 % and 93.75 % of the nuclear genetic material of B73 and B14, respectively [42], these two lines clustered into the same subpopulation as B73 and B14, respectively. This clustering was consistent with a report by Lorenz et al. [42]. Most importantly, the clustering analysis in this study identified new subpopulations that are represented by the following lines: GB, E28, Ye8112, PHG39, B37, A634, and Reid Yellow Dent. These lines correspond to the following heterotic groups: Chinese Landrace (GB) [19], Ludahonggu (E28) [41], PB (Ye8112, B37) [19] derived from modern U.S. hybrids, Commercial hybrid-derived lines (PHG39, A634) [10], and U.S. landrace (Reid Yellow Dent) [10], respectively. Of these groups, Chinese Landrace is mainly distributed in the northeast and southwest of China, and this variety originated from the North-American Mid-West and Mexican highlands, respectively [3]. These landraces yielded new subpopulations and have been widely used in maize-breeding programs [19]. For example, E28 is a representative line derived from crossing the landrace Ludahonggu with modified introduction lines according its pedigree [19]. Ye8112 was selected from the hybrid “8112”, which originates in the U.S. [41]. Some of the lines were derived from this line, such as Ye478 and 488, which were clustered in the heterotic group of PB [19, 41]. A634 was derived from the MN13 lineage [22], is highly utilized in U.S. hybrid maize breeding. This line constituted 4.2, 7.8, and 3.0 % of the total U.S. seed requirement in 1970, 1975, and 1979, respectively, and lots of lines were derived from A634 [13]. B37 is an important public line that was widely used to develop Pioneer hybrids during the 1980s [6]. The selection of a second cycle line from Pioneer hybrids resulted in new lines, which formed a subpopulation represented by B37. PHG39 is a representative inbred Maize Amargo germplasm line from which many protected corn lines have been developed. Furthermore, several important first cycle recombinant lines derived from PHG39 have been considered for commercial maize breeding [10]. These results provide maize breeders with more definitive information to effectively use historical genetic resources while maintaining the heterotic patterns necessary for hybrid breeding.

Genomic differentiation and putative functions

Genomic differentiation between subpopulations is a fundamental challenge in population genetics. Maize originated in tropical central-Mexico and rapidly spread to colder, temperate regions worldwide [32]. This diffusion caused maize to adapt to local environments by developing traits that allowed it to thrive in these environments, i.e., changes in FT and plant architecture. These changes allowed maize to reach maturity within different growing seasons. Some studies have documented the pair-wise F_ST between subpopulations while considering genomic differentiation [9, 43, 44]. Schaefer and Bernardo [22] reported an average pair-wise F_ST of 0.165 for one diverse panel of 284 maize inbreds; this value ranged from 0.054 between the A321 and Oh43 subpopulations to 0.325 between the Mo17 and B73 subpopulations. Romay et al. [12] found that most germplasms from classic breeding programs of the Corn Belt were closely related, with an average pair-wise F_ST of 0.04. However, the differentiation regions and putative function remain poorly understood. Moreover, the D_EST was also demonstrated as a measure genomic differentiation. This parameter relies on the genotypic rather than allelic number and is corrected for heterozygosity [27]; values close to zero indicate little differentiation, and values close to unity indicate nearly complete differentiation. Therefore, the D_EST was used in the present study to evaluate the genomic variation between the subpopulations, and the results of this analysis revealed strong differentiation among the subpopulations. This differentiation was attributed to the continuous fixation of target genomic regions within subpopulations and strong isolation between subpopulations during maize breeding practices. The pair-wise D_EST between the TS and the temperate subpopulations was 0.17 (Table 2), and 235 highly differentiated genomic regions were identified (Fig. 5). Most adaptive traits were selected and fixed during maize’s long evolution and adaptation from tropical to temperate climates [31]. This fixation caused the high genomic differentiation between TS and temperate germplasms (Table 1, Figs. 2 and 3). Interestingly, 85 strongly differentiated genomic regions with a D_EST exceeding 0.2 were identified between the TS and the temperate subpopulations. A genetic analysis showed that these 85 regions comprise 79 % of the genetic variation of this panel (Additional file 1: Figure S4). Of these regions, 15 were significantly associated with FT-related traits based on GWAS (Fig. 4d and Additional file 1: Figure S2). In addition, two significant QTNs were closely linked to zcn8 and vgt1 (Fig. 4d), which are involved in maize migration and adaptation from tropical to temperate climates [31]. Beyond that, 66 highly differentiated regions were located within the interval of plant architecture or FT-related QTL (Additional file 2: Table S3). In addition, 159 highly differentiated genomic regions were also identified between SS and NSS subpopulations, with a D_EST exceeding 0.16 (Fig. 5). Furthermore, 15 regions located within FT- or plant architecture-related QTL were also identified (Additional file 2: Table S3). This finding was consistent with the marked distance between SS and NSS (Figs. 1, 2 and 3). SS and NSS are two major heterotic groups used in U.S. breeding programs that are represented by the lines B73 and Mo17, respectively. Previous studies also reported a significant distance between these two groups [23]. Furthermore, other highly differentiated genomic regions between specific subpopulation pairs were also identified, and these regions were located within a number of QTLs associated with FT- or plant architecture-related traits mapped using different bi-parental populations (Additional file 2: Table S3). In total, 303 genomic regions with a high D_EST of more than 0.2 were detected, and these regions were located within 197 FT- or plant architecture-related QTLs. For example, the region containing the tag SNP PZE.108075114 was more differentiated between TS and the Temperate subpopulations (D_EST = 0.32), and this region was located within one FT-related QTL cluster that contained the flanking markers PHTi060 and bnlg1599. These results indicate genomic regions of interest for the formation of given subpopulations and provide new insight into the dissection of the genetic basis of complex traits.