A panel of 299 maize inbred lines, as reported by Yang et al. (2011), was planted in the Sanya Experimental Station of Beijing Academy of Agriculture and Foresting Sciences in Hainan, China in 2018. Irrigation, fertilization, and pest or disease control were performed according to Zhang et al. (2021).

The young developing kernels were dehydrated before X-ray tomography analyses. The dehydration processes were conducted according to Zhang et al. (2021) with slight modifications. Briefly, young kernels were harvested and immersed in formaldehyde–acetic acid fixative containing 50% ethanol, 5% glacial acetic acid, and 3.7% (v/v) formaldehyde immediately. After fixation, the samples were dehydrated through a sequential ethanol series from 50, 70, and 95–100% at a 24-h interval for each gradient. Then, samples were then transferred to tertiary butyl alcohol for 24 h before being frozen at −80°C overnight. Thereafter, the frozen kernels were freeze-dried in a freezer dryer machine (LGJ-10E, China). The mature and naturally dried kernels were scanned directly by X-ray μCT without the drying process. Three natural dried mature kernels harvested from the middle part, 1/3 of ear length in the ear central (Shen et al., 2018), of a representative ear were collected from each inbred line for phenotypical analysis.

A Skyscan 1172 X-ray μCT system (Bruker Corporation, Billerica, MA, United States) was used to acquire spatial information on the inner cavity of maize kernels of the 299 inbred lines with a unified setting, as described by Hou et al. (2019). Briefly, kernels were scanned with the imaging pixel of 13.5 μm in the 2K scanning mode (2000 × 2000 pixels). The sample was scanned over 180° rotation with an image taken at every 0.2°. For each sample, the X-ray μCT scan generates about 900 images. To further quantify the size (volume) of the kernel’s inner cavity, a series of 2D X-ray images was converted into a 3D image using the X-ray system integrated software Nrecon, which transferred the gray level of the raw image (16 gray level) into an 8-bit gray-level image ranging from 0 to 255, corresponding to white and black, respectively.

The principle of X-ray μCT analyses is based on the differences in X-ray attenuation arising from differences in material density and composition (Chawanji et al., 2012). The average attenuation of the sample is expressed in Hounsfield Unit (HU) or CT number (Donis-González et al., 2014). Thus, the differences in physical density or constituents are observed as changes in the HU or CT number. For example, the air has a low HU of −1,000, whereas a solid material with high density may have a HU value up to 3,000 (Donis-González et al., 2014). In the X-ray μCT image, the attenuation variance of different structural components within a sample was indicated by gray-level intensities, where the denser the region, the higher the attenuation, and the brighter the region appears on the image.

X-ray μCT analyses of grains from the maize inbred line revealed that the gray level of the cavity region in the grain was less than 50. Hence, the gray intensity of the X-ray image that was less than 50 was used as a threshold to calculate the EC size of the kernel. In this study, the batch function of CTAn software associated with the X-ray μCT system was employed for the high-throughput calculation of the cavity size for the 299 inbred lines. A 3D image processing software, ScanIP (Simpleware Ltd., Exeter, United Kingdom), was used to measure the volume of each cavity in maize grains as described by Hou et al. (2019).

The GWAS was performed as described by Lu et al. (2021) and Zhang et al. (2021) with minor modifications. Briefly, genotype data of the 299 inbred lines were obtained from Professor Yan Jianbing’s laboratory of Huazhong Agricultural University (URL¹). A total of 779,855 SNPs with a minimum allele frequency (MAF) greater than 0.05 and a call rate greater than 0.9 were used for GWAS analysis (Biscarini et al., 2016). A multi-locus random-SNP-effect mixed linear model tool (R package “mrMLM” version 4.0) coupled with the population structure (Q) and kinship (K) data was used to test the statistical association between the cavity trait and the genotypes. Six multi-locus GWAS methods (mrMLM, FASTmrMLM, FASTmrEMMA, ISIS EM-BLASSO, pLARmEB, and pKWmEB) were included in the “mrMLM” analyses. These six ML-GWAS methods were processed in two steps. First, each SNP on the genome was filtered with a p ≤ 0.5/N, where N is the total number of genome-wide SNPs. Then, all the SNPs that were potentially associated with the trait were included in a multi-locus genetic model further screened with a defeat p-value of 0.0002 to declare a significance of SNPs that were associated with a given trait. Only SNPs that were identified by more than two multi-locus GWAS methods were regarded as ‘‘Top’’ SNP. Annotation of SNPs was conducted by using ANNOVAR software against the maize B73 reference genome (B73 RefGen_v4) available in Ensembl Plants² and NCBI Gene database³.

The histological section of developing kernels staged at 8, 10, 12, and 14 DAP were prepared according to Zhao et al. (2018) with minor modification. The kernels were harvested freshly from the middle of the ear and immersed immediately in the formaldehyde–acetic acid fixative containing 50% (v/v) ethanol, 5% glacial acetic acid, and 3.7% (v/v) formaldehyde at 4°C overnight. Then, the samples were dehydrated through a gradient of ethanol, xylene, and embedded in paraffin. Samples were sectioned with a microtome at 8 μm in thickness (Leica Microsystems, Wetzler, Germany), stained with toluidine blue O (TBO), and pictured by Zeiss microscopy (Jean, Germany) with a CCD Camera.

JMP14 software was used in this study for data analysis. Detailed methods are captioned in the figure legend.

To investigate when and where the EC is formed, we performed a detailed X-ray μCT analysis on maize kernels starting from 8 DAP at 2-day intervals. Figure 1A represents a sagittal view of mature maize grain with the corresponding schematic 3D sectional view on the right (Figure 1B) from B73, a standard maize line from which the reference genome was produced (B73 RefGen_v4). Different tissues of the kernel were distinguished by different levels of the gray intensities, with the brightest and the darkest gray regions representing the embryo and the cavities, respectively (Figure 1A). Here, the embryo area represents a high level of attenuation indicating dense structure, while the cavities had the lowest attenuation, thus showing the darkest area.

As shown from the X-ray μCT imaging, the grain cavities consist of EC, embryo cavity, (EmC), and subcutaneous cavity (SubC, the gap between the pericarp and the endosperm), with EC being the largest among the three (Figure 1). It is worth noting that the EC described here is different from that of the EC in wheat. The latter is designated as a gap between the endosperm and the maternal tissue (Chateigner-Boutin et al., 2021) from which nutrients are taken up by the filial tissues. By contrast, the EC in maize grain represents a hollow space within the basal region of the endosperm (Figure 1A).

We then examined the temporal progression of EC formation by conducting X-ray μCT analyses on B73 maize kernels at stages 8, 10, 12, and 14 DAP. In the maize grain, the cellularized endosperm starts to accumulate starch from 8 DAP onward (Zhan et al., 2015). Consistent with this physiological feature is the increase in gray intensity in the endosperm from 8 to 14 DAP, reflecting the increasing accumulation of storage products and tissue density (Figures 1C–F). The overall dark appearance of the grain filial tissue at 8 and 10 DAP (Figures 1C,D) indicates low structural density, reflecting the low level of endosperm starch accumulation at these stages.

Surprisingly, in contrast to the overall increase of gray intensity in the endosperm as starch accumulates, a dark hollow space (cavity) appeared at 12 DAP in the upper part of the starch endosperm (Figure 1E), which enlarged further in the 14 DAP kernel (Figure 1F). To further confirm the timing and location of the emergence of this EC, the same set of maize kernels used for X-ray scanning was sectioned and stained with toluidine blue (Figures 1G–I). The histological analyses displayed an overall increase in endosperm cellular density from 8 to 14 DAP (Figures 1G–J). Most notably, the analyses revealed the EC was formed in the endosperm at 12 DAP and further expanded at the 14 DAP (Figures 1I,J). These observations are in agreement with that from the X-ray μCT analyses (Figures 1C–F). A magnified view of the indicated regions revealed the accumulation of starch granules in the endosperm cells at 10, 12, and 14 DPA but their absence at 8 DAP and the EC regions at 12–14 DAP (Figures 1G–J). Collectively, the X-ray μCT and histological analyses established that the EC was formed in the central starch endosperm at 12 DAP and enlarged at 14 DAP.

While the X-ray μCT system provides a solution to visualize the inner structure of the maize grain by applying a gray-level threshold (less than 50) to discriminate the cavity region from the surrounding non-cavity region, which enables the quantification of the total cavity size within a grain, it cannot specifically measure the volume of individual cavities such as EC due to the computing algorithm setting in the associated CTAn software. By incorporating the ScanIP software, we have recently developed an X-ray μCT-based micro-phenotype analysis process, which enables the measurement of individual cavities, although still a time-consuming process (Zhao et al., 2021). By analyzing grains from 11 different varieties differing in EC sizes revealed that EC is the predominant cavity in the maize kernel, which accounts for over 50% of the total cavity in the majority of the varieties (Zhao et al., 2021). This observation, together with the finding that the hardness of maize kernel is negatively correlated with inner cavity sizes (Guelpa et al., 2015), strongly indicates that the presence of EC reduces the hardness of maize kernel, thereby increasing its vulnerability to breakage.

We then performed a correlation analysis among these 11 varieties, which showed that the EC size closely correlated with the total cavity size of the entire kernel (sum of Ecs, EmC, and SubC), with R² = 0.8961 (Figure 2). Given that it is widely accepted in genetic studies on crops including maize to apply a description factor highly correlated with a target phenotype as a proxy for GWAS (e.g., Mantilla-Perez et al., 2020), we next used the total cavity volume of the entire kernel as a proxy for the high throughput calculation of EC size for the GWAS analyses as described below.

To identify potential loci and candidate genes linked to the development of EC, we used the data of total inner cavity size from the 299 inbred maize lines (Supplementary Table 1) as a proxy for the EC trait to perform GWAS. According to the population structure analysis conducted by Yang et al. (2011), the 299 maize panel was classified into four subgroups with 131 lines in the tropical–subtropical (TST) group, 14 in Stiff stalk (SS), 81 in the non-stiff stalk (NSS), and 73 in a mixed subpopulation. Overall, the inner cavity size of the maize kernel showed a wide variation in the association panel. The inner cavity size of the TST subgroup was significantly higher than the SS and the NSS subgroup (Figure 3).

The GWAS analysis was conducted using the multi-locus random-SNP-effect mixed linear models in R package “mrMLM” (version 4.0) together with six different association methods (mrMLM, FASTmrMLM, FASTmrEMMA, ISIS EM-BLASSO, pLARmEB, and pKWmEB) to balance the false positives and false negatives. Only the SNPs that were identified by more than two methods with p ≤ 6.4e-7 were considered as significantly associated with the target trait. A total of 10 SNPs were identified, of which 1 SNP was identified by five, 1 SNP by four, 2 SNPs by three, and 6 SNPs by two methods (Table 1). For example, the SNP chr2.S_7040126 on chromosome 1 was co-detected by five GWAS methods. Thereafter, these 10 unique SNPs were used to identify candidate genes, of which 9 SNPs were located in the intergenic regions, and 1 was situated in the intragenic region (Table 1). Here, a 100-kb window flanking 50-kb upstream and downstream of each intergenic SNP was defined as the candidate gene selection region (as shown in the schematic presentation of Figure 4). The selection of the 100-kb window size was based on a similar GWAS population in maize for the identification of candidate genes associated with a given locus (Li et al., 2013). Two flanking genes that were most close to the SNP were selected as the candidate genes. The intragenic SNPs were referred to as those located within a given candidate gene. For example, the chr2.S_13594672 is located in the exon region of the Zm00001d002477 gene, which was thus selected as the candidate gene.

In total, 19 candidate genes were identified through this exercise (Table 1). All the candidate genes were annotated according to the maize B73 reference genome (B73 RefGen_v4) available in Ensembl-Plants and NCBI Gene database. Among them, eight genes were annotated to proteins with unknown function, while the other 11 out of 19 candidate genes were annotated encoding specific proteins (Table 1). Moreover, data-mining of qTeller⁴ expression profiles revealed that 9 out of the 11 candidate genes were expressed during the EC formation (Figure 5). These include genes that encode proteins functioning as enzymes, transcription factors, and those used for gibberellin biosynthesis, as described below.

Here, candidate gene Zm00001d001164 encodes a glucuronosyltransferase, which functions in mitochondrial electron transport transferring electrons from succinate to ubiquinone. Zm00001d002166 encodes the eukaryotic translation initiation factor 3 subunit C (eif3c) that participates in most of the translation initiation processes. In Arabidopsis, mutation analysis of Ateif3c1 caused embryo development defects, leading to seed abortion (Raabe et al., 2019). Zm00001d045563 (Dwarf3) encodes cytochrome P450 enzyme, which functions in the early step of gibberellin synthesis (Winkler and Helentjaris, 1995). Zm00001d013005 encodes a ring-type E3 ubiquitin transferase. It functions in the final step of the ubiquitination process (Zhou and Zeng, 2017) and is critical to the plant immunity system (Park et al., 2016). Zm00001d013006 was annotated to encode a DNA gyrase subunit A (GyrA) of DNA topoisomerase, which is a key enzyme that controls the topological state of DNA (Champoux, 2001). DNA gyrase is the type II topoisomerases, essential for DNA replication and transcription (Champoux, 2001). Interestingly, one of the candidate genes, Zm00001d045097 was found to encode a multidrug resistance-associated protein 11 (MRP11), which belongs to the subfamily of ATP-binding cassette transporter (ABC-transporter) family. Previous studies revealed that the MRP proteins play important roles in plant detoxification by sequestrating toxic metabolites from the cytoplasm to the large central vacuole, thereby preventing cellular damage (Klein et al., 2006). For example, in maize, ZmMrp3 was localized to tonoplasts for anthocyanin transportation (Goodman et al., 2004). Zm00001d013173, which was annotated as a Phloem protein 2 (PP2) like gene, encodes a PHLOEM PROTEIN 2-LIKE A10 (PP2-A10). Despite little information available on PP2-A10 in maize, its Arabidopsis ortholog AtPP2-A10 was identified as a member of the Arabidopsis lectin family gene (Dinant et al., 2003; Naithani et al., 2021). Lectin family proteins possess at least one carbohydrate recognition domain (CRD), which enables reversible binding to carbohydrates. This is critical for maintaining sugar-related osmotic balance in the plant cell via endogenous lectin cycle, in which lectin binds to sugar to form a sugar-complex rendering them osmotically inert when sugar becomes surplus and release sugar for utilization when needed (Nonomura et al., 2020). Zm00001d010919 encodes a protein kinase superfamily protein, Pto kinase interactor1 with unknown biological function, whereas Zm00001d007403 encodes a chalcone synthase (CHS), which functions in the synthesis of maysin in maize, which is a C-glycosyl flavone found in maize silk tissue that confers resistance to corn earworm (Meyer et al., 2007).