Phenotyping for morpho-physiological, yield and its related traits:
The RIL population showed wide variation for the various morpho-physiological traits and yield related traits such as relative water content (RWC), normalized difference vegetation index (NDVI), canopy temperature depression (CTD), quantum yield (fv/fm), net CO2 assimilation rate (AN), stomatal conductance to water vapor (gs), transpiration rate (TR), leaf temperature (LT), anthesis-silking interval (ASI), plant height (PH), cob height (CH), cob length (CL), number of kernel rows (NKR), number of kernels per row (NKPR), cob weight (CW), total biomass (TB) and grain yield / plant (GY) under well-watered (WW) and water-deficit (WD) stress conditions (Supplementary Table S1a, b). All the traits were affected by water-deficit stress. Frequency distribution histograms of morpho-physiological and yield traits under WW and WD stress conditions were given in Fig. 1. All these traits showed near normal distribution except quantum yield which was negatively skewed under water-deficit stress conditions while stomatal conductance was positively skewed under well-watered conditions. Based on the significance values obtained using Kolmogorov–Smirnov test the traits CTD, NDVI, AN, LT, PH, CH, CL, NKPR, CW, GY and TB showed normal distribution under both well-watered and water-deficit stress conditions while the trait gs and TR showed normal distribution under WD stress conditions (Supplementary Table S2). This indicates that selected RILs captured the genetic variability of whole mapping population to be utilized for QTL identification for these traits. The descriptive statistics of mapping population under WW and WD stress for morpho-physiological and yield related traits including coefficient of variation (CV%), skewness, kurtosis and heritability were given in Supplementary Table S3. Moderate broad-sense heritability (H2) was observed for all the traits with a range of 0.45 – 0.72 in the population. A wide range of coefficient of variation was observed among all the traits.
The pooled ANOVA of trials across seasons (rainy and post rainy) showed significant interaction between seasons, treatment and genotype for most of the traits. For traits RWC, NKPR and TB, all interaction effects were significant except season × treatment (Supplementary Table S4). For the trait ASI, the interactions effect for treatment × genotype and season × treatment × genotype was not significant. The significant genotype, environment and their interaction variances for traits indicated that these traits were influenced by both genetic and non-genetic factors. Simple correlation coefficient analysis revealed significant positive correlations between NDVI and PH, CH, CW, GY and TB, AN and gs, TR, gs and TR, PH with CH and CL, CH and TB, CL and NKR, NKPR, GY and TB, NKR and NKPR, CW with GY and TB, GY with TB under both well-watered and water-deficit stress conditions. Similarly, significant positive correlations between QY and gs, TR, LT and TB under WD stress conditions was observed (Supplementary Table S5).
QTL mapping using high throughput SNP genotyping
The number of raw SNPs, polymorphic SNPs between the parents after filtering with MAF ≤ 0.05% and their distribution in various chromosomes, number of mapped SNPs in linkage map and average marker interval was presented in Table 1. The maximum number of markers (219) were detected on chromosome (chr) 2 and minimum number of markers (65) on chr 10. A total of 40 QTLs were associated with traits analysed, of which 26 were detected under well-watered and 14 under water- deficit stress conditions. Among these 14 major QTLs capturing ≥15% phenotypic variation were identified for the traits viz., RWC, NDVI, CTD, AN, gs, TR, PH, CH, CL, NKPR and CW. The LOD score values for these QTLs varied between 2.54 and 6.07 with phenotypic variation captured varied between 15.05 to 22.34% (Table 2). While, 26 minor effect QTLs capturing <15% phenotypic variation were identified for the traits viz., AN, ASI, CH, CTD, gs, GY, LT, NDVI, NKPR, NKR, QY, TB and TR. The LOD score values for these QTLs ranged between 2.52 and 3.37 with PVE% ranging between 6.75 to 14.91 (Supplementary Table S6). One major QTL for the trait RWC was detected on chromosome 9 with a LOD score of 2.71 and capturing 16.3% PVE. For the trait NDVI one major and minor QTL were detected on chr 2 having LOD scores of 3.93, 3.19 and capturing PVE % of 18.59, 14.02 respectively. Two minor QTLs were detected for the trait QY on chr 1 and 7 with a LOD score of 2.86, 2.81 and PVE% of 14.62, 14.47 respectively. For net CO2 assimilation rate (AN), a major and a minor QTL were detected on chr 3 (LOD score of 3.24, 2.52) capturing 17.24, 14.25 PVE% respectively. For the trait stomatal conductance to water vapor two major QTLs were detected on chromosome 6 and three minor QTLs one on chr 3 and two on chr 7. The LOD ranged from 2.54 to 4.45 and PVE% ranged from 9.69 to 18.21. For the trait transpiration rate one major QTL was detected on chromosome 7 (LOD: 6.07; PVE%: 21.47) and two minor QTLs one each on chr 1(LOD: 2.87; PVE%: 10.81) and 3 (LOD:2.72; PVE%: 10.11).
Among yield traits, one major QTL was detected for plant height on chr 8 (LOD score of 3.19 and PVE% of 17.98). For cob height one major QTL and two minor QTLs were detected on chromosome 1 with a LOD score of 3.49, 2.52, 2.90 and PVE% of 15.05, 13.84, 11.68 respectively. For total biomass, three minor QTLs, two on chr 2 and one chr 5 were detected with a LOD score of 2.53, 2.58, 2.55 and PVE% of 13.51, 11.50, 13.04 respectively. For cob weight one major QTL was detected on chr 2 under both WW and WD stress conditions. For grain yield two minor QTLs one each on chr 2 and 6 were detected (LOD score of 2.78, 2.54) capturing PVE% of 10.54, 11.33 respectively.
Identification of co-localized and QTL hotspot region
Out of the 40 QTLs, three co-localized QTLs were identified. For traits TR and ASI two QTLs (qTR1-1, qASI1-1) were co-localized at the marker interval rs128441140 - rs128842621 on chr 1 with LOD score of 2.87, 2.61 and capturing phenotypic variation of 10.81, 13.92 % respectively, while, gs and GY QTLs (qgs6-2 and qGY-6-1) were co-localized at the marker interval S6_89546385 - S6_179562549 on chr 6 with LOD score of 2.54 and capturing phenotypic variation of 15.86, 11.77% respectively. For the traits gs and TR, QTLs (qgs7-1, qTR7-1) were co-localized at the marker interval S7_166501967 - rs130671858 on chr 7 with LOD score of 2.73, 6.07 and capturing phenotypic variation of 11.3, 21.07 % respectively.
In the present study, we have also identified three QTL hotspots regions. The marker interval, S3_169283017- rs129386882 (15cM) on chromosome 3 encompassed three QTLs for the traits net CO2 assimilation rate and transpiration rate qAN3-1, qTR3-1 and qAN3-2 capturing 10.10 to 17.24 % phenotypic variation (Fig. 2). While, the marker interval, rs130916094 - rs849450935 on chromosome 8 (27cM) encompassed two QTLs for cob height and plant height qCH8-1 and qPH8-1 with LOD scores 2.78, 3.19 and the PVE% 14.9, 17.98% respectively (Fig. 3). The region between the marker interval, S2_213060766 - S2_14679066 on chromosome 2 (33cM) encompassed three QTLs for the traits cob weight, total biomass and grain yield qCW2-1, qTB2-1, qTB2-2 and qGY2-1 having LOD scores of 2.53 to 4.23 and capturing 10.54 to 22.33% phenotypic variation (Fig. 4).
Epistatic interaction among QTLs
Seventeen significant digenic epistatic QTLs for traits RWC, QY, AN, gs, TR, LT, ASI, CL, CW and TB were detected (Table 3; Supplementary Fig. S1a, b). The LOD score values ranged between 5.0 to 5.91 while the PVE% ranged between 11.42 – 37.43%. The negative epistatic values (add by add) indicated higher epistatic effect of the recombinant genotype than the parental genotype. The epistatic effect values of QTLs for traits TB, AN, RWC, CL, LT, and QY were negative while for traits gs, QY, ASI, TR, CL, AN and CW was positive. The genomic region on chromosome 2 between the marker interval of rs129243511-S2_229825946 showed epistatic interaction for CW on chromosome 8, located between the marker interval of S8_145945904 - S8_148216407 and this interaction contributed to 37.43% of phenotypic variation. The genomic region on chromosome 7 between the marker interval of rs131600615 - S7_50455390 showed epistatic interaction for RWC on chromosome 9, located between the marker S9_41587672 - rs131737404 contributing 33.03% of phenotypic variation. Two epistatic interaction were identified for AN. First the on chromosome 1 between marker S1_3798004 - rs128441140 showed epistatic interaction with region on chromosome 3, located between rs129555629 - S3_183844216. This interaction contributed to 17.60 % of phenotypic variation. Second, on chromosome 2 between markers, rs727228961- rs276685886 showed epistatic interaction for AN with chromosome 3, located between interval of S3_2276048 - S3_89193637 which contributed to 16.71% of phenotypic variation. The genomic region on chromosome 2 between the marker interval of rs131971876 -rs129196105 showed epistatic interaction for TR with the region on chromosome 7, located between S7_131791490 - S7_136261108 and contributed to 29.55% of phenotypic variation.
QTL-by-Environment Interaction Analysis
Using the MET (multi-environmental trials) module of ICIM a total of 78 QTLs were identified for different traits (Supplementary Table S7). Of these 32 QTLs were common to those identified by ICIM-ADD method. One QTL each for traits gs, ASI, GY and two QTLs for trait LT which were identified by ICIM-ADD method but were not detected by MET analysis. The MET analysis identified 46 QTLs for traits viz., RWC, QY, AN, CH, ASI, TB, GY, NDVI, LT, CW and PH which were not detected by ICIM-ADD method. The phenotypic variation captured by additive and dominance effects [PVE (A)] varied from 3.42-14.36 % and the PVE (A by E) (additive and dominance by environment effects for corresponding QTLs) varied from 0 - 7.92 %. The PVE (A by E) was significantly lower than PVE (A). For traits AN, gs, NDVI, NKPR, QY and TR significant differences were observed in the QTL × E interaction effect namely, the LOD (AbyE), PVE (AbyE) indicating the significant effect of the environment (Table 4).
Identification of candidate genes in the genomic region spanning QTLs
The major and minor QTLs associated with morpho-physiological and yield related traits, their chromosome position, locus ID, start and end position, SNP type and position, size of the locus, genes and their biological / molecular function was given in Table 5, Supplementary Table S8. The number of genes identified for the QTLs of various traits ranged between 1-3. The genes encoded for different traits belonged to the categories signal transduction (RWC - Zm00001eb395310 - Guanine nucleotide exchange factor SPIKE 1; GY - Zm00001eb297570 - Protein-serine/threonine phosphatase), transcription factors (NDVI - Zm00001eb118010 - G2-like-transcription factor 27; CTD - Zm00001eb325260 - AP2-EREBP-transcription factor 201; ASI - Zm00001eb295810 - NAC type transcription factor (NAC87), transporter activity (AN - Zm00001eb146040 - Chloride channel protein; gs - Zm00001eb324180 - Sugar carrier protein C; TR - Zm00001eb015510 - Phospholipid-transporting ATPase; CH - Zm00001eb363270 - Calcium-transporting ATPase; NKPR - Zm00001eb076560 - GDP-mannose transporter (GONST5), cell wall biosynthesis and organization (TR - Zm00001eb144960 - Lipoxygenase; TR - Zm00001eb145080 - Pectin acetylesterase), photosynthesis (gs, TR -Zm00001eb324240 - Chlorophyll a-b binding protein, chloroplastic) and Carbon utilization (PH - Zm00001eb359190 - Carbonic anhydrase; CH - Zm00001eb002270 - Glyceraldehyde phosphate dehydrogenase B1).