ORIGINAL RESEARCH article
Sec. Plant Breeding
Identification of a Dominant Chlorosis Phenotype Through a Forward Screen of the Triticum turgidum cv. Kronos TILLING Population
- 1John Innes Centre, Norwich, United Kingdom
- 2Department of Plant Sciences, University of California, Davis, Davis, CA, United States
- 3Institute of Experimental Botany, Centre of the Region Haná for Biotechnological and Agricultural Research, Olomouc, Czechia
- 4School of Biosciences, University of Birmingham, Birmingham, United Kingdom
Durum wheat (Triticum turgidum) derives from a hybridization event approximately 400,000 years ago which led to the creation of an allotetraploid genome. The evolutionary recent origin of durum wheat means that its genome has not yet been fully diploidised. As a result, many of the genes present in the durum genome act in a redundant fashion, where loss-of-function mutations must be present in both gene copies to observe a phenotypic effect. Here, we use a novel set of induced variation within the cv. Kronos TILLING population to identify a locus controlling a dominant, environmentally dependent chlorosis phenotype. We carried out a forward screen of the sequenced cv. Kronos TILLING lines for senescence phenotypes and identified a line with a dominant early senescence and chlorosis phenotype. Mutant plants contained less chlorophyll throughout their development and displayed premature flag leaf senescence. A segregating population was classified into discrete phenotypic groups and subjected to bulked-segregant analysis using exome capture followed by next-generation sequencing. This allowed the identification of a single region on chromosome 3A, Yellow Early Senescence 1 (YES-1), which was associated with the mutant phenotype. While this phenotype was consistent across 4 years of field trials in the United Kingdom, the mutant phenotype was not observed when grown in Davis, CA (United States). To obtain further SNPs for fine-mapping, we isolated chromosome 3A using flow sorting and sequenced the entire chromosome. By mapping these reads against both the cv. Chinese Spring reference sequence and the cv. Kronos assembly, we could identify high-quality, novel EMS-induced SNPs in non-coding regions within YES-1 that were previously missed in the exome capture data. This allowed us to fine-map YES-1 to 4.3 Mb, containing 59 genes. Our study shows that populations containing induced variation can be sources of novel dominant variation in polyploid crop species, highlighting their importance in future genetic screens. We also demonstrate the value of using cultivar-specific genome assemblies alongside the gold-standard reference genomes particularly when working with non-coding regions of the genome. Further fine-mapping of the YES-1 locus will be pursued to identify the causal SNP underpinning this dominant, environmentally dependent phenotype.
Polyploidisation events underpin plant evolution and have been suggested to be key drivers of innovation, particularly within the angiosperms (Soltis and Soltis, 2016). All angiosperm species, including important crops such as wheat, rice, and maize, carry signatures within their genomes of ancient whole genome duplication (WGD) events that occurred within their lineage, such as the monocot-specific duplication, τ (Paterson et al., 2012). These polyploidisation events lead to the presence of multiple copies of genes which previously carried out the same function. It has been proposed that, following WGD, the resulting diploidisation of the genome leads to neo-functionalization or sub-functionalisation of gene copies derived from the original WGD (Dodsworth et al., 2016; Clark and Donoghue, 2018). The diploidisation process reduces the redundancy present within the genome by minimizing the number of genes with duplicate functions.
However, unlike rice and maize, wheat has also undergone two more recent allopolyploidisation events, where inter-species hybridizations bring together the chromosomes of each parent, creating a hybrid species with higher ploidy. The first event, approximately 400,000 years ago, occurred when two wild grasses hybridized to produce a tetraploid grass (wild emmer) which would go on to be domesticated as durum wheat (Triticum turgidum) (Dubcovsky and Dvorak, 2007; Borrill et al., 2019). The second polyploidisation event occurred more recently, only 10,000 years ago, when the tetraploid emmer hybridized with another diploid wild grass, leading to a hexaploid species which was then domesticated as bread wheat (Triticum aestivum). Unlike in ancient WGDs, these polyploidisation events have occurred relatively recently, such that most wheat genes are present as homoeologous duos or triads in durum and bread wheat, respectively, and may often have redundant functions (Ramírez-González et al., 2018).
A direct result of this homoeolog redundancy is that the inheritance of many traits in polyploid wheat tend to be quantitative, with multiple homoeologous loci contributing partly to the phenotype (Borrill et al., 2019; Brinton and Uauy, 2019). The phenotypic consequences of mutations in single homoeologs in wheat can be broadly classified into three categories — dominant (e.g., VRN1), whereby the mutant allele leads to a complete change in phenotype akin to mutations in diploids (Yan et al., 2003); additive (e.g., NAM, GW2), whereby mutants in each homoeolog lead to a partial change in phenotype which becomes additive as mutations are combined (Avni et al., 2014; Pearce et al., 2014; Borrill et al., 2018; Wang et al., 2018); and full redundancy (e.g., MLO), whereby the single and double mutants are similar to wildtype individuals, and only the full triple mutant (in hexaploid wheat) leads to significant phenotypic variation (Acevedo-Garcia et al., 2017). The presence of homoeolog redundancy, therefore, can hinder the use of forward genetic screens in polyploid wheat.
Therefore, beyond its status as an important crop, tetraploid durum wheat can provide a useful system to reduce the redundancy inherent in polyploid wheat. New advances in wheat genomics resources are increasing the speed and resolution with which we can now map loci corresponding to quantitative traits (Uauy, 2017). Recently gold-standard reference genomes for wheat were released, based on the hexaploid landrace Chinese Spring (International Wheat Genome Sequencing Consortium, 2018) and the tetraploid cultivar Svevo (Maccaferri et al., 2019). Additional wheat cultivars from across the globe are being sequenced as part of the wheat 10+ pan-genome project (10+Wheat Genomes Project, 2016). Crucially, this also includes durum wheat cultivar Kronos, which was used in the development of an in silico TILLING population (Krasileva et al., 2017). This mutant resource contains over 4 million chemically induced point mutation variation that can be rapidly accessed for a gene of interest through Ensembl Plants (Vullo et al., 2017).
An additional challenge when working in wheat is the sheer size of the genome, approximately 16 Gb in hexaploid and 11 Gb in tetraploid wheat. This is particularly important when designing sequencing strategies of mutant populations or individuals for mapping-by-sequencing. Various reduced representation methods exist for subsampling the wheat genome. These include gene-based methods through exome capture (Mamanova et al., 2010; Krasileva et al., 2017) or sequencing a specific gene family, as in R-gene enrichment sequencing (RenSeq) (Jupe et al., 2013; Steuernagel et al., 2016). However, these methods are less successful in obtaining variant information from non-coding regions due to their focus on genic regions. This is particularly important in the case of dominant phenotypes, which are often due to variations in regulatory regions that are not within the gene body (Yan et al., 2004; Fu et al., 2005; Borrill et al., 2015), although not exclusively (Simons et al., 2006; Greenwood et al., 2017). Methods do exist, however, to facilitate subsampling of the wheat genome while still retaining information from non-coding regions. In particular, chromosome flow sorting reduces the size of the genome by isolating an entire chromosome which can then be sequenced (Doležel et al., 2012). Other techniques (implemented in rice) include skim sequencing, which uses low coverage to obtain information about deletions or duplications, as well as SNPs, across the genome (Huang et al., 2009).
The new genomic resources and techniques for use in wheat now allows the more in-depth study of traits at the genetic level. One of these traits is chlorosis, characterized by the degradation of chlorophyll pigments in the leaves of the wheat plant. Chlorosis can be a symptom of disease, such as for yellow (stripe) rust and barley yellow dwarf virus (BYDV) (Kimura et al., 2016), as well as of many different nutrient deficiencies, such as nitrogen, magnesium, and potassium (Snowball and Robson, 1991). The presence of chlorosis, therefore, is often dependent on environmental conditions, such as the level of nutrients present in the soil. As a result, while chlorosis is not a desirable phenotype for breeding purposes, understanding the genetics underpinning chlorosis in different circumstances can increase our understanding of complex processes involved in nutrient homoeostasis and disease resistance.
Here, we used the Kronos TILLING population as a case study to identify and fine-map a novel locus in a tetraploid background (Krasileva et al., 2017). We performed a forward screen of the Kronos TILLING population for lines that exhibited late or early senescence phenotypes. From this set, we identified a line that segregated for a dominant chlorosis phenotype and was consistent across multiple years of field trials in the United Kingdom. We used mapping-by-sequencing to define the dominant phenotype as a single Mendelian locus on chromosome 3A, which we called Yellow Early Senescence-1. Using exome capture and chromosome flow-sorting to subsample the large wheat genome, we utilized the new RefSeqv1.0 hexaploid reference genome (International Wheat Genome Sequencing Consortium, 2018) alongside an assembly of the durum cultivar Kronos to identify SNPs across the region of interest. Following this, we mapped the Yellow Early Senescence-1 locus to 4.3 Mb, containing 59 high-confidence genes.
Materials and Methods
TILLING Population Screen
The initial screen of the sequenced Kronos TILLING population (N = 951 M4 lines) was carried out on un-replicated single 1 m rows (Supplementary Figure 1A), sown in November 2015 at Church Farm, Bawburgh (52°38′N 1°10′E). Note that all John Innes Centre (JIC) trials were sown at Church Farm, but in different fields at the farm in each year. Lines were sown in numerical order (i.e., line Kronos0423 was followed by Kronos0427). For simplicity, TILLING lines will be referred to as KXXXX throughout the manuscript (i.e., Kronos0423 as K0423). Wild-type controls (cvs. Kronos, Paragon, and Soissons) were sown randomly throughout the population. Rows were phenotyped for senescence as detailed below. Following scoring, independent M4 individuals from 10 mutant lines with early flag leaf and/or peduncle senescence and 11 mutant lines with late flag leaf and/or peduncle senescence were crossed in the glasshouse to wild-type Kronos (Supplementary Table 1). The F1 plants were then self-pollinated to obtain F2 seed (Figure 1E). For three mutant lines (K0331, K3085, and K3117) we recovered insufficient F2 seeds and hence these populations were not pursued further. All original mutant lines described are available through the JIC Germplasm Resources Unit1.
Figure 1. A premature yellowing phenotype from the Kronos TILLING population segregates as a single dominant locus. F2 populations of the K2282 Kronos TILLING line grown at the JIC in 2016 showed an early yellowing phenotype (A) Pigment content was measured in the yellow mutant plants (F2M) compared to the wild-type plants (F2WT) (B; n = 3 per genotype) and was also quantified using SPAD (C; n = 153 F2M, n = 61 F2WT). The yellow group (F2M) senesced significantly earlier than the late bulk (F2WT) (D; n = 148 F2M, n = 56 F2WT). Scoring of the plants demonstrated that the F2 population was segregating 3:1 for the yellow trait, indicative of a dominant single locus where heterozygous plants, such as the F1 generation, also present the mutant phenotype (E; numbers are combined for both populations). F2M and F2WT refer to plants which are yellow and green, respectively, and which derive from the F2 population (see bottom of E), while WT and MP refer to Kronos WT plants or M4 K2282 plants, respectively (see top of E).
F2 populations of the selected TILLING lines (backcrossed to cv. Kronos) were sown at Church Farm in March 2016 and grown as described previously (Harrington et al., 2019b). Briefly, individual F2 seeds were hand-sown in 6 × 6 1 m2 grids, leaving approximately 17 cm between each plant (Supplementary Figure 1B). In total, we sowed 31 F2 populations representing 18 distinct TILLING mutant lines. For K2282, two F2 populations were sown, K2282-A and K2282-B, and phenotyped. Seeds from both K2282 populations were taken forward for further field trials. 119 kg/ha of Nitrogen was applied to the trials during 2016.
In 2017 and 2018, the F3 seed from the K2282 F2 plants that were either heterozygous across the identified region on chromosome 3A or contained recombination within the mapped interval were grown. In 2017, we selected 30 lines from the K2282-A population and 8 lines from the K2282-B population. F3 seed from these 38 lines were sown in a randomized block design, replicated between 1 to 4 times depending on seed availability. Each experimental unit consisted of a 1 m2 plot that contained three 1 m rows of a single lines, separated from each other by ∼17 cm (Supplementary Figure 1A). The primary tillers of 12 individual plants from each row were tagged before heading. In 2018, 374 individual seeds derived from 16 F2 plants completely heterozygous across the SH467/SH969 region were hand-sown into a 1 m2 grid (Supplementary Figure 1B) and scored as in 2016. In each year, tissue for genotyping was sampled from the tagged plants (2017) or each individual plant (2018). Senescence phenotyping was carried out as detailed below. Precipitation data for the JIC field trials was obtained from a weather station at 52°37′ 52.29″ N, 1°10′ 23.57″ E. Nitrogen was applied to the field as fertilizer in both years, 147 kg/ha in 2017 and 124.5 kg/ha in 2018.
Based on the 2017 genotypic information, nine individual F3 lines genotyped as fully homozygous mutant (N = 4) or homozygous wild-type (N = 5) across the initial mapped region (from marker SH467 to SH969) were selected. Plants from these nine genotypes were sown in 2018 in 1 m2 plots (two double 1 m rows separated by approximately 33 cm; Supplementary Figure 1C) and replicated 3 times in a complete randomized design. Wild-type Kronos and M5 seed from the K2282 line were also sown as controls (N = 3). Two tillers in each row were tagged at heading and used for SPAD readings and genotyping. Senescence phenotyping was carried out as detailed below.
Putative Candidate Gene Phenotyping
In 2018, TILLING lines containing mutations in TraesCS3A02G414000 were grown in the field at JIC. Twenty-four seeds were hand sown per line in single, unreplicated 1 m rows (Supplementary Figure 1A). The TILLING lines were selected using the in silico TILLING browser implemented through Ensembl Plants (based on reanalysis of data from Krasileva et al., 2017). Lines containing frameshift and deleterious missense mutations (SIFT = 0) were identified, and three were selected for testing in the field based on mutation location and germination success of the TILLING line.
Davis, California Trial
91 F3 lines from population K2282-A and 55 F3 lines from population K2282-B were sown at the University of California Field Station near Davis, California (38° 31′ N, 121° 46′ W) in November 2016. Lines were selected if the F2 parent contained recombination within the SH467/SH969 region or was fully heterozygous across the region. In addition, seed from F2 parents completely mutant or wild-type across the region (N = 12 each) were also selected. Lines were sown in a complete randomized design, as double 1 m rows each separated by an empty row (as in the JIC 2018 trial; Supplementary Figure 1C). Eight individual plants were tagged per row at heading for plants derived from heterozygous parents, to allow genotyping and scoring of individual plants. At least two plants per double row were also tagged and sampled to verify the genotype of the completely wild-type or mutant lines. Heading and visual senescence was scored as in 2018 at the JIC, detailed below. 225 kg N/hectare was applied to the trial (as ammonium sulfate), half before planting and half on March 31st (Z30 stage). The trial was treated with appropriate fungicides to prevent stripe rust (Puccinia striiformis f. sp. tritici). Precipitation data for the Davis trial was obtained from the Davis, California weather station (38° 32′ 07″ N, 121° 46′ 30″ W).
F3 plants derived from mutant or wild-type F2 parents, genotyped across the SH179-SH969 region, were pre-germinated on damp filter paper for 48 h at 4°C in the dark. The seedlings were sown into P96 trays with 85% fine peat and 15% horticultural grit. Plants were transplanted to 1 L pots at the 3-leaf stage. The pots contained either (a) Petersfield Cereal Mix (Petersfield, Leicester, United Kingdom), (b) Horticultural Sand (J. Arthur Bower’s, Westland Horticulture), or (c) Soil taken from the Church Farm site used for JIC field trials (Bawburgh, United Kingdom). Plants sown into sand were also supplied with 100 mL of Hoagland solution every 3 days (Hoagland and Arnon, 1950). K2282 mutant and wild-type F3 plants were also tested under low water conditions in each of the three soil conditions listed above. Under the low water conditions, the plants were watered once weekly, and additionally to maintain a soil volumetric water content of approximately 20%, as measured with the Decagon GS3 sensor (ICT International, Armidale, Australia). Three plants of each genotype were treated in each condition. Plants were visually phenotyped for chlorosis onset, determined as a visual yellowing of the main flag leaf (see Figure 1A for a visual example).
Visual Color Scoring
Plants were scored as either yellow or green based on their visual appearance in the field. Those which phenocopied the wild-type Kronos plants were scored as green, while those which phenocopied the TILLING mutant parent (MP) were scored as yellow; the two categories are illustrated in Figure 1A. These initial visual scores were then corroborated with SPAD and chlorophyll scoring; see methods below.
Plants were scored for senescence across the different field trials as detailed previously (Harrington et al., 2019b). Briefly, when scoring individual plants, all phenotyping was carried out on the main tiller, tagged upon heading. Heading was scored at Zadoks growth stage 57, when the spike was 25% emerged (Zadoks et al., 1974). Flag leaf senescence was scored for the main tiller when 25% of the flag leaf showed visual yellowing and tissue death (necrosis) from the tip. Senescence of the main peduncle was scored when the top 3 cm were fully yellow. When scoring rows of the same genotype, all stages were scored across the entirety of the row. Rows were considered to have reached heading when 75% of the main spikes reached Zadoks growth stage 57. Leaf senescence was similarly scored when 75% of the flag leaves were yellowing and necrotic across 25% of the leaf, from the tip. Peduncle senescence was scored when the top inch of 75% of the peduncles were completely yellow.
Alongside visual scoring, we utilized the SPAD-502 meter (Konica Minolta, Osaka, Japan) to obtain non-destructive chlorophyll content readings. For measurements of individual plants (2016, 2017, 2018) eight readings were taken along the flag leaf on each side of the midrib and averaged to obtain a final reading which was considered the SPAD score for that biological replicate. For measurements of rows (2018), the two tagged tillers were both measured in the same way, and the average of their measurements was taken as the SPAD reading for that biological replicate.
Chlorophyll content was measured directly from sampled leaf tissue in 2016 and 2018 at JIC. In 2016, flag leaf tissue was sampled at heading (N = 3 per genotype); in 2018 flag leaf tissue was sampled at anthesis and the third leaf was sampled at the third leaf stage (Zadoks 13–14), approximately 24 days before anthesis (Mutant, N = 8; Wild-type, N = 10). In 2016, one leaf was sampled per individual plant and was treated as an independent biological replicate. Similarly, in 2018, one leaf was sampled per row, and treated as an independent biological replicate. Three 1 cm2 discs were extracted from each leaf, one at the base of the leaf, one in the middle, and one from near the leaf tip. Chlorophyll was extracted as described previously (Wellburn, 1994); briefly, the discs of tissue were soaked in N, N-Dimethylformamide (analytical grade, Sigma-Aldrich, United Kingdom) for 48–64 h until all pigment was completely removed from the leaf tissue. Pigment content was then quantified as previously described (Wellburn, 1994).
Leaf Mineral Content
Mineral content was taken from leaf tissue samples (2018). Leaf samples of approximately 0.2 g were dried and ground to a fine powder before digestion with 2 mL nitric acid (67–69%, low-metal) and 0.5 mL hydrogen peroxide (30–32%, low-metal) for 12 h at 95°C. Samples were then diluted 1:11 in ultrapure water before analysis with ICP-OES (Vista- PRO CCD Simultaneous ICP-OES; Agilent). Calibration was carried out using standards of Zn, Fe, and Mg at 0.2, 0.4, 0.6, 0.8, and 1 mg L–1 and Mn and P at 1, 2, 3, 4, and 5 mg L–1.
Thin sections of flag leaves were cut using a razor from mutant and wild-type plants in 2018 were imaged using a Leica MZ16 light microscope (Meyer Instruments, Houston, United States; N = 3 per genotype).
Bulked Segregant Analysis
Individual plants with green and yellow phenotypes from the K2282 F2 populations sown at the JIC in 2016 were selected for bulked segregant analysis. DNA from plant tissues, sampled at seedling stage, was extracted using the QIAGEN DNeasy Plant Mini Kit. The quality and quantity of the DNA was checked using a DeNovix DS-11 Spectrophotometer, Qubit (High Sensitivity dsDNA assay, Q32854, Thermo Fisher), and by running a sample of the DNA on an agarose gel (1%) to visualize the high molecular weight DNA. Four bulks were assembled by pooling DNA from plants which had been scored as either “yellow” or “green” (K2282-A, N = 75 for yellow, N = 16 for green; K2282-B, N = 33 for yellow, N = 22 for green). Equal quantities of DNA from the individual plants were pooled into each bulk to minimize bias.
Library preparation and sequencing was carried out at the Earlham Institute (Norwich, United Kingdom) as follows. DNA quality control was carried out using the High Sensitivity Qubit assay, before library preparation was carried out with a KAPA HTP Library Prep Kit. Size selection was carried out using Beckman Coulter XP beads, and DNA was sheared to approximately 350 bp using the Covaris S2 sonicator. Four libraries were produced, one for each bulk detailed above, which were barcoded and pooled. Five cycles of PCR were carried out on the libraries before carrying out exome capture.
Hybridization to the wheat NimbleGen target capture, previously described in Krasileva et al. (2017), was carried out using the SeqCapEZ protocol v5.0, with the following changes: 2.8 μL of Universal Blocking Oligos was used, and the Cot-1 DNA was replaced with 14 μL of Developer Reagent. Hybridisation was carried out at 47°C for 72 h in a PCR machine with a lid heated to 57°C.
The library pool was diluted to 2 nM with NaOH and 10 μL transferred into 990 μL HT1 (Illumina) to give a final concentration of 20 pM. This was diluted further to an appropriate loading concentration in a volume of 120 μL and spiked with 1% PhiX Control v3 before loading onto the Illumina cBot. The flow cell was clustered using HiSeq PE Cluster Kit v4, utilizing the Illumina PE_HiSeq_Cluster_Kit_V4_cBot_recipe_V9.0 method on the Illumina cBot. After clustering, the flow cell was loaded onto the Illumina HiSeq2500 instrument following the manufacturer’s instructions. The sequencing chemistry used was HiSeq SBS Kit v4. The library pool was run on two lanes with 125 bp paired end reads. Reads in bcl format were demultiplexed using the 6 bp Illumina index by CASAVA 1.8, allowing for a one base-pair mismatch per library, and converted to FASTQ format by bcl2fastq.
Chromosome Flow-Sorting and Sequencing
Seeds from the original K2282 M5 mutant line were used for the chromosome sorting and sequencing to ensure all parental SNPs were included. Suspensions of intact mitotic chromosomes were prepared from synchronized root tip meristems according to Vrána et al. (2000). To achieve better discrimination of individual chromosomes by flow cytometry, GAA microsatellite loci were fluorescently labeled by FISHIS (Giorgi et al., 2013) using FITC-labeled (GAA)7 oligonucleotides as described (Vrána et al., 2016). Chromosomal DNA was then stained by 4′,6-diamidine-2′-phenylindole (DAPI) at final concentration 2 μg/ml and the chromosome suspensions were analyzed by FACSAria SORP II flow sorter (BD Biosciences, San Jose, United States) at rates of 1000–2000 particles/s. Bivariate flow karyotypes DAPI vs. GAA-FITC were obtained and individual populations were flow sorted to identify the population representing chromosome 3A and to estimate the extent of contamination by other chromosomes (Supplementary Figure 2). Briefly, 2000 chromosomes were sorted onto a microscopic slide and evaluated by fluorescence microscopy after FISH with probes for GAA microsatellite and Afa-family repeat (Kubaláková et al., 2002). Three batches of 30,000 copies of chromosome 3A corresponding to ∼50 ng of DNA each were then sorted into PCR tubes containing 40 μl sterile deionized water. Chromosomal DNA was purified and amplified by Illustra GenomiPhi V2 DNA amplification Kit (GE Healthcare, Piscataway, United States) according to Šimková et al. (2008).
Library preparation and sequencing were carried out at Novogene. DNA integrity was confirmed on 1% agarose gels. A PCR-free library preparation was carried out, using the NEBNext Ultra II DNA Library Prep Kit for Illumina, following manufacturer’s instructions. Libraries were sequenced using a HiseqX platform, generating 150 bp paired end reads.
Sequencing Alignments and SNP Calling
For the bulked segregant analysis, the raw Illumina reads were aligned to the Chinese Spring reference genome, RefSeqv1.0 (International Wheat Genome Sequencing Consortium, 2018), using bwa-mem (v 0.7.5) with the default settings (-k 20, -d 100) (Li, 2013). Alignments were sorted, indexed, and PCR duplicates removed using SAMtools (v 1.3.1) (Wysoker et al., 2009), and SNPs were called using freebayes (v 1.1.0, default settings) (Garrison and Marth, 2012). Depth of coverage was calculated using the exome capture size detailed previously (Krasileva et al., 2017; Supplementary Table 2). Following SNP calling, we then filtered the original output to obtain only SNPs that were previously called in the K2282 parent line (Krasileva et al., 2017) using an original script available online to convert SNP coordinates to the RefSeq v1.0 genome2. The relative enrichment of each SNP in the yellow and green bulks was visualized across the wheat genome using the Circos package (Krzywinski et al., 2009). To do this, the AO/DP ratio was calculated for each of the green and yellow bulks, where AO is the number of reads with the mutant, or “alternate,” allele at the position in question, and DP is the depth of reads at that position. A Δ value was calculated by subtracting the AO/DP ratio for the green bulk from the yellow bulk. A schematic of the pipeline is provided in Supplementary Figure 3.
Following flow-sorting of chromosome 3A, reads were aligned to both RefSeq v1.0 and the Kronos assembly. We obtained access to the draft Kronos assembly produced at the Earlham Institute, which was assembled using the methods previously described (Clavijo et al., 2017a, b). The Kronos assembly is available in advance of publication from Grassroots Genomics3. In both cases, the alignment was carried out with bwa-mem (v 0.7.5; default settings -k 20, -d 100) (Li, 2013). Illumina reads from the wild-type Kronos assembly were aligned to RefSeq v1.0 using hisat (v 2.0.4, default settings with -p 8) (Kim et al., 2015). In all cases, files were sorted, indexed, and PCR duplicates removed with SAMtools (v 1.3.1) (Wysoker et al., 2009). For alignments to RefSeq v1.0, depth of coverage across part 2 of chromosome 3A was calculated using genomic windows of 1 Mb (Supplementary Table 2). Depth of coverage was not calculated for the complete Kronos alignment, as the scaffolds are not associated with a chromosome. SNPs were called on the respective alignments using freebayes (v 1.1.0) at default settings in all cases. BCFtools (Wysoker et al., 2009) was used to filter the SNPs based on quality (QUAL ≥ 20), depth (DP > 10), zygosity (only homozygous), and EMS-like status (G/A or C/T SNPs).
For the alignment against the Kronos genome, SNPs were manually filtered to remove those in regions of very high SNP density. We then identified scaffolds from the Kronos genome which fall within the YES-1 locus in the Chinese Spring RefSeq v1.0 genome using BLAST (v2.2.30; default parameters with -max_hsps 1 and -outfmt 6) (Altschul et al., 1990) against the gene sequences annotated within that region, using the v1.1 gene annotation. All further analysis of the SNP data for mapping and marker design focused solely on the 32.9 Mb YES-1 region. This retained 18 scaffolds which contained high quality SNPs within the YES-1 region.
Varietal SNPs between Kronos and Chinese Spring were identified by aligning the raw wild-type Kronos reads against RefSeqv1.0 (as detailed above). These varietal SNPs were then removed from the SNPs called from the alignment of the chromosome 3A reads to RefSeqv1.0. These SNPs were also manually curated to remove any which fell in regions of unreasonably high SNP density. A schematic of this workflow is provided in Supplementary Figure 4.
KASP Marker Genotyping
Markers were designed for the identified SNPs predominantly using the PolyMarker pipeline (Ramirez-Gonzalez et al., 2015b). Those not successful in PolyMarker were designed manually to be homoeolog specific. Markers were run on the recombinant populations using KASP genotyping, as previously described (Ramirez-Gonzalez et al., 2015a). Markers specific to K2282 are listed in Supplementary Table 3. Markers used for NAM-A1 genotyping were previously published (Harrington et al., 2019b).
Appropriate statistical tests for all data analyses were carried out and are detailed explicitly in the section “Results”. When needed, adjustments for false discovery rate were carried out using the Benjamini–Hochberg adjustment. This is referred to in the results as “adjusted for FDR.” All statistics were carried out in R (v3.5.1) (R Core Team, 2018), and data was manipulated using packages tidyr (Wickham and Henry, 2018) and dplyr (Wickham et al., 2019). Graphs of phenotyping and expression data were produced using ggplot2 (Wickham, 2016) and gplots (Warnes et al., 2019), respectively.
A Forward Screen of the Kronos TILLING Population Identifies a Line Segregating for a Dominant Chlorosis Phenotype
951 M4 lines of the Kronos TILLING population (Krasileva et al., 2017) were grown at the JIC in 2015 and scored for flag leaf and peduncle senescence timing. Ten lines showed early senescence phenotypes, while 11 showed late senescence phenotypes relative to Kronos wild-type (Supplementary Figure 5 and Supplementary Table 1). We developed F2 populations for these 21 lines crossed to wild-type Kronos. In 2016 the F2 mapping populations for 18 of these 21 lines were grown at JIC, and again scored for the senescence. From these populations, two showed significantly delayed peduncle senescence; K1107, with delayed peduncle senescence present in two independent F2 populations, and K2711, with delayed peduncle senescence in one of two F2 populations (Supplementary Figure 6). These two lines both contained mutations in the NAM-A1 gene, known to be a positive regulator of senescence (Uauy et al., 2006). The presence of the NAM-A1 mutation was sufficient to account for the variation in peduncle senescence timing found within the F2 populations for both K1107 and K2711 (Tukey’s HSD, p < 0.01, Supplementary Figure 7), indicating that the NAM-A1 SNPs were causal. The effect of the NAM-A1 mutations was followed up separately (Harrington et al., 2019b).
Based on the data from the 2016 field trials, we identified a single line, K2282, which showed a significant deviation in the timing of flag leaf senescence onset between the F2 population and the wild-type controls (p < 0.001, Kolmogorov-Smirnov test, adjusted for FDR; Supplementary Figure 6). Two F2 populations derived from K2282, K2282-A, and K2282-B, both showed earlier senescence compared to the wild-type controls. This phenotype, however, did not appear to be typical of a leaf senescence mutant. Although leaf senescence (scored based on leaf-tip necrosis) was indeed earlier in the K2282 populations, by anthesis the leaf tissue of individual plants was already highly chlorotic (Figure 1A). Quantification of chlorophyll levels confirmed that the yellow F2 individuals from both populations contained significantly less pigment than green F2 individuals (p < 0.05, Student’s t-test, Figure 1B). We also observed that the chlorosis phenotype predominated in the interveinal regions in the yellow plants, leading to a characteristic striated phenotype (Supplementary Figure 8).
We scored the K2282 F2 populations for chlorosis as a binary trait; i.e., plants were scored as yellow or green [see Figure 1A for an image of yellow (MP/F2M) and green (WT) flag leaves]. We confirmed that our visual scoring of the plants corresponded to the true chlorotic phenotype using non-destructive measurements of relative chlorophyll units. This identified a significant reduction in chlorophyll in the yellow (F2M) plants compared to the green (F2WT) plants, as expected (p < 0.001, Student’s t-test; Figure 1C). After classifying the F2 population into the green (F2WT) and yellow (F2M) groups, we found that the yellow group had significantly earlier leaf senescence (when scored to include necrotic symptoms) than the green group (p < 0.001 Kolmogorov–Smirnov test, Figure 1D). The yellow group showed a slight but significant delay of approximately 2 days in heading (p < 0.001, Wilcox test; Supplementary Figure 9) compared to the green group. This suggests that the early senescence observed is not due to the yellow group flowering earlier than the green group. The segregation of the chlorotic phenotype within the two populations was not significantly different from a 3:1 yellow to green ratio (X2, p = 0.07; Figure 1E), consistent with the trait being underpinned by a single dominant locus, hereafter referred to as Yellow Early Senescence 1 (YES-1).
The YES-1 Locus Maps to the Long Arm of Chromosome 3A
To map the trait, we carried out bulked segregant analysis on the two independent populations, K2282-A and K2282-B. A diagram of the analysis pipeline used is provided in Supplementary Figure 3. Following library preparation and exome capture, reads were aligned against the RefSeqv1.0 genome (International Wheat Genome Sequencing Consortium, 2018) and SNPs were called (Supplementary Table 2). To reduce the number of false SNP calls, we initially filtered the SNPs to only include those previously identified in the original M2 TILLING line (Krasileva et al., 2017). We recovered 1,548 SNPs out of the 3,060 SNPs present in the original K2282 M2 line which was sequenced. We expected to recover fewer SNPs than those identified in the original TILLING line as SNPs that were initially heterozygous in the M2 generation, may have been lost in the following two generations. Similarly, ∼50% of heterozygous mutations present in the M4 line crossed to wild-type Kronos to produce the F2 population would also have been lost.
We initially focused our analysis on the K2282-A population and calculated the ratio of the mutant (alternate) allele over total depth of coverage (AO/DP) at each SNP location in the yellow and green bulks (Figure 2A, inner track). From this, we then calculated the Δ value representing the enrichment of the mutant allele in the yellow bulk compared to the green bulk (Figure 2A, outer track). The segregation ratio seen in the field suggested this was a dominant single locus trait. Hence, we assumed that the yellow bulk would contain individuals homozygous or heterozygous for the causal mutant allele, while the green bulk should only contain homozygous wild-type plants. As a result, the AO/DP value should approach 0 in the green bulk, and 0.66 in the mutant bulk, and thus have a Δ value of 0.66. Using a conservative limit of 0.5 for the Δ value (gray line, outer track of Figure 2A), we identified only one region, on chromosome 3A, that was enriched for the mutant allele (Figure 2B). This result was consistent with that obtained from mapping carried out on the second population, K2282-B (Supplementary Figure 10).
Figure 2. Bulked segregant analysis identifies the YES-1 locus on chromosome 3A. Exome capture was carried out on yellow and green bulks from K2282 × Kronos F2 populations grown at the JIC in 2016. The K2282-A yellow bulk (yellow line, inner track; smoothed to a moving average of 4) and green bulk (green line) were scored at each SNP locus identified for enrichment of the mutant allele. The level of enrichment in the green bulk was subtracted from that of the yellow bulk to obtain the Δ value (outer track; smoothed to a moving average of 4). A high Δ value, indicative of a region enriched for mutant alleles within the yellow bulks, was identified on the long arm of chromosome 3A (B; smoothed to a moving average of 4). Markers designed on known TILLING SNPs within this region mapped the YES-1 locus to a 32.9 Mb interval within the F2 population (C). Green bars indicate wild-type calls, while yellow bars indicate mutant or heterozygous calls. The relevant phenotype of each recombinant group is shown to the left of the panel. The numbers of individual plants that fell into each recombination interval are shown to the right. The chromosome scale in (A) is given in Mb.
To validate this mapping, we developed KASP markers for the SNPs within and surrounding the region of interest (Figure 2C and Supplementary Table 3). Mapping of the individual F2 plants which were used to perform the exome capture confirmed the location of the region of interest on the long arm of chromosome 3A. Using the recombination events within this region and requiring at least two independent F2 plants to define the mapping interval, we narrowed the YES-1 region to between markers SH179 and SH969, a region of 32.9 Mb in the RefSeq v1.0 genome containing 345 genes (RefSeq v1.1 gene annotation) (Figure 2C).
Leaf Chlorosis Precedes Anthesis but Is Inconsistent Across Environments
To further characterize the phenotype, individual lines which were genotyped as completely mutant or wild-type across the YES-1 region were grown at the JIC in 2018. The mutant lines contained less chlorophyll A, B, and carotenoid pigment as early as the 3rd leaf stage (Zadoks 13–14) (Student’s t-test, p < 0.01; Figure 3A). This difference was increased at anthesis (Student’s t-test, p < 0.005), at which stage there was a larger spread in pigment content within the mutant lines than the wild-type lines. Chlorophyll content, measured with SPAD units, was also monitored across the development of the plants, from 14 days before anthesis to 39 days post-anthesis. SPAD readings were consistently lower in the mutant lines up to 24 days post-anthesis (p < 0.01, Pairwise Wilcoxon Rank Sum adjusted for FDR). The chlorotic phenotype remained highly visible on the leaves of the mutant plants, compared to wild-type (shown at 20 DPA, Figure 3C). In both wild-type and mutant lines, the level of chlorophyll in the flag leaf peaked at approximately 6 DPA (Figure 3B). No significant decline in SPAD units was observed in the wild-type plants until 24 DPA (p < 0.01, Pairwise Wilcoxon Rank Sum adjusted for FDR). In contrast, the mutant plants contained significantly less chlorophyll at 18 DPA compared to the peak at 6 DPA (p < 0.01, Pairwise Wilcoxon Rank Sum adjusted for FDR). Despite this earlier onset of senescence, the mutant lines continued to lose chlorophyll until the final stage of the time course (39 DPA), in line with the wild-type plants. We also found that the chlorosis phenotype is associated with significant decreases in leaf mineral content, with chlorotic leaves containing less magnesium at the 3rd leaf stage, and less of all four measured minerals at anthesis (Mg, Fe, Zn, p < 0.05; Mn, p = 0.05; Supplementary Figure 11).
Figure 3. The YES-1 locus causes lower chlorophyll levels before anthesis and earlier onset of senescence. The early chlorosis phenotype was recapitulated in the JIC 2018 field trials. Pigment content in the mutant lines is significantly lower at the third leaf stage (Zadoks 13–14, 24 days before anthesis) and becomes more extreme by anthesis (A; ∗∗p < 0.01; ∗∗∗p < 0.001 Student’s t-test). Relatively chlorophyll content, as measured with a SPAD meter, is significantly decreased in the mutant lines before anthesis, and remains significantly lower until 29 days post-anthesis (B; ∗∗p < 0.01; ∗∗∗p < 0.001 Pairwise Wilcoxon Rank Sum, adjusted for FDR). The yellowing phenotype in the leaves were clear in the field at 20 DPA (C).
Mutant and wild-type lines were also grown at UC Davis during the summer of 2016. Unlike in the United Kingdom, no chlorosis or senescence phenotype was observed either through visual or SPAD scoring (Supplementary Figure 12). This suggested that the causal locus underpinning YES-1 was environmentally dependent. Given the similarity between the interveinal chlorosis phenotype observed in the YES-1 mutant plants to that seen in plants with varying forms of nutrient deficiency (Snowball and Robson, 1991) and the decrease in leaf mineral content seen in the mutant plants (Supplementary Figure 11), we hypothesized that the environmental variation in phenotype may be due to nutrient content in the soil. To test this, F3 plants fully mutant across the YES-1 region were grown under glasshouse conditions in three soil types: standard glasshouse cereal mix, soil taken from the JIC field site in 2017, and horticultural sand supplemented with nutrient-replete Hoagland solution. However, none of the three conditions tested recapitulated the yellowing phenotype observed in the United Kingdom field (Supplementary Figure 13). This was surprising given the consistency of the phenotype at the JIC field site across four different fields during four successive field seasons (2015–2018).
We then investigated weather-related environmental variation across the two field sites and across years. We obtained rainfall and temperature data from Davis, CA, for the 2016–2017 growing season, and from the JIC field site for the 2016, 2017, and 2018 growing seasons. The trials carried out in California in 2017 received substantially more rainfall between sowing and heading than in any of the JIC trials (Supplementary Table 4 and Supplementary Figure 14). This suggested that perhaps reduced rain levels were correlated with the appearance of the mutant yellow phenotype. However, attempts to recapitulate the yellowing phenotype in the glasshouse through reduced watering of plants was also unsuccessful, as no early chlorosis or senescence was observed under different watering conditions (Supplementary Figure 13).
Fine-Mapping Reduces the YES-1 Locus to 4 Mb on Chromosome 3A
To identify further SNPs within the YES-1 locus, we purified chromosome 3A from the K2282 mutant by flow cytometry sorting. However, as the population of 3A chromosomes partially overlapped with the population of 7A chromosomes on a bivariate flow karyotype DAPI vs. GAA-FITC (Supplementary Figure 2), flow-sorted fractions comprised 80% of chromosome 3A and 20% of chromosome 7A as determined by microscopic observation. For sequencing, three batches of 30,000 chromosomes (∼50 ng) were flow-sorted and subsequent DNA amplification of three independent samples resulted in a total of 4.51 μg DNA.
Following sequencing, reads were mapped against the A and B genomes of the wheat RefSeq v1.0 genome (Supplementary Figure 4). 60.38% of reads aligned to chromosome 3A while 25.37% aligned to chromosome 7A, consistent with the expected contamination. The remaining reads (14.25%) mapped against the rest of the genome. We obtained on average 82X coverage across chromosome 3A, using genomic windows of 1 Mb.
In order to maximize our ability to discover novel SNPs in the YES-1 region, we carried out a simultaneous approach to SNP discovery utilizing both the Chinese Spring reference genome as well as the draft Kronos assembly, as depicted in Supplementary Figure 4. In brief, paired end sequencing of the K2282 mutant chromosome 3A was used to obtain high-quality SNPs outside of the previously captured exome. We used the Kronos assembly to identify SNPs in non-coding regions that are less conserved between the Kronos and Chinese Spring cultivars. In tandem, we took advantage of the contiguity of the RefSeq v1.0 genome which facilitated the identification of high-quality SNPs in and around all genes within the YES-1 locus.
Reads from the mutant K2282 chromosome 3A were mapped against the draft Kronos assembly and were filtered for homozygous, EMS-like SNPs, passing minimum quality and depth thresholds. To obtain only SNPs that fell within the physical region encompassed by the YES-1 locus, we carried out a BLAST between the Kronos scaffolds which contained SNPs and the Chinese Spring gene sequences within part of the YES-1 region. Conducting a BLAST against gene sequences within the YES-1 region, rather than the entirety of the region, reduced the number of scaffolds that mapped to the YES-1 region due to shared repetitive sequences rather than true synteny. Based on recombination seen in individual plants, we focussed on a region spanning 3 Mb upstream of marker SH179 and 3 Mb downstream of marker SH838, approximately 16 Mb in size. Within this region, we identified 18 unique Kronos scaffolds which both contained SNPs and at least one gene found in the RefSeqv1.0 YES-1 physical interval (Supplementary Table 5). 26 of the genes within the YES-1 region in Chinese Spring were identified (out of 345 total) within these 18 Kronos scaffolds. Genes that were not identified in the Kronos scaffolds may fall in scaffolds that contained no high-quality SNPs, may be split across multiple scaffolds, or may be absent from the Kronos genome. The SNPs within these scaffolds were manually curated, to exclude any regions that contained an unexpectedly high density of SNPs, leaving a final list of 16 scaffolds containing high-quality SNPs (Supplementary Table 5). These scaffolds covered an approximately 2.2 Mb of genome sequence, containing SNPs which fell mostly in non-coding regions (Supplementary Table 5). The SNPs underlying markers SH838 and SH179, initially identified in the exome capture data, were also recovered in the Kronos genome, validating the use of this method. KASP primers were designed for a subset of the SNPs and were used, together with the previous phenotypic data, to map YES-1 to a 6.6 Mb region between markers SH123480 and SH59985 (Figure 4A).
Figure 4. The YES-1 locus fine-maps to a 4.3 Mb region containing 59 genes. Markers were developed for novel SNPs identified in non-coding regions. We used phenotypic data from JIC 2017 and 2018 field trials to classify recombinant lines as green or yellow, indicated to the right of the panel. (A). These markers mapped YES-1 to a 4.3 Mb interval between markers SH044 and SH59985. Expression data for the 59 high-confidence genes in the region (B) from developmental time course data (Ramírez-González et al., 2018) highlights gene expression in root (yellow, top), leaf/shoot (orange, second from top), spike (green, second from bottom) and grain (purple, bottom) tissues across developmental stages. Genes mentioned in the text are highlighted by an asterisk (TraesCS3A02G412900 to TraesCS3A02G413200; OsSAG12 orthologs), a circle (TraesCS3A02G410800; Tryptophan Decarboxylase 2) or an inverted triangle (TraesCS3A02G414000; putative magnesium transporter). The full list of genes is provided in physical order in Supplementary Table 7.
To obtain more markers across the region, we also called SNPs against the Chinese Spring reference. To account for varietal SNPs between Kronos and Chinese Spring, we aligned raw reads from wild-type Kronos to the RefSeq v1.0 genome (Supplementary Figure 4). Using a subset of reads, we obtained a coverage of approximately 30X across chromosome 3A. SNP calling was then carried out against RefSeq v1.0 to obtain a list of varietal SNPs between Chinese Spring and Kronos. A total of 968,482 homozygous SNPs with quality greater than 20 and depth of coverage greater than 10 were identified across the second part of chromosome 3A, encompassing YES-1.
SNPs were then called between the K2282 mutant chromosome 3A reads and the RefSeq v1.0 genome. The set of SNPs was filtered for quality and depth and to exclude the varietal SNPs identified above. Following this filtering, a total of 7,153 SNPs were identified between markers SH123480 and SH969, a region of approximately 30 Mb. This was substantially more SNPs than would be expected from the known mutation density of 23 mutations/Mb for the Kronos TILLING lines (Krasileva et al., 2017). However, SNP density across the region was highly irregular which we hypothesized was due to mismapping and spurious SNP calling in repetitive regions.
To reduce the impact of repetitive regions on SNP calling, we extracted SNPs only from regions encompassing 1 Kb up and downstream of annotated genes within the YES-1 region. Following manual curation of SNPs, we identified a set of 15 SNPs that were located near genes within the annotated region (Supplementary Table 6). Of these SNPs, three were located in gene bodies (including the known TILLING SNP SH838), while the remainder were intronic (5) or fell in the promoter (5) or 3′ UTR (2). Note that some SNPs are in sufficient proximity to two gene models to be counted twice. Of the mutations in the coding region, SH838 and SH858 are both missense variants with low SIFT scores (0.00 and 0.03, respectively), while SH567 is a synonymous mutation. We designed markers based on these new SNPs and based on the JIC 2017 and 2018 phenotypic data we mapped YES-1 to a 4.3 Mb interval, between markers SH044 and SH59985 (Figure 4A).
Genes Within the Region
Within this region we identified 59 high-confidence genes based on the RefSeq v1.1 gene annotation (Supplementary Table 7). Using developmental time course data from two wheat varieties (Chinese Spring and Azhurnaya) (Borrill et al., 2016; Ramírez-González et al., 2018), we found that 25 genes within the region are expressed above 0.5 transcripts per million (TPM) in at least one stage of leaf or shoot tissue during development, consistent with our observation of a leaf-based phenotype (Figure 4). Of these genes, 18 were expressed above 0.5 TPM in leaf and shoot tissue during both vegetative and reproductive stages (Supplementary Table 7). This set of genes includes a putative magnesium transporter, TraesCS3A02G414000 (Gebert et al., 2009), which contains a missense mutation in the first exon of the gene which is predicted to be highly deleterious (SIFT = 0). This is the only gene within the 4 Mb region that contains a coding-region SNP, however, no chlorosis phenotype was observed for three additional lines with mutations in this gene (Supplementary Table 8). Within the 59 high-confidence genes, five genes were found to have senescence-related functions in their closest rice orthologues. A set of four tandem duplicated genes, TraesCS3A02G412900 to TraesCS3A02G413200, are orthologues to the rice gene OsSAG12-1, a negative regulator of senescence (Singh et al., 2013). A fifth gene, TraesCS3A02G410800, is orthologous to Tryptophan Decarboxylase 2, a rice gene that causes higher serotonin levels and delayed leaf senescence when over-expressed (Kang et al., 2009). All five genes with senescence-related phenotypes are lowly expressed or non-expressed across the set of tissues and developmental stages considered. However, the majority of the genes within the region are un-annotated and lack orthologous copies in either rice or Arabidopsis.
Here we have fine-mapped a region causing a dominant, environmentally dependent early-chlorosis phenotype. We have taken advantage of the recently released genetic and genomic resources for wheat to increase our ability to identify SNPs de novo in a Kronos TILLING mutant line. We have shown how the use of cultivar specific genome assemblies can be used to increase the ability to identify high-quality SNPs in non-genic regions which are often relatively less conserved between varieties than coding sequences.
Induced SNP Variation Can Lead to Novel Dominant Phenotypes
Many of the critical domestication alleles in polyploid wheat are derived from dominant mutations (Borrill et al., 2015; Uauy et al., 2017). This includes genes with critical variation in flowering time and free-threshing alleles resulting from dominant mutations (Yan et al., 2004; Fu et al., 2005; Simons et al., 2006; Greenwood et al., 2017). In wheat, the high level of redundancy between homoeologous genes adds to the importance of identifying dominant alleles to develop novel traits. Dominant alleles have retained their importance in modern breeding programs, underpinning the Green Revolution via the dominant dwarfing Rht-1 allele (Peng et al., 1999; Borrill et al., 2015). Most traits selected for in modern breeding programs, however, lack standing variation of dominant alleles in both the modern breeding pool and in older wheat landraces and progenitors. Instead, forward screens for phenotypes of interest typically identify multiple quantitative trait loci (QTL) that each contribute toward a small portion of the desired phenotype. These more complex effects, often caused by loss-of-function mutations, are inherently more difficult to identify due to the need to acquire/combine mutations in both or all homoeologous copies of a gene to attain a clear phenotypic effect (Borrill et al., 2015, 2019).
Here we have demonstrated that novel dominant alleles can be identified in chemically mutagenized TILLING populations (Krasileva et al., 2017). Forward screens of the TILLING population are most likely to identify novel traits caused by dominant mutations, given the low likelihood of obtaining simultaneous mutations in multiple homoeologous copies of the same gene. Indeed, the fact that mutations in NAM-A1 underpinned the only other senescence phenotype identified during this forward screen underscores this. The B-genome homoeolog of NAM-A1 is non-functional in Kronos; as a result, a single mutation in the A-homoeolog equates to a complete null and was sufficient to show a strong and consistent phenotype (Pearce et al., 2014; Harrington et al., 2019b).
The dominant phenotype identified in the K2282 line was particularly clear in that individual plants could unambiguously be scored for a binary green/yellow trait. However, we suggest that the TILLING population is equally well suited for forward genetic screens to identify novel dominant alleles governing other phenotypes. Recently, the Kronos TILLING population was used to identify a line which contained a deletion of Rht-B1, the partially dominant dwarfing allele (Mo et al., 2018). Here we have identified a novel dominant allele with no previously characterized genes located within the candidate region. This highlights the potential for novel dominant alleles to be identified in populations with induced variation, such as the Kronos and Cadenza TILLING populations (Krasileva et al., 2017).
The Use of Cultivar-Specific Assemblies Facilitates the Identification of Non-genic SNPs
A complication of working with dominant induced variation, however, is that dominant mutations may often act through changes to regulatory elements. Variation in the promoter and intron sequence of the flowering time gene VRN1 underpins the transition from winter to spring growth habit in wheat and barley (Yan et al., 2004; Fu et al., 2005). More recently, CRISPR editing has been used in tomato to edit the promoter region of various yield-related genes, leading to a high level of variation in trait morphology (Rodríguez-Leal et al., 2017). These results, amongst many others, highlight the importance of non-coding regions in regulating agronomically relevant traits. However, many reduced representation methods focus on enrichment of coding regions (Borrill et al., 2019). Such methods of genome complexity reduction, therefore, are less likely to contain the information needed to identify a dominant causal SNP in a regulatory region. Compounding this difficulty is the fact that non-coding regions of the genome are typically less conserved between cultivars. As a result, SNP identification against the reference variety may fail to identify critical SNPs or, conversely, identify a large number of spurious SNPs.
We have shown here that the draft Kronos assembly can instead be used, alongside non-biased methods of genome size reduction (e.g., chromosome flow sorting), to identify cultivar-specific SNPs in non-coding regions of the genome. We started by calling SNPs against scaffolds of the Kronos assembly, obtaining a large amount of SNP variation between the wild-type and mutant lines. Once we had this data, we then positioned the scaffolds which contained SNPs against the reference genome, identifying the SNPs which were located within our region of interest (here YES-1). This approach overcame two of the main drawbacks to using the reference genome and the Kronos assembly. On one hand, the reference genome would be expected to have different sequence content to another variety, such as Kronos, limiting its utility for SNP identification. On the other hand, unlike the gold-standard reference genome, the Kronos assembly does not have long-range assemblies needed to obtain positional information for SNPs. Long-range contiguous assemblies of additional cultivars, such as the recently published Svevo genome (Maccaferri et al., 2019), will greatly improve this current limitation. Until then, using variety-specific genomes, such as those being produced by the 10+Wheat Genomes project, alongside the highly contiguous reference genomes will facilitate the identification of non-genic SNPs.
Variability in Phenotype Points to an Environmentally Dependent Causal Locus
The early chlorosis and senescence phenotype caused by the YES-1 locus was consistent across 4 years in field trials at the JIC. However, mutant lines showed no evidence of a chlorotic phenotype when grown in Davis, CA. Comparison of rainfall and temperature patterns between the years and locations highlighted the fact that the plants received a high level of rainfall in Davis before flowering, substantially more than that received in any of the years at the JIC (Supplementary Table 4). This was due to the highly unusual wet winter that occurred in California in 2016/2017, with an average rainfall of 781 mm across the state from October 2016 to March 2017 (NOAA National Centers For Environmental Information, 2017). This suggested initially that the chlorosis response may be a response to higher water stress, yet we were unable to recapitulate the phenotype when grown in the glasshouse under different watering conditions. It is also possible that differences in temperature and photoperiod may influence the presence of the chlorosis phenotype. In particular, given the longer growing season in Davis, CA compared to at the JIC, as well as the difference in latitude, the photoperiod experienced by the Davis field trial would vary substantially from that experienced in the United Kingdom. However, as we also failed to recapitulate the yellow chlorosis phenotype in the glasshouse in the United Kingdom, under natural photoperiod in the spring and summer season of 2017, this does not seem to be the most likely source of environmental variation affecting the phenotype.
We also considered whether the phenotype was due to variation in soil nutrient content. The presence of a missense mutation within the coding region of a putative Mg2+ transporter (Gebert et al., 2009) highlighted this as a promising candidate gene. Similarly, the observed interveinal chlorosis phenotype (Supplementary Figure 8) is reminiscent of that characteristic of a magnesium deficiency (Snowball and Robson, 1991). However, we failed to recapitulate the phenotype when grown in the glasshouse using soil taken from the field at JIC, where the three trials were grown that showed a clear chlorotic phenotype (Supplementary Figure 13). Compounding this, we found that three Kronos TILLING lines which contained other SNPs within the transporter gene sequence did not show the same chlorotic phenotype (Supplementary Table 8). This included lines with both missense and premature stop codon mutations which lacked the exon containing the identified SNP in K2282. This implies that, if the magnesium transporter were the cause of the YES-1 phenotype, the specific missense mutation present in K2282 has a unique ability to cause a dominant change in function. As the transporter is predicted to function in a pentamer (Payandeh et al., 2013), it is possible that the mutation could be sufficient to prevent the pentamer to function effectively once formed, but not sufficient to prevent the mutant monomer from being incorporated into the pentamer. In this way it may be possible that plants heterozygous for the mutation show an equally strong phenotype as homozygous mutants as incorporation of a mutant monomer disrupts completely the function of the hexameric complex. This hypothesis could be tested in the future using Cas9-driven base editing in wheat to recapitulate the exact mutation in an independent background (Zong et al., 2017).
An alternative possibility is that a separate SNP located in a regulatory region may be acting either on the identified magnesium transporter, or on a separate, currently uncharacterised gene. Few dominant chlorosis phenotypes have previously been reported in the literature. A dominant chlorosis phenotype was previously reported in Brassica napus, however, this phenotype disappeared after budding (equivalent to heading in cereals) unlike here, where the yellowing phenotype became increasingly strong post-heading (Wang et al., 2016). In wheat, a Ygm (yellow-green leaf color) mutant has been identified with a semi-dominant phenotype where the heterozygous plants are an intermediate yellow-green color between the wild-type and homozygous mutant plants (Wu et al., 2018). This phenotype is underpinned by abnormal chloroplast development and is associated with differential expression of genes involved in chlorophyll biosynthesis and carbon fixation, amongst other traits. Further work to fine-map the YES-1 locus will hopefully shed light on the specific causal SNP underpinning the environmentally dependent chlorosis phenotype observed here, as well as on mechanisms governing dominant traits in polyploid wheat.
The raw reads from the exome capture and the flow-sorting experiments have been deposited on the SRA (PRJNA540141). The Kronos assembly is available from http://opendata.earlham.ac.uk/Triticum_turgidum/. Codes used for coordinate conversion and creation of the Circos plots (Figure 2 and Supplementary Figure 9) are stored on Github (https://github.com/Uauy-Lab/K2282_scripts). The original phenotypic data from the screen of all TILLING lines has been deposited on Dryad (doi: 10.5061/dryad.g3r3hp7).
CU and PB conceived the study. SH, PB, and NC carried out the field trials and phenotyping. SH carried out the mapping. MK and JD flow-sorted chromosome 3A, determined the purity in flow-sorted fractions and amplified chromosomal DNA for sequencing. SH and CU wrote the manuscript. All authors have read and approved the final manuscript.
This work was supported by the United Kingdom Biotechnology and Biological Sciences Research Council (BBSRC) through the Designing Future Wheat (BB/P016855/1) and GEN (BB/P013511/1) ISPs and an Anniversary Future Leader Fellowship to PB (BB/M014045/1). SH was supported by the John Innes Foundation. MK and JD were supported from ERDF project “Plants as a tool for sustainable global development” (No. CZ.02.1.01/0.0/0.0/16_019/0000827). NC acknowledges the funding received from Fulbright and Comisión Nacional de Investigación Científica y Tecnológica (CONICYT) Becas-Chile 72111195. This research was also supported in part by the NBI Computing infrastructure for Science (CiS) group through the HPC resources.
Conflict of Interest Statement
The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.
We would like to acknowledge A. Smith and M. Banfield (John Innes Centre) for helpful input during this project. We thank J. Vrána, M. Kubaláková and R. Šperková (Institute of Experimental Botany, Olomouc) for assistance with chromosome flow sorting and chromosome DNA amplification. We also thank J. Dubcovsky (University of California, Davis) for hosting the field trials at the UC Davis in 2016–2017. We also acknowledge the help of the field and glasshouse teams at the JIC and the UC Davis, as well as the help of G. Chilvers at the University of East Anglia ICP-OES platform. Finally, we thank B. Clavijo and the wheat pan-genome team at the Earlham Institute for allowing the use of the Kronos assembly before final publication. This manuscript has been released as a Pre-Print at https://www.biorxiv.org/content/10.1101/622076v1 (Harrington et al., 2019a).
The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fpls.2019.00963/full#supplementary-material
- ^ www.seedstor.ac.uk
- ^ https://github.com/Uauy-Lab/K2282_scripts; SNPs from K2282 were obtained from www.wheat-tilling.com
- ^ https://opendata.earlham.ac.uk/opendata/data/Triticum_turgidum/
10+Wheat Genomes Project (2016). The Wheat ’Pan Genome’. Available at: http://www.10wheatgenomes.com/ (accessed January 25, 2019).
Acevedo-Garcia, J., Spencer, D., Thieron, H., Reinstädler, A., Hammond-Kosack, K., Phillips, A. L., et al. (2017). mlo-based powdery mildew resistance in hexaploid bread wheat generated by a non-transgenic Tilling approach. Plant Biotechnol. J. 15, 367–378. doi: 10.1111/pbi.12631
Borrill, P., Harrington, S. A., Simmonds, J., and Uauy, C. (2018). Identification of transcription factors regulating senescence in wheat through gene regulatory network modelling. bioRxiv 456749. doi: 10.1104/pp.19.00380
Clavijo, B. J., Garcia Accinelli, G., Wright, J., Heavens, D., Barr, K., Yanes, L., et al. (2017a). W2rap: a pipeline for high quality, robust assemblies of large complex genomes from short read data. bioRxiv 110999.
Clavijo, B. J., Venturini, L., Schudoma, C., Accinelli, G. G., Kaithakottil, G., Wright, J., et al. (2017b). An improved assembly and annotation of the allohexaploid wheat genome identifies complete families of agronomic genes and provides genomic evidence for chromosomal translocations. Genome Res. 27, 885–896. doi: 10.1101/gr.217117.116
Doležel, J., Vrána, J., Šafář, J., Bartoš, J., Kubaláková, M., and Šimková, H. (2012). Chromosomes in the flow to simplify genome analysis. Funct. Integr. Genomics 12, 397–416. doi: 10.1007/s10142-012-0293-0
Fu, D., Szücs, P., Yan, L., Helguera, M., Skinner, J. S., Von Zitzewitz, J., et al. (2005). Large deletions within the first intron in Vrn-1 are associated with spring growth habit in barley and wheat. Mol. Genet. Genomics 273, 54–65. doi: 10.1007/s00438-004-1095-4
Garrison, E., and Marth, G. (2012). Haplotype-based Variant Detection from Short-read Sequencing. Available: https://ui.adsabs.harvard.edu/abs/2012arXiv1207.3907G (accessed July 01, 2012).
Gebert, M., Meschenmoser, K., Svidová, S., Weghuber, J., Schweyen, R., Eifler, K., et al. (2009). A root-expressed magnesium transporter of the Mrs2/Mgt gene family in Arabidopsis thaliana allows for growth in Low-Mg2+ environments. Plant Cell 21:4018. doi: 10.1105/tpc.109.070557
Giorgi, D., Farina, A., Grosso, V., Gennaro, A., Ceoloni, C., and Lucretti, S. (2013). FISHIS: fluorescence in situ hybridization in suspension and chromosome flow sorting made easy. PLoS One 8:e57994. doi: 10.1371/journal.pone.0057994
Greenwood, J. R., Finnegan, E. J., Watanabe, N., Trevaskis, B., and Swain, S. M. (2017). New alleles of the wheat domestication gene Q reveal multiple roles in growth and reproductive development. Development 144:1959. doi: 10.1242/dev.146407
Harrington, S. A., Cobo, N., Karafiátová, M., DoleŽel, J., Borrill, P., and Uauy, C. (2019a). Identification of a dominant chlorosis phenotype through a forward screen of the Triticum turgidum cv. Kronos TILLING population. bioRxiv 622076.
Harrington, S. A., Overend, L. E., Cobo, N., Borrill, P., and Uauy, C. (2019b). Conserved residues in the wheat (Triticum aestivum) Nam-A1 Nac domain are required for protein binding and when mutated lead to delayed peduncle and flag leaf senescence. bioRxiv 573881.
International Wheat Genome Sequencing Consortium (2018). Shifting the limits in wheat research and breeding using a fully annotated reference genome. Science 361:eaar7191. doi: 10.1126/science.aar7191
Jupe, F., Witek, K., Verweij, W., Sliwka, J., Pritchard, L., Etherington, G. J., et al. (2013). Resistance gene enrichment sequencing (RenSeq) enables reannotation of THE NB-LRR gene family from sequenced plant genomes and rapid mapping of resistance loci in segregating populations. Plant J. 76, 530–544. doi: 10.1111/tpj.12307
Kang, K., Kim, Y.-S., Park, S., and Back, K. (2009). Senescence-induced serotonin biosynthesis and its role in delaying senescence in rice leaves. Plant Physiol. 150, 1380–1393. doi: 10.1104/pp.109.138552
Kimura, E., Bell, J., Trostle, C., Neely, C., and Drake, D. (2016). Potential Causes of Yellowing During the Tillering Stage of Wheat in Texas. Texas A&M AgriLife Extension. Available: https://agrilifecdn.tamu.edu/texaslocalproduce-2/files/2018/07/Potential-Causes-of-Yellowing-During-the-Tillering-Stage-of-Wheat-in-Texas.pdf (accessed June 06, 2019).
Krasileva, K. V., Vasquez-Gross, H. A., Howell, T., Bailey, P., Paraiso, F., Clissold, L., et al. (2017). Uncovering hidden variation in polyploid wheat. Proc. Natl. Acad. Sci. U.S.A. 114:E913. doi: 10.1073/pnas.1619268114
Krzywinski, M., Schein, J., Birol, İ, Connors, J., Gascoyne, R., Horsman, D., et al. (2009). Circos: an information aesthetic for comparative genomics. Genome Res. 19, 1639–1645. doi: 10.1101/gr.092759.109
Kubaláková, M., Vrána, J., Číhalíková, J., Šimková, H., and Doležel, J. (2002). Flow karyotyping and chromosome sorting in bread wheat (Triticum aestivum L.). Theor. Appl. Genet. 104, 1362–1372. doi: 10.1007/s00122-002-0888-2
Li, H. (2013). Aligning Sequence Reads, Clone Sequences and Assembly Contigs with BWA-MEM. Available: https://ui.adsabs.harvard.edu/abs/2013arXiv1303.3997L (accessed March 01, 2013).
Maccaferri, M., Harris, N. S., Twardziok, S. O., Pasam, R. K., Gundlach, H., Spannagl, M., et al. (2019). Durum wheat genome highlights past domestication signatures and future improvement targets. Nat. Genet. 51, 885–895. doi: 10.1038/s41588-019-0381-3
Mamanova, L., Coffey, A. J., Scott, C. E., Kozarewa, I., Turner, E. H., Kumar, A., et al. (2010). Target-enrichment strategies for next-generation sequencing. Nat. Methods 7, 111–118. doi: 10.1038/nmeth.1419
Mo, Y., Howell, T., Vasquez-Gross, H., De Haro, L. A., Dubcovsky, J., and Pearce, S. (2018). Mapping causal mutations by exome sequencing in a wheat Tilling population: a tall mutant case study. Mol. Genet. Genomics 293, 463–477. doi: 10.1007/s00438-017-1401-6
Pearce, S., Tabbita, F., Cantu, D., Buffalo, V., Avni, R., Vazquez-Gross, H., et al. (2014). Regulation of Zn and Fe transporters by the Gpc1 gene during early wheat monocarpic senescence. BMC Plant Biol. 14:368. doi: 10.1186/s12870-014-0368-2
Peng, J., Richards, D. E., Hartley, N. M., Murphy, G. P., Devos, K. M., Flintham, J. E., et al. (1999). ‘Green revolution’ genes encode mutant gibberellin response modulators. Nature 400, 256–261. doi: 10.1038/22307
Ramírez-González, R. H., Borrill, P., Lang, D., Harrington, S. A., Brinton, J., Venturini, L., et al. (2018). The transcriptional landscape of polyploid wheat. Science 361:eaar6089. doi: 10.1126/science.aar6089
Ramirez-Gonzalez, R. H., Segovia, V., Bird, N., Fenwick, P., Holdgate, S., Berry, S., et al. (2015a). Rna-Seq bulked segregant analysis enables the identification of high-resolution genetic markers for breeding in hexaploid wheat. Plant Biotechnol. J. 13, 613–624. doi: 10.1111/pbi.12281
Rodríguez-Leal, D., Lemmon, Z. H., Man, J., Bartlett, M. E., and Lippman, Z. B. (2017). Engineering quantitative trait variation for crop improvement by genome editing. Cell 171, 470–480.e8. doi: 10.1016/j.cell.2017.08.030
Šimková, H., Svensson, J. T., Condamine, P., Hřibová, E., Suchánková, P., Bhat, P. R., et al. (2008). Coupling amplified Dna from flow-sorted chromosomes to high-density Snp mapping in barley. BMC Genomics 9:294. doi: 10.1186/1471-2164-9-294
Simons, K. J., Fellers, J. P., Trick, H. N., Zhang, Z., Tai, Y.-S., Gill, B. S., et al. (2006). Molecular Characterization of the major wheat domestication gene Q. Genetics 172, 547–555. doi: 10.1534/genetics.105.044727
Singh, S., Giri, M. K., Singh, P. K., Siddiqui, A., and Nandi, A. K. (2013). Down-regulation of OSSAG12-1 results in enhanced senescence and pathogen-induced cell death in transgenic rice plants. J. Biosci. 38, 583–592. doi: 10.1007/s12038-013-9334-7
Steuernagel, B., Periyannan, S. K., Hernández-Pinzón, I., Witek, K., Rouse, M. N., Yu, G., et al. (2016). Rapid cloning of disease-resistance genes in plants using mutagenesis and sequence capture. Nat. Biotechnol. 34, 652–655. doi: 10.1038/nbt.3543
Uauy, C., Distelfeld, A., Fahima, T., Blechl, A., and Dubcovsky, J. (2006). A NAC GENE regulating senescence improves grain protein, zinc, and iron content in wheat. Science 314, 1298–1301. doi: 10.1126/science.1133649
Uauy, C., Wulff, B. B. H., and Dubcovsky, J. (2017). Combining traditional mutagenesis with new high-throughput sequencing and genome editing to reveal hidden variation in polyploid wheat. Annu. Rev. Genet. 51, 435–454. doi: 10.1146/annurev-genet-120116-024533
Vullo, A., Allot, A., Zadissia, A., Yates, A., Luciani, A., Moore, B., et al. (2017). Ensembl Genomes 2018: an integrated omics infrastructure for non-vertebrate species. Nucleic Acids Res. 46, D802–D808. doi: 10.1093/nar/gkx1011
Wang, W., Simmonds, J., Pan, Q., Davidson, D., He, F., Battal, A., et al. (2018). Gene editing and mutagenesis reveal inter-cultivar differences and additivity in the contribution of TAGW2 homoeologues to grain size and weight in wheat. Theor. Appl. Genet. 131, 2463–2475. doi: 10.1007/s00122-018-3166-7
Warnes, G. R., Bolker, B., Bonebakker, L., Gentleman, R., Huber, W., Liaw, A., et al. (2019). gplots: Various R Programming Tools for Plotting Data. R package version 126.96.36.199.Google Scholar
Wellburn, A. R. (1994). The spectral determination of chlorophylls a and b, as well as total carotenoids, using various solvents with spectrophotometers of different resolution. J. Plant Physiol. 144, 307–313. doi: 10.1007/s11120-013-9834-1
Wickham, H., François, R., Henry, L., and Müller, K. (2019). dplyr: A Grammar of Data Manipulation. R package version 0.8.0.1.Google Scholar
Wickham, H., and Henry, L. (2018). tidyr: Easily Tidy Data with ‘spread()’ and ‘gather()’ Functions. R package version 0.8.2.Google Scholar
Wu, H., Shi, N., An, X., Liu, C., Fu, H., Cao, L., et al. (2018). Candidate genes for yellow leaf color in common wheat (Triticum aestivum L.) and major related metabolic pathways according to transcriptome profiling. Int. J. Mol. Sci. 19:1594. doi: 10.3390/ijms19061594
Yan, L., Helguera, M., Kato, K., Fukuyama, S., Sherman, J., and Dubcovsky, J. (2004). Allelic variation at THE VRN-1 PROMOTER region in polyploid wheat. Theor. Appl. Genet. 109, 1677–1686. doi: 10.1007/s00122-004-1796-4
Yan, L., Loukoianov, A., Tranquilli, G., Helguera, M., Fahima, T., and Dubcovsky, J. (2003). Positional cloning of the wheat vernalization gene Vrn1. Proc. Natl. Acad. Sci. U.S.A. 100, 6263–6268. doi: 10.1073/pnas.0937399100
Keywords: durum wheat, genomics, senescence, chlorosis, bulked-segregant analysis, TILLING, mapping-by-sequencing
Citation: Harrington SA, Cobo N, Karafiátová M, Doležel J, Borrill P and Uauy C (2019) Identification of a Dominant Chlorosis Phenotype Through a Forward Screen of the Triticum turgidum cv. Kronos TILLING Population. Front. Plant Sci. 10:963. doi: 10.3389/fpls.2019.00963
Received: 09 May 2019; Accepted: 10 July 2019;
Published: 24 July 2019.
Edited by:Agata Gadaleta, University of Bari Aldo Moro, Italy
Reviewed by:Filippo Maria Bassi, International Center for Agricultural Research in the Dry Areas (ICARDA), Morocco
Shichen Wang, Texas A&M University, United States
Copyright © 2019 Harrington, Cobo, Karafiátová, Doležel, Borrill and Uauy. This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner(s) are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.
*Correspondence: Cristobal Uauy, email@example.com