Non-escaping frost tolerant QTL linked genetic loci at reproductive stage in six wheat DH populations

Jingjun Zhng, MD Shhiul Islm, Yun Zho, Msoo Anwr, Zi Alhr, Moyun Sh,Ronghng Yng, Angl Juhsz, Guixing Tng, Jinshng Chn, Hng Liu, Ynji Jing,Shngnn Zhi, Xin Hu, Junkng Rong, Yingqun Zhng, Yo Qin, Qir Liu, Zitong Yu,Yujun Zhng, Sgh Blot, Mirz Dowl, Soni Arin, Nnit Roy, M Rs Mllik,M Atik Us Si, Shnji Rhmn, Nigrin Sultn, Srh Al-Shikh Ahm, Chris Floris,Ki Chn,Drshn Shrm,Nthn Hight,Bn Biulph,Miqin Lu,Jorg E.Myr,Wujun M,*

a Australian-China Joint Centre for Wheat Improvement, Western Australian State Agriculture Biotechnology Centre, Food Futures Institute, Murdoch University, South Street,Murdoch, WA 6150, Australia

b College of Agriculture & Food, Zhejiang Agriculture & Forestry University, Lin’an 311314, Zhejiang, China

c SAGI West, Faculty of Science and Engineering, Curtin University, Bentley, WA 6102, Australia

d Department of Primary Industries and Regional Development, Western Australia, 3 Baron-Hay Court, South Perth, WA 6151, Australia

e Australian Grain Technologies, 12656 Newell Highway, Locked Bag 1100, Narrabri, NSW 2390, Australia

f Ag RD&IP Consult P/L, Nerang, QLD 4211, Australia

Keywords:Anthesis-related genes Wheat DH populations Reproductive stage Non-escaping frost tolerance QTL

ABSTRACT Reproductive stage frost poses a major constraint for wheat production in countries such as Australia.However, little progress has been made in identifying key genes to overcome the constraint.In the present study, a severe frost event hit two large-scale field trials consisting of six doubled haploid (DH)wheat populations at reproductive stage (young microspore stage) in Western Australia, leading to the identification of 30 robust frost QTL on 17 chromosomes.The major 18 QTL with the phenotype variation over 9.5% were located on 13 chromosomes including 2A, 2B, 2D, 3A, 4A, 4B, 4D, 5A, 5D, 6D, 7A, 7B and 7D.Most frost QTL were closely linked to the QTL of anthesis,maturity,Zadok stages as well as linked to anthesis related genes.Out of those,six QTL were repetitively detected on the homologous regions on 2B,4B,4D,5A,5D,7A in more than two populations.Results showed that the frost damage is associated with alleles of Vrn-A1a, Vrn-D1a, Rht-B1b, Rht-D1b, and the high copy number of Ppd-B1.However, anthesis QTL and anthesis related genes of Vrn-B1a and TaFT3-1B on chromosomes 5B and 1B did not lead to frost damage, indicating that these early-flowering phenotype related genes are compatible with frost tolerance and thus can be utilised in breeding.Our results also indicate that wild-type alleles Rht-B1a and Rht-D1a can be used when breeding for frost-tolerant varieties without delaying flowering time.

1.Introduction

Frost can cause significant grain yield and quality losses in wheat crops [1].In spring,when plants sense the gradual increase of temperature and their development proceeds beyond the jointing stage,both winter and spring types show considerable sensitivity to low temperatures (0-12 °C) and frost (<0 °C), particularly to short chilling and frost events at night [2].The early-flowering phenotype in modern wheat and barley cultivars has resulted in significant grain yield and quality losses from frost damage [3].Frost damage during reproductive stage can lead to a multiplicity of symptoms, including dead stems, floret and spikelet abortion,and empty shells along the spikes, thus significantly reducing the seed number per spike.In spring wheat cultivation in Australia,late sowing plants and long season varieties can escape the low temperatures and get less impact,which is called as a frost escaping mechanism.However, those plants most likely will face drought and heat stresses during the grain filling period, resulting in significant yield losses.Farmers and breeders are looking forward to frost tolerant varieties with early maturity.Therefore,short season varieties with non-escaping frost mechanism are desirable in Australia.

In a previous abiotic research, during the reproductive stage,male organs increase the sensitivity dramatically from the start of meiosis to the break-up of the tetrad.A single anther needs approximately 24 h in that event[4].The significant grain number loss was happened during the period of 8-17 days before anthesis when exposing to drought stress [5].Meiosis (10 days before anthesis)is the most sensible stage to abiotic stress[6].At meiosis,male sterility occurs under non-freezing temperatures below 10°C in cereals [7,8].

Many previous wheat frost tolerance studies focused on the vegetative development stages [9-14].Genetic segregation for vegetative frost tolerance (or susceptibility) have been reported on chromosomes 5A, 5B, 5D, and 7B in wheat, and 5H in barley.The frost tolerance QTLFr-A1andFr-B1on chromosomes 5A and 5B were closely linked with the vernalisation genesVrn-A1aandVrn-B1a[11,15].Frost resistance 2(FR2) genes are in control of delayed heading.In combination withVRN1showed reproductive frost tolerance[16],it is not clear yet whether these frost tolerance loci are identical withVrn-A1aandVrn-B1a[17].Similarly, it remains unclear whether the frost locus on chromosome 7B is influenced byVRN3(VrnB4) [9].In barley, three doubled haploid(DH) populations were used to identify reproductive stage frost tolerance QTL[3].Two major QTL were identified on chromosomes 2H and 5H.The QTL for frost-induced floret sterility and grain damage overlapped with the anthesis QTL on theVrn-H1locus in all three DH populations,whereas in two of the populations the floret sterility QTL on 2H was not close to the growth QTL orPpd-H1,but close to theearliness per segene (Eps 2), the cold-regulated gene(Cor14b)and the barley low-temperature gene(Blt14)loci.It seems that the frost tolerance QTL are closely associated with the vernalisation genes in both vegetative and reproductive stages,although the effects kick in at different growth stages.

Regarding the anthesis related genes, three vernalization gene groups includingVRN1,VRN2, andVRN3have been well studied[18].Three homologous copies of theVRN1gene, known asVrn-A1a,Vrn-B1aandVrn-D1awere mapped on chromosome 5A, 5B and 5D, respectively [19-24].TheVRN2gene was mapped on the distal region of 5AL which was a repressor of flowering.WhenVRN2is down-regulated by vernalization, the gene expression ofVRN1was promoted [19,25].The mutated and dysfunctionalVRN2resulted in spring lines [25].VRN3, similar to theFlowering Locus T(FT) gene, encodes a Rapidly Accelerated Fibrosarcoma(RAF) kinase inhibitor-like protein [18].The mutated lines withVrn-B3flowered considerably earlier than the recessivevrn-B3allele [26].One of theVrn-B3genes was mapped on 7BS (VRNB3) [26].The results of Yan et al.[26] implied that theVRN2modulated the quantitative levels ofFT(directly or indirectly) and the absence ofVRN2function showed no effect toVRN1andVRN3mutations.Lately,another gene for developing spring growth habitVRN-D4was identified in the short arm and close to the centromeric region of chromosome 5D [27,28].VRN-D4is a homologous gene ofVrn-A1.

Photoperiod insensitivity provides wheat with the ability to flower in short day as well as in long day conditions.The genes(PPD1)involved in this process arePpd-A1,Ppd-B1,andPpd-D1(formerlyPpd3,Ppd2,andPpd1),located on 2A,2B,and 2D,respectively[29-31].The homologous genesPpd-A1andPpd-B1showed less effect on flowering in short days than thePpd-D1[29].The insensitivity ofPpd-A1is greater thanPpd-B1[32], andPpd-A1tends to increase thousand grain weight and yield whilePpd-B1seems associated with high kernel number through increasing spikelet number[33].In recent study,a heading time QTL detected on 2B in durum wheat and explain 26.2% of the phenotypic variations.The early heading QTL correspond to higher copy number ofPpd-B1[34].

The influence ofearliness per se(eps) genes tend to influence developmental rate at a much lower level as compared to vernalization and photoperiod.Numerousepsand the related flowering-time QTL in wheat have been mapped to chromosomes 1DL, 2B, 3A, 4A, 4B, and 6B [35-37].Epsloci are associated with spikelet number and size, thereby affecting wheat yield [38].An ortholog to theArabidopsis thaliana LUX ARRHYTHMO/PHYTOCLOCK1(LUX/PCL1) gene was identified asEps-3Amin einkorn wheat(Triticum monococcumL.).Lines containingEps-3Amshowed a distorted circadian clock,spikelet number variation and temperature sensitivity [39].In wheat, homologs to ArabidopsisEarly flowering 3(ELF3)gene have been identified on chromosome group 1[40],namely TaELF3-1AL,TaELF3-1BL,and TaELF3-1DL.The gene ELF3 was identified as a candidate gene of Eps-Am1.Wheat lines harbouring ELF3 showed flowered earlier and less spikelets per spike,and stronger photoperiod sensitivity,which indicate the significant epistatic interaction with Ppd1 [41,42].Eps-D1 deletion reduced the total expression of TaELF3 indicating TaELF3-1DL is the major isoform of gene TaELF3 [43].On average, lines harbouring allele TaELF3-1DLb headed two days earlier compared with those holding TaELF3-1DLa [40].

Apart from the well-known anthesis-related genes mentioned above, on 3A, 3B and 3D, there is a set ofT.aestivumGIGANTEA(TaGI) encoding genes, whose products interact with FLAVINBINDING, KELCH REPEAT, and F-BOX 1 (FKF1) domains to form a complex regulating photoperiod-dependent flowering by regulating CONSTANS (CO) expression [44].

ASOC1(Suppressor of Overexpression of CO 1)-like gene on chromosome 4DL,WSOC1, was reported to influence flowering time in wheat[30].A gene for wheat vegetative to reproductive transition on the chromosome 7 group,TaVRT-2, interacts withVRN1andVRN2and regulates the floral transition [45].Three short-day flowering-time genes on 1B, includingflowering locus T3(TaFT3-B1), WUSCHEL-like (TaWUSCHELL-B1) andTARGET OF EAT1

(TaTOE1-B1)have been cloned[46],with the early-flowering function forTaFT3-B1having been validated.A set of heading-date genes (TaHD1) identified on 6A, 6B and 6D, are regulated by long-day conditions and the circadian clock,directly affecting vernalisation genes under long-day conditions.Its mutants showed a delayed flowering response in a long-day environment [47,48].

It is difficult to screen frost tolerance in the field at the reproductive stage, since trials need to be hit by a natural frost event at the right developmental stage, which is purely a matter of chance.To increase the chance, trials consisting of a wide range of genotypes usually need to be planted on multiple sowing dates[3,49].The type of damage to the affected plants will also be influenced by the weather on days leading to the frost event and after.This goes a long way in explaining why untangling the complex genetic basis of frost tolerance under controlled conditions in glasshouses and cold chambers has been of limited utility to breeders.

Our 2018 large-scale field trials encountered such a chance event exactly at the right stage,which provides a valuable resource for QTL detection for frost damage.Six DH populations, namely,Bethlehem/Westonia (BW); Gregory/Bethlehem substitution line 7AS (G7A); Spitfire/Bethlehem substitution line 7AS (Sp7A).Spitfire/Bethlehem(SpB);Spitfire/Mace(SpM);and Suntop/Bethlehem substitution line 3BL(St3B),representing genetically divergent origins were impacted by frost in two distinct environments.All populations suffered considerable frost damage.Genes and markers contributing to the frost tolerance or susceptibility phenotype were identified.Frost tolerance segregation patterns in the current study provide valuable genetic information that can be used in wheat frost tolerance breeding.

2.Materials and methods

2.1.Plant materials

This study was conducted using six DH populations derived from the following crosses:Bethlehem/Westonia (BW); Gregory/Bethlehem substitution line 7AS(G7A); Spitfire/Bethlehem substitution line 7AS (Sp7A).Spitfire/Bethlehem (SpB); Spitfire/Mace(SpM); and Suntop/Bethlehem substitution line 3BL (St3B).Mace,Spitfire,Suntop,Westonia and Wyalkatchem are mid to early flowering varieties and produce good yield (note:Mace is a highyielding cultivar) (Tables S1 and S2).Bethlehem and Bethlehem substitution lines 3BL (B_3B) and 7AS (B_7A) flower early and are considered as frost sensitive wheat lines producing medium yield(Tables S1 and S2).Gregory, Tungsten and Yitpi are lateflowering and are considered frost tolerant varieties and produces medium to high yield in high rainfall seasons (Tables S1 and S2).The population sizes of the six DH populations were 105, 327,304, 168, 222, and 350 for BW, G7A, Sp7A, SpB, SpM, and St3B,respectively (Table S3).The wheat varieties Wyalkatchem, Tungsten and Yitpi(Tables S1 and S2)were included as internal controls in the 2018 field trials.A total of 171 historical lines were used to investigate the proportion of mutatedVrn-A1a,Vrn-B1a,Vrn-D1a,Rht-B1bandRht-D1bgenes (described in Table S4).

2.2.Field and glasshouse experiments

In 2017, individual lines from the six DH populations were planted as 1-m2plots in Katanning (Kat), Wongan Hills (WH) and South Perth, while four DH populations, namely BW, SpM, SpB,and Sp7A were also planted in 4-litre pots in a glasshouse at Murdoch University.

In 2018,the field experiments with the six DH populations were conducted at three locations across Australia representing distinct environments including 1716 plots at Narrabri in New South Wales, 1884 plots at Muresk and 2316 plots at Williams in Western Australia.The majority of the DH lines(greater than 96%)were replicated two times at one or more than one of the three locations.Eight parental lines and three control varieties were also utilised in each of the field experiments.The partially replicated experiments were designed using DiGGer in R [50].Since no frost were recognized in Narrabri, the Zadok data recorded in Narrabri were used to validate the flowering QTL in this study.

To identify the young microspore stage, the parental lines,together with other varieties, were grown in Muresk in 2019.

Field-grown plants were rainfed under standard agronomic management practices, whereas adequate water and fertilizer were provided to plants grown in pots in the glasshouse [51].

2.3.Growth stage measurement

Given the importance of establishing the plant growth stage during the frost events,anthesis date and plant growth stages during ear emergence were recorded for all populations in the glasshouse, in Wongan Hills and Katanning in 2017, and in Williams and Muresk during the 2018 trials (Fig.S1).The days to anthesis were calculated by the anthesis date minuses the sowing date.In the Narrabri trials in 2018, the Zadok stages around 50-60 were recorded.In Williams trials 2018,days to maturity were estimated as follows for each plot except the St3B population:on October 31,based on the plant maturity performance,the plant maturity scoring was entered as 0 (mature) to 15 (based on an estimate of 15 days to maturity).The days to maturity were the estimated maturity date minuses the sowing date (June 4, 2018).

2.4.Temperature recording

In 2018, probably because of high rainfall and low soil temperatures in August, the average anthesis time was delayed by 15 to 20 days in the WA field trials.Flowering time ranged from September 25 to October 12 in Williams,and from September 13 to October 6 in Muresk.The lowest temperatures (-1.1 and 0.3 °C,respectively) were recorded from whether station beside the field trials on September 14, 15 and 16 in Williams, and September 15 and 16 in Muresk (Fig.1a; Table S5).The low temperatures(≤2 °C) lasted 975 and 300 min in Williams and Muresk, respectively.The low temperatures during heading and anthesis caused frost damage in all six DH populations.The days between frost and anthesis were calculated by the average days to anthesis of each line minus the days to frost (September 15, 2018).

2.5.Frost damage measurement

The highest frost damage was to spikes (Fig.1b).On the field,frost damage was scored visually by a single person using a 0-15 scale (0:no frost damage; 1-2:a proportion of 0.5%-5% of spikes damaged by frost within a plot (95%-99% of grain remaining);3-4:5%-20% damaged (80%-95% of grain remaining); 5-6:20%-40% damaged (60%-80% of grain remaining); 7-8:40%-60%damaged (40%-60% of grain remaining); 9-10:60%-80%damaged(20%-40%of grain remaining);11-12:80%-90%damaged(10%-20%of grain remaining);13-15:90%-100%damaged(0-10%of grain remaining) (Table S6; Fig.S2).

To evaluate the reliability of field scoring of frost damage,a set of randomly selected lines were assessed for the floret sterility and the data was compared with the field visual scoring results.For validating the visual scoring,ten heads were picked from each plot around 20 days after anthesis for frost assessment.Three parental lines (Bethlehem, Mace, and Westonia) and two control lines(Wyalkatchem, Yitpi) from three different anthesis time windows were used, and 7 to 14 plots for each line were randomly selected for floret sterility assessment (Fig.S3a).The two outside florets of each spikelet in a spike were used for sterility calculation.The sterility percentage equals to the sterile florets divided by the total florets in each spike,then times 100.The average sterility of the ten heads of each plot formed the whole plot frost damage level.The correlation coefficient between visual frost damage and floret sterility was 0.79(P<0.01)in randomly selected parental and control lines (Fig.S3b).

2.6.Plant height measurement

Plant height was measured using three main tillers of three representable plants in each plot.The plant height of G7A and St3B populations was used for QTL analysis in this study as the semidwarf genes were segregating in these two populations (Fig.S4).

2.7.Meiosis stage measurement

In 2019 field trials in Muresk,plants were sampled in 2-3 days interval since the auricle distance was 1 cm and the main tillers with the same auricle distance as sampled plants were tagged.The sampled plants were stored above ice and the spikes were dissected and photographed with a ruler.The anthers stained by acetocarmine solution (45%) were placed on a light microscopy(400× magnification) and the pollen developmental stages were captured.The days between meiosis and anthesis were calculated by the days to anthesis minus the meiosis date.

Fig.1.Frost occurred 9-18 days before anthesis.(a) Daily temperature below 2 °C presented in Muresk and Williams sites (temperature was recorded every 15 min).In Muresk,the low temperature appeared on the 15th and 16th of September 2018 while it mainly occurred on the 14th,15th, 16th in September, shorter time period on the 11th,13th,17th,18th,25th of August,and 6th of September in Williams.(b)Different types of frost damage occurred in both Williams and Muresk sites.A1,Frost impacted on the lower part of spikes and the top of peduncles;B1 and B2 showed Bethlehem frost impacted spikes and seed settings,respectively;W1 and W2 showed Westonia frost impacted spikes and seed settings, respectively.(c)Data of days to anthesis and days between last frost event and anthesis on each parental line and control variety in two population areas of SpB and BW in Muresk (M) and Williams (W), respectively.Numbers represent days for different statistical classes.The vertical bars represent SE.(d)Number of days between the frost event and anthesis on different parental lines and control varieties.Empty column means the missing line in the assayed area.Values with the same letter are not different at P = 0.05.

2.8.Linkage map construction

Genomic DNA was extracted from a single plant for each DH line and their parents [52].SNP linkage maps were constructed for six DH populations, including two linkage maps of 90K SNP and four linkage maps of 12K Targeted Genotyping-By-Sequencing (tGBS).SNP genotyping was performed using an Infinium iSelect assay on an Illumina iScan instrument according to the manufacturer’s protocols (Illumina, San Diego, CA).SNP clustering and genotype calling was performed using GenomeStudio v2011.1 software (Illumina, San Diego, CA) with the custom genotype-calling algorithm described by Cavanagh group [53].Identical lines were detected and removed using non-metric multidimensional scaling(MDS)of genetic dissimilarity using software from Numerical Taxonomy System (NTsys) v2.2 and Plymouth Routines in Multivariate Ecological Research (PRIMER v6) [54,55].Lines with large proportion of missing values on SNP genotyping were also removed, together with distorted markers and doublecross markers.Most co-segregating markers were made redundant and removed from the genetic map.As results, a fine-map of BW population was constructed using 77 lines and 2387 SNP markers;and so as the maps of G7A by 218 DH lines and 3592 SNP markers,Sp7A by 191 DH lines and 2367 SNP markers, SpB by 94 DH lines and 2570 SNP markers,SpM by 188 DH lines and 2235 SNP markers, and St3B by 185 DH lines and 1924 SNP markers.(Table S3).

Gene-based markers were used in conjunction with SNPs to account for the fact that in spring wheat, anthesis is highly controlled by the vernalisation genes [56] and seed number per spike is strongly associated with the wild-type ofRht-D1agene [51].Markers forRht-B1b,Rht-D1b,Vrn-A1a,Vrn-B1a, andVrn-D1awere used to construct the final maps of the six DH populations using Map Manager [57] and the QTL mapping package R/qtl [58].Primers forVRN1(Vrn-A1a,Vrn-B1a, andVrn-D1a), andRht-B1,Rht-D1were as described [51,56].Primers forTaELF3-1A,TaELF3-1BandTaELF3-1D,VRN-D4genes were used for PCR amplifications in all parental lines.The gene-specific markers were adopted from Alvarez et al.[41], Wang et al.[40], and Kippes et al.[28].

2.9.QTL mapping

As proposed for additive, dominant and epistatic QTL mapping in biparental populations [59], we employed inclusive composite interval mapping (ICIM) for QTL detection, using IciMapping V4.1(http://www.isbreeding.net).A LOD score of 2.5 was used as significance of QTL detection.Permutations were set to 1000 at a significant level of 0.05.The inclusive composite interval mapping addition (ICIM-ADD) method was selected for QTL mapping [60].

2.10.TaqMan assays for Ppd1-2B copy number determination

Significant anthesis and frost QTL were detected in thePpd1-2Bregion.As variable copy number of thePpd1-2Bgene among parental lines was suspected, TaqMan assays were conducted to determine copy number in those lines.Based on a published protocol[61],20 μL PCR reactions were set up including 10 μL ddPCR Supermix for Probes (Bio-Rad), 0.4 μL probe plus forward and reverse primers (10 μmol L-1), 5 μL DNA (10 ng μL-1), and 4.6 μL RNase/DNase-free water.The primer and probe sequences forPpd1-2Bgene andTaCO2internal control were as those published by Díaz et al.[61].PCR cycling parameters were 95°C for 15 min;40 cycles of 95°C for 15 s;and 60°C for 60 s.Ppd1-2Bcopy number was analysed based on the ratios of absolute copy numbers against theTaCO2control.

2.11.Statistical analysis

Linear mixed models (LMM) were fitted with ASReml-R [62] in the analysis of the frost and growth stage phenotypic traits,where the variance parameters in the mixed model are estimated using the residual maximum likelihood(REML) procedure[63].Residual diagnostics were performed to examine the validity of the model assumption(normality and homogeneity of variance).The best linear unbiased predictions (BLUPs) were used for the phenotypic traits.

Phenotypic data were analyzed by multivariate analysis of variance (MANOVA) using the general linear model implemented in IBM SPSS statistics 24 (https://www.ibm.com/au-en/products/spss-statistics).Wilks Lambda was used as the multivariate test statistic.Post-hoc Tukey’s Multiple Range tests were used to identify significant groupings.Pearson correlations of the parameters investigated were analysed by SPSS software using the BLUP values across environments.The broad-sense heritability was calculated through R-studio.

3.Results

3.1.Frost impact levels associated with growth stages in DH parental and control lines

In the 2018 field trials sown at Williams and Muresk(WA),low temperatures during heading and anthesis caused frost damage to all six DH populations(Fig.1a).The most visible frost damage was on spikes (Fig.1b).Probably because of the different developmental stages between main spikes and tiller spikes, the levels of frost damage to spikes of the same wheat variety in the whole plot show large variation, such as ‘‘white spikelets”- spikelets appear white and dead; ‘‘half-cut” and ‘‘bald-pointed” spikes- only the lower part of the ear has grains;‘‘toothed ears”-grains set in several spikelets.Since flowering times were different not only between field trial locations but also between sowing areas within the same location, the representative anthesis data for the parental and control lines in the areas sown to Spitfire/Bethlehem (SpB) and Bethlehem/Westonia (BW) are presented in Table S2 and Fig.1c, d.On average, anthesis times for cultivars Gregory, Tungsten and Yitpi were significantly later than those for the other parental cultivars at both locations.The cultivar Suntop flowered two and five days later than the rest of cultivars in Williams and Muresk, respectively.The cultivar Bethlehem and its substitution lines B_3B and B_7A were the earliest flowering lines, while Westonia, Spitfire and Mace flowered two to three days later.The cultivar Wyalkatchem was one and two days later than Mace and Westonia.In Muresk, the anthesis time of Suntop showed large differences(three days) between the SpB and BW sown areas (119 and 122),which may lead to significantly different levels of frost damage.

In Williams,the earlier-flowering lines of Bethlehem, B_3B and B_7A,were impacted the most by frost whereas the late-flowering lines Gregory, Tungsten and Yitpi were not affected (Fig.2a, b;Table S7).Among the cultivars with a similar flowering window,such as Spitfire, Mace and Wyalkatchem, Westonia was the most susceptible to frost,whereas Mace showed the highest level of tolerance.The same trend was observed in Muresk.Spitfire displayed frost tolerance in Williams but was susceptible in Muresk.Its anthesis window was very close to that of Mace in the SpB sown area and to Wyalkatchem in the BW area.Suntop showed high frost tolerance in Williams but was susceptible in Muresk, where its flowering window was shorter.Floret sterility also showed the same levels of frost impact in selected parental and control lines belonging to the three anthesis windows (Fig.S3a, b).

According to the date of the frost event (September 14, 15 and 16, 2018), the time of impact was between 9 and 18 days before the average anthesis time.The 2019 Muresk field trial data for Suntop at 16 days before anthesis (DBA) are presented in Fig.S3c-g.The data show that wheat plants were at the flower development and anther differentiation stages at 16 DBA,according to Bonnett’s description[64],and some flowers were at meiosis stage.Based on Koonjul et al.[4],lines are most vulnerable to abiotic stress during spikelet development and meiosis stages.The 2018 frost events incidentally hit the most susceptible meiosis stages of our trial.

3.2.Phenotype differentiations across DH populations and environments

Different frost impact levels were observed in six DH populations.On average, the frost impact levels in BW (8.0), Sp7A (8.0),SpB (8.1) and St3B (6.2) were higher than that in G7A (5.0) and SpM (3.7) populations.The impact levels in Muresk (8.6) were almost doubled in Williams (4.4) (Table S8).The coefficients of variation were significantly higher in Williams compared with Muresk.The high heritability in frost impact was observed (0.71-0.84) in all six populations across two environments.Apart from the frost impact, plant height showed the highest heritability(0.96), followed by Zadok stages (0.79-0.95), anthesis (0.52-0.82)and maturity (0.4-0.87).The significantly positive correlation levels of frost impact (P< 0.01) between replicates and environments further demonstrated that the traits were under genetic control (Table S9).

3.3.Rht,VRN and related anthesis gene segregations in DH populations

Fig.2.Sensitivities of parental lines to frost impact associated with parental line maturity in different population areas in Muresk and Williams.(a)Sensitivities in G7A,Sp7A and SpB area in both Muresk(M)and Williams(W).(b)Sensitivities in SpM,BW and St3B area in both Muresk(M)and Williams(W).The vertical bars represent SE;values with the same letter are not different at P=0.05;*represents no replicate;empty column means the missing line in the assayed area.(c)Summarized table for parental lines with or without VRN1, Rht, ELF3-1 Da and VRN-D4 genes; Ppd-B1(Ppd1-2B) copy number ratio and the gene segregated populations.

A screen of a suite of phenology genes in the DH lines identified 18 haplotypes segregating forVrn-A1a, Vrn-B1a, Vrn-D1a,Rht-B1b,Rht-D1b,ELF-1 Da,VRN-D4andPpd-B1(Figs.2c,S5a-e).In DH populations, SpB and Spitfire × B_7A (Sp7A) segregated forVrn-A1a,whileVrn-D1asegregated in Spitfire/Mace (SpM) and BW, andVrn-B1ain the BW,SpB,Sp7A and St3B populations.ForRhtgenes,only Suntop and Gregory carry both theRht-B1band the wild-typeRht-D1aalleles.The other six parental lines have the opposite genotype for these alleles.Therefore,Rht-B1bandRht-D1bsegregated in two populations of Suntop × B_3B (St3B) and Gregory × B_7A (G7A).No amplification was detected by primer markers of geneTaELF3-1AandTaELF3-1B.Fragments of 709 bp amplified by primers of TaELF-1DL-F5/R4 presented in Gregory,Mace, Suntop, Westonia and Chinese Spring which hold alleleTaELF3-1DLa(delaying heading two days) whereas no amplification was in Bethlehem, B_3B, B_7A and Spitfire (Fig.S5d).The results indicate theTaELF3-1Dsegregations in BW, SpM, G7A and St3B population (Fig.2c).Fragments of 1.4 kb and 1.28 kb were amplified by primer pairs of VRND4-ins.F4/R3 and VRND4-ins2.F1/R1 in both upstream and downstream ofVRND4insertion,respectively, and showed that Bethlehem, B_3B, B_7A, Gregory,Mace and Spitfire harbourVRN-D4while the Suntop and Westonia were negative(Fig.S5d),indicating theVRN-D4segregations in the populations of St3B and BW(Fig.2c).Copy number differences for thePpd-B1allele were detected in parental lines and the segregations in SpB, Sp7A, G7A, and St3B populations were expected(Fig.S5e).

3.4.Cross population verified frost QTL and related potential candidate genes

Thirty individual frost QTL were detected in the six DH populations.Out of these, there were 18 major QTL that were each responsible for a phenotypic variation greater than 9.5% and were distributed across 13 chromosomes(2A,2B,2D,3A,4A,4B,4D,5A,5D, 6D, 7A, 7B, and 7D).Most frost QTL were closely linked to the QTL for anthesis and physiological maturity Zadok stages and to anthesis-related genes.The QTL results for individual DH populations are included in Table S10.

A large proportion(83%)of the detected frost QTL were showed up consistently in two or more DH populations,except for the QTL on SpM_1A, SpM_1D, SpM_3D, Sp7A_6D, and G7A_3B (Fig.S6).A frost QTL was detected on the homologous region on the short arm of chromosome 2A in the SpB and SpM populations, while the other two frost QTL distant from each other in SpM were closely linked to the anthesis and physiological maturity loci on the long arm of chromosome 2A (Table S10; Fig.S7).The photoperiod(Ppd)allelePpd-A1is located on the short arm of 2A(Table S11).In the SpB and SpM populations, thePpd1-2A(2A:36.9 Mb) allele is located between 2AM5988 (2A:31.7 Mb) and 2AM67517 (2AL 40.8 Mb) and tightly linked to the frost QTL on 2A.The frost QTL on the 2A long arm in the SpM population is a new locus segregating for frost.

A highly significant frost QTL (the highest LOD score 9.2) was detected on homologous regions of the 2B chromosome in four populations of Sp7A, SpB, G7A and St3B, and the additive effects were all contributed by Bethlehem and its substitution lines,while another QTL (LOD score 3.9) in SpM was also closely associated with the homologous region, with the phenotype contributed by Mace (Table S10; Fig.3).APpd1-2Bgene (gene bank number DQ885765) was used to identify the physical map location.The potential location ofPpd1-2Bwas on the short arm of 2B at 56.2 Mb, which was close the marker 2BM46469 (2B:58.8 Mb)in Sp7A, SpB and G7A, and between markers of 2BM29812 (2B:53.4 Mb) and 2BM73250 (2B:77.2 Mb) in St3B (Table S11).Strikingly, the location ofPpd1-2Bwas in the frost and anthesis QTL region in those four populations, and the QTL results matched the potential segregation ofPpd-B1gene copy number (Fig.2c;Fig.S5e),which points toward a contribution ofPpd-B1copy number to the phenotype.This QTL corresponds to the newly identifiedPpd1-2Blocation in this study.

On linkage groups 2D3 and 2D2,frost QTL were detected in both SpM and Sp7A, contributed by Mace and Spitfire, respectively.These QTL were associated with anthesis time, as anthesis QTL were closely linked to the homologous region of 2D(Fig.S8).GenePpd1-2Dwas located at 2D:33.9 Mb from the upstream on the short arm of 2D.SNP M22365 (2D:18.2 Mb) on Sp7A_2D1 and SNP M15795 (2D:17.3 Mb) on SpM_2D1 were close toPpd1-2Don the physical map.The frost and anthesis QTL on SpM_2D3 and Sp7A_2D2 tend not to be associated with thePpd1-2Dlocus,since thePpd-D1alocus does not expect to segregate in those populations(Fig.S5b).This frost QTL may imply a new phenology gene on 2D, associated with anthesis time.

On chromosome 3A,the frost QTL in SpB population was located on the short arm, whereas in St3B it was located on the long arm.The QTL on the short arm overlapped with anthesis QTL (Fig.S9).TaGIgene sequences (AF543844) were isolated [65] and were located on 3A (3A:84.1 Mb), 3B (3B:117.9 Mb), and 3D (3D:71.9 Mb).One frost QTL in SpB_3A was located between SNP markers 3AM63098 (3A:56.4 Mb) and 3AM34922 (3A:107.7 Mb)(Table S11).TaGI-3Ais most likely the candidate gene contributing to this frost QTL and the associated anthesis QTL.The 3A long arm QTL in the St3B population was located between 3AM9626 (3A:711.0 Mb) and 3AM39004 (3A:729.7 Mb), which are not associated with theTaGIgene.The Genbank number of KF769443 was used to identify the location ofEps-3Amon the physical map.Eps-3Am(3A:740.1 Mb) was located between the markers 3AM10770(3A:739.3 Mb) and 3A73079 (3A:741.2 Mb), about 10 cM away from the frost QTL (Table S11).

Significant frost QTL were consistently detected on the long arm of chromosome 4A in the SpM, G7A, and St3B populations(Fig.S10).Interestingly, the phenotypes were all attributable to the male parents of Mace, B_7A and B_3B.The frost QTL in SpM and G7A populations were located on the anthesis QTL region,and the late anthesis QTL were contributed by the female parents,Spitfire and Gregory.In other words, the early-flowering phenotype contributed by male parents led to frost susceptibility.The homologous gene of the DELLA protein (rht1-D1a; AJ242531) on 4A is TraesCS4A02G271000 (4A:582.4 Mb), which is located between markers 4AM76744 (4A:575.0 Mb) and 4AM43375 (4A:597.6 Mb), and above the common marker 4AM77169(Table S11).Interestingly, the newly identifiedrht1-4Alocus in the current study was located 20 cM away from the minor frost QTL regions on the short arm of 4A in SpM and G7A populations.The homologous gene sequence ofWSOC1-4Dis TraesCS4A02G320300, on chromosome 4A:608.8 Mb, near SNP marker 4AM77381(4A:612.1 Mb),which is 30-40 cM away from the frost QTL on the long arm regions of 4A in the G7A and St3B populations(Table S11).The potentialWSOC1-4Alocus is further away from the frost QTL in the SpM population.

Significant frost QTL were detected on theRht-B1bregion on chromosome 4B in G7A and St3B populations segregating forRht-B1b(Fig.4a;Table S11),which were contributed by the female parents, Gregory and Suntop, and both harbour theRht-B1ballele.No flowering QTL were detected in those regions.The plant height QTL with the LOD values over 14 were overlapped with the frost QTL in the semi-dwarf genesRht-B1bregion and were contributed by male parents B_7A and B_3B containing the wild-typeRht-B1aon 4B.Wide-type ofRht-B1aattributed the average height phenotype of 20%-30% and the semi-dwarf geneRht-B1bcontributed to frost damage QTL on 4B in G7A and St3B populations (Table S10;Fig.4a).

Fig.3.Frost QTL were overlapped with anthesis QTL on the homologous region on 2B in five populations of Sp7A, SpB,G7A, St3B and SpM.The location of Ppd1-2B in High Confidence (HC) 1.0 physical map was indicated.

Fig.4.Significant frost QTL were detected on the Rht1(Rht-B1b)and Rht2(Rht-D1b)regions.(a)Frost and plant height QTL on Rht1(Rht-B1b)region in the populations of G7A and St3B while the frost QTL in BW and SpM were on the distal downstream and upstream, respectively.A new gene of WSOC1-4B may close to the frost QTL on 4B in BW population.(b) Significant frost and plant height QTL were detected on the Rht2 (Rht-D1b) region in the populations of G7A and St3B.

A significant frost QTL with a LOD score of 16.0, attributable to Mace, was repeatedly detected on the distal region of chromosome 4B short arm in the SpM population.The QTL region was about 80 cM away from the recessiveRht-B1agene loci.The contributing gene for this QTL is unknown.A minor frost QTL was detected in the BW and G7A populations on the homologous region of the long arm of chromosome 4B.In the BW population,the frost QTL were closely linked to the anthesis QTL.Interestingly, frost QTL in the BW population was in the 4BM4887(4B:646.6 Mb) region, which is only 6.4 Mb away from theWSOC1-4Band the anthesis QTL, close to 4BM8859 (4B:660.7 Mb) (Table S11).These results suggest thatWSOC1-4Bmight be a new gene influencing flowering time and contributing to segregation for frost.

Likewise, significant frost QTL were detected on theRht-D1bgene region on 4D in the G7A and St3B populations segregating for theRht-D1bgene, contributed by B_7A and B_3B, which harbour theRht-D1ballele (Fig.4b; Table S11).The plant height QTL with high LOD value over 14 were located onRht-D1bregion and contributed 20%-44%of the total phenotype by Gregory and Suntop containingRht-D1aon 4D.The wide-type ofRht-D1acontributed to plant height whereas the semi-dwarfRht-D1battributed to the frost damage QTL on 4D in G7A and St3B populations (Table S10; Fig.4b).

Frost QTL were consistently detected on the homologous region of chromosome 5A in SpM, Sp7A and SpB populations, with the highest LOD score of 7.7 found in Sp7A(Fig.5).The frost QTL were in theVrn-A1aregion, attributed to B_7A and Bethlehem in the Sp7A and SpB populations (Table S11), respectively, whereas it was contributed by Spitfire in the SpM population.The frost QTL in the Sp7A population overlapped with the anthesis QTL contributed by Spitfire.Both B_7A and Bethlehem host aVrn-A1amutant and flower early.The results further indicate that frost damage was caused by the early flowering genotypes of B_7A and Bethlehem.A minor frost QTL in the St3B population was on the distal region of the short arm 5A, away from the anthesis QTL, and not associated with the flowering genes or their related phenotypes.

In the case of chromosome 5D,frost QTL were detected in SpM,BW and St3B populations(Fig.6).The largest frost QTL(the highest LOD score 14.2) on 5D was in theVrn-D1aregion in the SpM population while in the BW population it was close toVrn-D1a(Table S11).The frost phenotype in the SpM and BW populations was contributed by Spitfire and Bethlehem, which harbourVrn-D1a,whereas for the St3B population in the Williams trial,the contributed phenotype was attributed to Suntop.Several anthesis and maturity QTL were detected in theVrn-D1aregion, indicating that frost damage is closely associated with plant growth stages.Another vernalisation related gene,VRN-D4,which is a homologue ofVrn-A1a, originated from a large segment of chromosome 5A inserted into the short arm of 5D[27,28].According to the physical map ofVRN-D4(5D:193 Mb),the gene location should be close to marker 5DM45442 (5D:278 Mb) on the St3B _5D map, while it is next to 5DM62708 (5D:154 Mb) on the BW map (Table S11).The frost QTL may stem from theVRN-D4locus as theVRN-D4was segregating between the parental lines of Bethlehem and Westonia,Suntop and B_3B, respectively.

Significant frost QTL were detected on chromosome 7A in SpM,SpB and BW populations(Fig.7),with the highest LOD score of 10.The phenotype in the SpB and BW populations was contributed by Bethlehem.Frost QTL in the SpM population were located in three positions, mainly attributing to Mace.According to the physical map ofTaVRT-2-7A(7A:128.8 Mb), the closest SNP marker is 7AM1849 (7A:128.5 Mb), close to 7AM75587, which was tightly linked to the frost QTL in the SpB and BW populations, while in the SpM population the location was closely linked to the anthesis QTL on 7A(Table S11).One frost QTL was detected in the SpM and SpB populations on the homologous distal region on chromosome 7A.This corresponds to a new segregating locus for frost impact on 7A, which was not associated with the flowering genes or their related phenotypes.

For chromosome 7B, frost QTL were detected in the SpM and BW populations (Fig.S11a).However, the frost QTL regions were distinct from each other.The frost QTL in the SpM population,contributed by Mace, was on the terminal region of the short arm of chromosome 7B, away from the anthesis QTL, whereas in the BW population,the frost QTL were closely linked to the QTL for anthesis and maturity, and were contributed by Westonia.The anthesis QTL in the SpM and BW populations were located on the homologous region on 7B.Another anthesis QTL detected exclusively in the SpM population was located on the distal region of the long arm of 7B.In the SpM population, the physical map ofVRN-B3(7B:9.7 Mb) was close to marker 7BM76084 (7B:6.7 Mb), which may contribute the frost phenotype in SpM (Table S11), while the location was 30 cM away from the frost QTL in the BW population.Based on the physical map,TaVRT-2-7B(7B:90.1 Mb) is located between 7BM53206 (7B:64.7 Mb) and 7BM10089(7B:115.2 Mb) on the BW map (Table S11), closely linked to the QTL for anthesis and frost,suggesting the involvement of this locus in frost damage in the BW population.In the SpM population,TaVRT-2-7Bis not close to the frost QTL.

A QTL with a LOD score of 4.6 also appeared on the homologous regions on chromosome 7D in the Sp7A and SpB populations(Fig.S11b).The frost phenotype was contributed by B_7A and Bethlehem and was highly associated with anthesis QTL.TaVRT-2-7D(7D:128.9 Mb) was close to the upper SNP marker 7DM76171 (7D:112.0 Mb) in Sp7A population and distant to the lower SNP marker 7DM42766 (7D:182.6 Mb) in both the Sp7A and SpB populations (Table S11).This gene tends to contribute the significant frost QTL on 7D.

3.5.Anthesis QTL not segregating for frost tolerance

No frost damage QTL were detected on chromosome 1B.However, anthesis and maturity QTL were present on the central chromosome regions in the G7A and SpM populations and on the distal region of 1B short arm in the SpB population (Fig.8).GeneTaFT3-1B(1B:581.4 Mb) was tightly linked to SNP marker 1BM77588(1B:581.6 Mb), which is next to the anthesis QTL region(1BM54518) in the G7A population (Fig.8; Table S11).In the SpM population, theTaFT3-1Blocus was possibly located above 1BM73284 (1B:662.9 Mb), on the anthesis QTL region(Table S11).On the physical map,TaELF3-1Bwas located on the distal region of the long arm in these three populations.Two other anthesis-related genes,TaWUSCHELL-1B(1B:53.3 Mb)andTaTOE1-1B(1B:59.1 Mb), were close to markers 1BM42781 (1B:56.8 Mb)and 1BM47932 (1B:58.7 Mb) on the Chinese Spring physical map(Table S11).These two genes were located 25 cM away from the maturity QTL in the SpB population.It is not clear whether the QTL was influenced by them.

Another chromosome without frost damage QTL was 5B, even though anthesis QTL are present in the St3B, SpB, Sp7A and SpM populations, andVrn-B1awas segregating in the St3B, SpB and Sp7A populations (Fig.9).The anthesis QTL was on or closely linked to theVrn-B1aallele in the St3B and SpB populations(Table S11), and to the distal regions of either the long arm in the Sp7A population or the short arm in the SpM population.

3.6.Allele frequency of VRN1 and Rht genes in historical lines

Fig.5.Significant frost QTL were overlapped with anthesis QTL where Vrn-A1a(VPA)tightly linked in Sp7A and SpB populations.The frost QTL on SpM_5A was about 25 cM below the Vrn-A1a and 90 cM above VRN2 loci while it was on the upstream on St3B_5A.

According to the frost trial results, except forVrn-B1a, the mutatedVrn-A1a,Vrn-D1aandRhtalleles were highly associated with the frost QTL and contributed to the frost damage phenotype.Since global warming accelerates terminal drought severity, early flowering genes with early maturity phenotype have been gradually introduced into modern varieties through wheat breeding for avoiding late dry season.In the meantime, the proportions ofRhtgenes were also increased.The utilization of these gene may increase the frost sensitivity in modern varieties.Therefore, it is interesting to see the allele frequency of these mutated genes in historical varieties.In the current study, a set of 171 commercial wheat cultivars, released between 1890 and 2015, were surveyed for theVRN1andRhtallele frequencies.ForVRN1andRht,recessive alleles made up 95%-100%of varieties released before 1967,except forvrn-B1(15%) (Fig.S12).After 1968,Vrn-A1awas maintained at around 61% whereasVrn-B1adropped from 85% to 50% and remained at that level thereafter.TheVrn-D1aallele appeared mainly in the varieties released after 2000.TheRht-B1ballele was presented in a large proportion of varieties released after 1968 (56%-69% of lines) while theRht-D1ballele gradually increased to near 41% in recent years.

Fig.6.Significant frost QTL were overlapped with anthesis QTL in SpM where Vrn-D1a(VrnD)tightly linked in SpM.The Frost QTL were linked to Vrn-D4 gene in BW and St3B populations on chromosome 5D.

4.Discussion

4.1.The most sensitive growth stage to frost

In previous drought studies in wheat, the young microspore stage of pollen development, before anthesis, appears to be the most sensitive to mild water stress[4,66].It has been reported that in rice cold-induced pollen sterility at the young microspore stage had effects comparable to those of drought stress[2].Cold-tolerant lines at the microspore stage are also tolerant to drought stress.The same mechanism was observed in sorghum [67,68].

The young microspore stage in wheat is the time when the auricle distance(AD)between flag leaf and penultimate leaf is 5-8 cm,around 10 days before anthesis [5,66].In our 2019 wheat trials at the Muresk site,the young microspore stage of Suntop was 16 days before anthesis and the AD was 6 cm.In 2018, the frost events at the Williams and Muresk locations occurred 9-18 days before the average anthesis time, which was the most vulnerable growth stage for wheat to endure low temperatures.Although the lowest temperature in Muresk remained above 0 °C, plants still suffered from the sudden temperature drop.According to previous studies[69,70], cold stress induces ABA accumulation in rice anthers,which represses anther cell wall invertase activity.This in turn hampers sugar transport from the tapetum to the pollen and pollen sterility occurs.It is interesting that cold and drought stresses share the same pathways in inducing pollen sterility.

4.2.Non-escaping mechanism related frost tolerance

Since the reproductive stage is the period most sensitive to frost, flowering time becomes very critical.Ideally, late flowering is helpful in avoiding frost events, which happens mostly during early spring.This is,however,contradictory to the needs of avoiding terminal drought, which requires early maturity [71].Being able to combine early flowering with reproductive frost tolerance is thus highly desirable for wheat breeding.In the previous study,frost QTL were detected in the same region as theVrn-A1aallele on chromosome 5A [10,72].It has been reported that plants withVRN1copy showed normal flowering but reducing frost tolerance[73].In our study, significant frost QTL were detected close to theVrn-A1aallele in both the Sp7A and SpB populations at the Williams location.In the SpM population grown in Muresk,significant frost QTL were detected close to theVrn-D1aallele,where the QTL for anthesis and maturity were located for all three environments,which further validated the tight connections between frost damage and the early anthesis allelesVrn-A1aandVrn-D1a.

Fig.7.Significant frost QTL were linked to TaVRT-2-7A gene in SpB and BW populations.One of the frost QTL in SpM were linked anthesis QTL on 7A while TaVRT-2-7A closely linked two anthesis QTL.

Fig.8.Anthesis QTL on 1B in G7A,SpB and SpM populations.The locations of published anthsis related genes of TaFT3-1B,TaWUSCHELL-1B and TaTOE1-1B in physical map in Chinese Spring and closest SNP markers are indicated.

Fig.9.Anthesis QTL detected on 5B in St3B, SpB, Sp7A and SpM populations. Vrn-B1a (VrnB) was segregating in St3B, SpB and Sp7A populations.

On the other hand,theVrn-B1aallele segregated in the BW,SpB,Sp7A,and St3B populations,with anthesis QTL closely linked to theVrn-B1agene in the SpB and St3B populations.Strikingly, no frost damage QTL was detected on chromosome 5B in any of the six populations studied at the two locations.The link between frost and anthesis was not observed on 5B.In our previous study on days to anthesis, the contribution of winter allelevrn-B1was less than that ofvrn-D1andvrn-A1[56].In the historical lines studied, the frequency of spring allelesVrn-A1aandVrn-D1ain recent varieties amounted to 0.62 and 0.40, respectively, while the frequency ofVrn-B1adiminished from 0.85 to 0.52 over time.An increased proportion ofVrn-A1aandVrn-D1amay increase the risk of frost damage.The most popular variety, Mace, only possesses theVrn-B1aallele and shows good frost tolerance.This clearly demonstrates that the early flowering phenotype induced byVrn-B1ais associated with higher frost tolerance.In wheat, VRN1 determines the most natural variation in flowering.The functional VRN1 proteins are responsible for the apical meristem transition from the vegetative to reproductive phase.However, VRN1 protein are not essential for wheat flowering [74].Although the underlying mechanism is unknown, it can be speculated that the metabolites associated with theVrn-B1afunctional network may contribute to the increased tolerance [75].Because of the global warming, the earlier flower lines without frost damage are highly demanded[76].Our results revealed that theVrn-B1aallele can be utilised in breeding for frost tolerance.

Another potential anthesis gene associated with the nonescaping frost mechanism isTaFT3-B1on 1B.GeneTaFT3-1Bwas very tightly linked to the anthesis QTL on 1B in G7A and SpM populations.The study of Zikhali et al.[46]showed thatTaFT3-B1deletion lines was associated with late flowering while increased the gene copy number was related to early flowering.TaFT3-B1gene was suggested to promote flowering.Two maturity QTL were detected on 1B in SpB populations, where two potential anthesis genesTaWUSCHELL-B1andTaTOE1-B1were closely linked to.In Maize,TARGET OF EAT1was also named asTaSRR1-B1,TaWUSCHELL-B1andTaTOE1-B1, respectively.Similar toTaFT3-B1,TaTOE1-B1had higher gene expression during short days and showed early flowering function [46].Nevertheless, no frost QTL were detected on 1B in those three populations in two locations.It implies that the markers or the copy numbers for theTaFT3-1Bare able to be used for the purposes of shortening the days to anthesis and reducing frost damage risk.

4.3.Recessive Rht-B1a and Rht-D1a genes associated with frost tolerance

Significant frost QTL were detected in the proximity of theRht-B1bandRht-D1bgene regions in the G7A and St3B populations, in which the two genes were segregating.In both populations, the frost QTL were contributed by the mutated type ofRht-B1bandRht-D1b,the semi-dwarf alleles that are associated with short plant height as the plant height QTL were contributed by the wide-type ofRht-B1aandRht-D1a(recessiverhtgenes).This indicates that the mutated or short plant genotypes are more prone to frost.On the other hand,the wild typesRht-B1aandRht-D1aare relatively frost tolerant.

Wheat life cycle can be divided into two phases:one is the stem elongation from jointing to anthesis and another is the grain filling stage from anthesis to maturity.The former stage is more important for yield as the number of fertile florets at anthesis will be determined during this stage [77].The pre-anthesis phase may be more sensitive to photoperiod or temperature, and different levels of hormones(gibberellin, auxin or cytokinin) [78-81].Since several dwarfing genes are associated with GA biosynthesis or signalling,the possible effect on spike development may exist[82,83].It is well known that lines withRht-B1bandRht-D1balleles have reduced sensitivity to gibberellic acid (GA), and GA-responsive growth is therefore repressed [84].This may affect GA levels among the different combinations of dwarfing genes,and thus lead to varying flowering times.It has been reported that GA can increase the transcription level ofSOC1including a MADS-box gene, promoting flowering [47].In a study of ethylene effects in the GA-GID1-DELLA signalling pathway,ethylene leads to a reduction of GA, which delays floral induction[84].Nevertheless, previous research has shown that GA-insensitive dwarfing genes have no effect on spikelet primordia at the shoot apex nor on the number of leaves and internodes [85].In the current study, no flowering QTL were detected in theRht-B1bandRht-D1bgene regions.This suggests that the regulatory function of theRht-B1bandRht-D1bgenes on anthesis is minor.On the other hand, the variations in plant canopy architecture linked to dwarfing genes,for example,height, tillering, ear structure etc.may affect the frost tolerant levels.Selecting taller or wild-typeRht-B1aandRht-D1aalleles for frost tolerance in breeding programs will not cause delayed flowering time or maturity,thus theseRht-B1aandRht-D1aalleles can be recommended for frost tolerance breeding.

In previous studies, the semi-dwarfing allelesRht-B1bandRht-D1bdid not lead to higher grain yields in drought environments[51,86].The QTL analysis showed that the recessiveRht-B1aandRht-D1aalleles are associated with higher grain yields while the semi-dwarfing alleles are associated with lower grain yields [51].Our results further demonstrated the disadvantages of theRht-B1bandRht-D1balleles in frost-prone environments.From 1967 to 2015, the frequencies of theRht-B1bandRht-D1balleles have increased from 0 to 0.6 and 0.4 in Australian wheat varieties,respectively.This is in line with the first green revolution in that the dwarfing genes are utilised widely.However, a recent study has revealed that the short plant phenotype has become a bottleneck for wheat grain yield improvement [51,86].Achieving high grain yield by increasing plant height has become a newly established breeding approach.Our study indicates that eliminating the dwarfing allelesRht-B1bandRht-D1bin breeding can increase frost tolerance in wheat.

4.4.High copy number of Ppd1-2B increase the sensitivity to frost

Day length (photoperiod) affects plant flowering time.QTL for frost and the time of anthesis and maturity were repeatedly detected on chromosomes 2A, 2B and 2D.ThePpd1-2Aallele was tightly linked to the frost QTL on 2A.Remarkably, four out of five frost QTL (next to anthesis QTL) were detected in the proximity ofPpd1-2Bin five populations.Due to thePpd1-2Bcopy number differences between parental lines of these four DH populations,it is likely that the significant frost QTL were contributed by the segregation ofPpd1-2Bgene copy numbers,with a high copy number ofPpd1-2Bshortening the days to anthesis and inducing early flowering[34],thus contributing to the frost QTL.Ppd1-2Ddid not segregate in our six populations and no strong linkage to the frost QTL on 2D was detected.

4.5.Escaping mechanism related frost tolerance

Besides theVrn-A1aandVrn-D1a, the following early flowering related genes are associated with the frost damage QTL and the recessive flowering genes (flower late genes) correspond to frost tolerance, which led a frost escaping mechanism.GeneVRN-B3was mapped to the short arm of chromosome 7B (7B:9.7 Mb)[26].In the SpM population,VRN-B3is tightly linked to a 7B short arm frost QTL that mapped on the top region (7BM55522) of the chromosome.This implies thatVRN-B3contributes to the frost phenotype in the SpM population.In addition, theTaVRTgenes on group 7 appeared to be tightly linked to the QTL of anthesis and frost.TheTaVRTgene is independently modulated by photoperiod and vernalization.It can inhibitVRN1activity through binding to the CArG box of theVRN1promoter in vivo.After vernalization,bothTaVRT2andVRN2functions are repressed,which causesVRN1accumulation for the transition from vegetative to reproductive phase [87].The functions ofTaVRTgenes seem to be critical in the SpM, Sp7A, SpB, and BW populations.

Numbers of genes in photoperiod pathway also regulate anthesis.TaGIlocated on 3A, 3B and 3D, regulate CO which then mediates photoperiodic flowering[44,65].In previous studies,TaGIwas found to be initiated by photoperiod and then expressed in both vegetative and reproductive tissues, suggesting a function in anthesis regulation [44,65].Coincidentally, in our study, the frost and anthesis QTL, and maturity QTL were localized on theTaGIregions on 3A and 3B in SpB and G7A populations, respectively.This indicates thatTaGIfunctions are both in anthesis and frost.

Apart from the vernalisation and photoperiod genes,epsgenes are the third factor regulating heading and anthesis.Oneepslike QTL identified on the long arm of 1DL, namely,TaELF3, held the function to theT.monococcum Eps-Am1locus [37,43].The homologs ofTaELF3were also identified in 1AL and 1DL [40].One significant frost damage QTL contributed by Spitfire holdingTaELF3-1DLb(heading two days earlier) on the long arm of 1D in SpM population were tightly link to the homologous genes ofTaELF3on 1D.It indicates theepsgene functions in heading and contributes to the frost impact in the location of 1D.

Heading related geneTaHD1also designated as CONSTANS2(CO2), is the homolog of CO in wheat, and the competitor ofVRN2[47,48].TaHD1are identified on 6A,6B and 6D.In our study,one of the frost QTL together with the anthesis QTL were on the region of TaHD1-6D physical map location, implying that the TaHD1 gene contributes to the phenotype.

In summary, six DH populations planted at two different locations in 2018 were hit by a severe frost event at the critical reproductive phase.Evaluation of frost damage showed that the plant growth stage most susceptible to low temperature(<2°C)was during the young microspore stage (10-18 days before anthesis).Out of the 30 frost QTL detected, 18 major QTL were mapped onto 13 chromosomes.Most frost QTL overlapped or were closely linked to the QTL for anthesis and maturity Zadok stages as well as to anthesis-related genes.However, the frost tolerance contribution by these QTL clearly stems from the late flowering alleles,illustrating the frost escape mechanisms and indicating that they are not useful in wheat breeding for water-limited environments.The mutatedVrn-A1a,Vrn-D1a,Rht-B1b, andRht-D1balleles and a high-copy numberPpd-B1allele contributed significantly to frost damage.Nevertheless, QTL or genetic factors outside the escape mechanisms were detected in the current study.Anthesis QTL were repeatedly detected in the proximity of theVrn-B1aregion and on chromosome 1B, whereas no frost QTL were detected on these two chromosomes.These striking results strongly imply that theVrn-B1aandTaFT3-1Balleles on 5B and 1B should be utilised in breeding for frost tolerance,as the early-flowering phenotype associated with these two genes is frost tolerant.Meanwhile,the recessive non-dwarfing allelesRht-B1a(rht1) andRht-D1a(rht2)associated with normal plant height could also be used in breeding for reproductive frost tolerance without delaying the flowering time.

CrediT authorship contribution

Jingjuan Zhang, MD Shahidul Islam, Wujun Maproposed the frost damage investigation.Jingjuan Zhangvisually scored the frost damage, implemented the molecular work and field work,analyzed the data and wrote the paper.Ben Biddulphsupervised frost damage phenotyping;Wujun MaandJorge E.Mayercontributed to the substantial paper writing and revisions.Junkang Rongsurveyed the parental lines and provided Bethlehem and CSL seeds for DH population construction;Kefei Chendesigned the field trials and statistically adjusted phenotype data.Zaid Alhabbar, Masood Anwar, andRongchang Yangparticipated in the DH population construction;Jingjuan Zhang,Yun Zhao,Hang Liu, Maoyun She, Mirza Dowla, Sonia Afrin, andNandita Roy

constructed the genetic linkage map;Jingjuan Zhang, MD Shahi-dul Islam,Wujun Ma,Yun Zhao,Masood Anwar,Zaid Alhabbar,Maoyun She, Rongchang Yang, Angela Juhasz, Guixiang Tang,Jiansheng Chen, Hang Liu, Yanjie Jiang, Shengnan Zhai, Xin Hu, Junkang Rong, Yingquan Zhang, Yebo Qin, Qier Liu, Zitong Yu, Sadegh Balotf, Yujuan Zhang, Mirza Dowla, Sonia Afrin,Nandita Roy, Md Resad Mallik, Md Atik Us Saieed,Shanjida Rahman, Nigarin Sultana, Sarah Al-Sheikh Ahmed,Chris Florides, Darshan Sharma, Nathan Height, Ben Biddulph,andMeiqin Lucontributed to the DH population seed bulk up and field trial work (planting, monitoring, onsite phenotyping,sampling, and harvesting) and recorded the flowering data.

Declaration of competing interest

The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.

Acknowledgments

This work was supported by Murdoch University and the Australia Grains Research & Development Corporation (GRDC)(grant number UMU00048), the Department of Primary Industries and Regional Development (DPIRD), Western Australia, and Kalyx Australia Pty Ltd.We thank Dr.Dean Diepeveen and his colleagues in DPIRD, and Drs.Hugo.Alsono-Cantabrana and Rowan Maddern from the GRDC for their support and assistance.We thank Mrs.Sue Broughton from the DPIRD for constructing the six DH populations used in this study.

Appendix A.Supplementary data

Supplementary data for this article can be found online at https://doi.org/10.1016/j.cj.2021.02.015.