Marker-Assisted elimination of drought-susceptible accessions in upland rice breeding

Abstract

Breeding for water-deficit tolerance is fundamental to guarantee the sustainability of upland rice production, mainly due to the possibility of an increased frequency of drought episodes due to climate change. This work aimed to identify single nucleotide polymorphism (SNP) markers, derived from RNA sequencing (RNA-Seq), genome-wide association study (GWAS) and candidate genes from Arabidopsis, with potential for use in marker-assisted selection (MAS) for drought tolerance. RNA-Seq and GWAS were efficient in identifying useful SNP markers from the data obtained from three years of field experiments for 175 upland rice accessions, which were sequenced using 32 genes by Capture-Seq. Three genes were equally able to generate SNP markers that discriminated 95% of the 20 most drought susceptible accessions in the joint analysis of the experiments. The elimination of the genotypes with the unfavourable SNP allele reduced the initial number of accessions to one third, and transferring this result in a breeding routine, would enable to conduct smaller experiments per target location, increasing the precision and reducing the cost of the drought phenotyping.

Introduction

Drought is one of the major problems affecting crop yields, particularly for upland rice (Oryza sativa L.), mainly when cultivated in areas that depend exclusively on rainfall, resulting in low average grain yield (Serraj et al., 2011). Drought in rainy areas is expected to increase due to climate change, which can cause losses in productivity, especially if it occurs in the reproductive stage (Fedoroff et al., 2010; Kumar et al., 2014). In Brazil, the upland rice cultivated area has decreased by 57% in the last decade due to area competition with corn and soybean, and due to the lower yields (1.8 ton/ha) compared to lowland rice (7.0 ton/ha), contributing only about 38% of Brazilian rice production (Barbosa Filho and Yamada 2002; Heinemann et al., 2015). However, rice production in the lowlands can decrease due to reduced water availability around the world and, in this scenario, upland rice is a great alternative. To be feasible, future upland rice cultivars should combine the characteristics of resistance to drought in couple with high productivity when cultivated under conditions of water availability (Serraj et al., 2009; Fageria et al., 2014).

In upland rice, genetic variation for drought tolerance is available in gene banks (Pantalião et al., 2016), and the useful variability must be better exploited in order to be incorporated into the breeding programs. (Serraj et al., 2011; Kumar et al., 2014). However, the development of drought tolerant cultivars, combining high productive potential and drought tolerance is difficult, since drought is a quantitative trait, with low heritability, environmental interaction, and confounding effect with other abiotic stresses (Swamy and Kumar, 2013). Despite these issues, a great deal of effort has been carried out to obtain drought tolerant lines and cultivars, focusing in experiments in target environments, and the selection of many genotypes throughout the seasons and years (Degenkolbe et al., 2013; Swamy and Kumar, 2013).

Much information on the structural and functional genome of rice is publicly available, and such information can be used as a starting point for studies on various traits of interest. Genetic mapping, physiological characterization and expression analyses are available to potentially identify genes that have major effects on grain yield under drought treatment (Swamy et al., 2013). In addition to the development of drought tolerant cultivars by breeding programs, an important strategy would be to use recent scientific advances in biotechnology to improve this trait (Serraj et al., 2011). The complex genetic nature of drought tolerance, which hampers classic breeding, and poor understanding of the physiological and molecular mechanisms associated with this trait, make also the molecular breeding a difficult task.

Genes related to drought tolerance in rice have recently been obtained via mapping of QTLs (Dixit et al., 2014) and association mapping (Pantalião et al., 2016) by means of single nucleotide polymorphism (SNP) markers, functional genomics, microarrays (Parida et al., 2012; Borah et al., 2017), RNA-Seq (Huang et al., 2014; Silveira et al., 2015) or simply by using candidate genes homologous to those identified in other species (Shaar-Moshe et al., 2015; Abreu et al., 2016). The cost of high-resolution genotyping, as well as sequencing, has fallen year by year, which has opened the possibility of identifying genes related to drought tolerance via whole-genome resequencing (Jain et al., 2014). Nevertheless, the selection of drought-tolerant material based on these markers is still not a routine practice in marker-assisted selection (MAS) because of their inconsistency due to large QTL × E interaction, and no clear evidence exists that these markers have been successful in improving yield under drought stress in rice (Bernier et al., 2009; Swamy and Kumar, 2013). Most efforts to improve the drought tolerance of rice through MAS have been devoted to the use of low-resolution QTLs represented by large chromosomal segments, which only capture a small proportion of the total phenotypic variation (Babu, 2010). The possible reasons for this could be the inconsistency of QTL effects when considering multiple environments, and the detection of small effect QTLs on grain yield (Serraj et al., 2011; Dixit et al., 2014).

In the present work, a set of 175 upland rice genotypes (tropical japonica) were evaluated in three years of field experiments under water deficit, and also sequenced for a series of genes previously identified as being related to drought response. This work aimed to develop a strategy to select SNPs useful for a marker-assisted selection for drought tolerance in rice.

Material and Methods

Plant material

A set of 175 upland rice (Oryza sativa ssp. japonica) accessions from the Embrapa Rice Core Collection (Abadie et al., 2005) was evaluated, being 83 Brazilian landraces (acronym LRC), 47 improved varieties from Brazil (IVB) and 45 improved varieties from other countries (IVC), in addition to three checks (BRS Esmeralda, BRS Soberana and BRSGO Serra Dourada) (Supplemental Table 1).

Table 1: List of the 32 drought-related genes identified by RNA-Seq (1 to 12), GWAS (13 to 26) and candidate genes (27 to 32), sequenced by Capture-Seq.

Field trials assessing drought stress

The experiments were conducted in three years (2010, 2011 and 2012) at the Drought Phenotyping Experimental Station, located in Porangatu, Goias State, Brazil (49o06′W, 13o18′S, 396 m of altitude). For 2010 and 2011, the experimental design was Federer augmented blocks (FABs), and for 2012, it was randomized complete blocks, with two replications. The plot size consisted of four rows of three meters with a density of 60 seeds per meter. The experiments were monitored by tensiometers at a depth of 15 cm to maintain suitable soil water conditions (- 0.025 MPa). The drought was applied 30 days after emergence, when irrigation was reduced to approximately 50% until the end of the experiment. After the complete physiological maturity of the grains, the rice plants were manually harvested. Statistical analysis of grain yield data (kg/ha) was performed using the lme4 package from R software version 3.0.1 (The R Foundation for Statistical Computing, 2015) using the mixed model analysis procedure. Coefficient of variation (CV), accuracy and heritability for yield were estimated according to Resende and Duarte (2007). Estimates of the variance components were obtained by the REML (restricted maximum likelihood) method, according to Patterson and Thompson (1971), and the genetic values of everyone were estimated using the empirical best linear unbiased predictor. The genotype stability was estimated by the harmonic mean of the genotypic values (HMGV), and the adaptability was estimated by the relative performance of the genotypic values (RPGV), whereas the harmonic mean of the relative performance of genotypic values (HMRPGV) evaluated both the adaptability and stability of the genotypes (Resende, 2004).

DNA extraction and capture sequencing probe design

The DNA of 175 upland rice genotypes was obtained from the leaf tissue of seedlings using the Biopur Mini Spin Plant Extraction Kit (www.biometrix.com.br). The DNA samples were quantified on a Nanodrop spectrophotometer (ThermoScientific) and diluted to a concentration of 40 ng/μL. The integrity of the samples was assessed on 1% agarose gel. The probe set used for Capture-Seq were designed according Neves et al. (2013) for 32 genes (Table 1) whose selection was based on studies with RNA-Seq (Silveira et al., 2015), a genome-wide association study (Pantalião et al., 2016) and candidate genes (Abreu et al., 2016).

Capture sequencing analysis

The 175 upland rice genotypes were sequenced using the Capture-Seq methodology described in Neves et al. (2014). The DNA samples were sent to Rapid Genomics (USA) for capture sequencing (Capture-Seq). For each sequenced gene, a fragment of 120 nucleotides was selected and complimentary biotinylated RNA probes where synthesized. Prior to the enrichment step, Illumina-comparable libraries where prepared for each sample (Bentley et al., 2008). An input of 500 ng of DNA was used and fragmented to an average of 300 bp. Libraries were constructed by end-repairing the sheared DNA, A-tailing, ligating standard adapters to each molecule and enriching the libraries for ligated fragments with PCR amplification. The libraries were multiplexed prior to hybridization and 750 ng of pooled libraries were combined with the probes to enrich for the genes of interest following Sure Select Target Enrichment protocol (Agilent Technologies). The library pool was sequenced on a HiSeq 3000 platform (Illumina) according to the manufacturer's recommendations, with 100 cycles with paired-end fragments. Sequencing data was converted to fastq files using Illumina’s bcl2fastq software on default settings.

Identification of drought-related SNPs

The sequences from Capture-Seq analysis were trimmed by removing the low-quality bases, using the parameter phred < 20 at the 3' ends of each sequence in the Trimmomatic tool. MosaikAligner and Freebayes were used to align and call SNPs, respectively. To identify markers related to tolerance or susceptibility to drought, grain yield data from the 2012 experiment were used, since the genotypes were evaluated with repetition. The SNP profiles of the 20 genotypes with the highest yield and 20 genotypes with the lowest yield were compared.

Results

Drought phenotyping

The coefficient of variation, accuracy and heritability for yield for the 2010 experiment was 14%, 0.85% and 0.62%, respectively; for 2011, these values were 11%, 0.77% and 0.51%; and for 2012, they were 16%, 0.86% and 0.61%. Although the experimental design of the 2010 and 2011 trials was in FABs, these values were very similar to the values found in the 2012 experiment, which used randomized blocks, enabling direct comparisons among genotype and yield performances. The F test of the joint analysis of variance showed that the effects of year and year × genotypes were significant (p<0.001) for grain yield. Accessions above the yield average for the three field experiments were found for all groups of accessions (LRC, IVB and IVC). This means that the selection of drought-tolerant materials should not necessarily be prioritized in landraces, which are known to have unwanted traits for rice commercial cultivars, such as long cycle, tall size, and lower grain quality for trade.

Capture sequencing

The 175 individuals sequenced for the 32 genes generated a total of 276 SNPs. The set of genes with the highest number of SNPs identified were derived from GWAS, with 134 SNPs (mean of 9.6 SNPs per gene), the same average obtained by the set of SNPs from RNA-Seq. The set derived from candidate genes showed only 26 SNPs (19 of these SNPs came from LOC_Os06g08450) and the lowest average of SNP number per gene (4.3) (Table 1). The gene with the highest number of SNPs identified was LOC_Os10g31780 (35 SNPs, from RNA-Seq), whereas the gene with the highest SNP density was LOC_Os06g40570 (13 bp/SNP, from GWAS).

Drought-related SNPs

Considering all 276 SNPs found for 32 genes, 15 SNPs were informative (which were discriminant for the 20 drought-tolerant and 20 drought-susceptible accessions) from GWAS series (9 from LOC_Os08g13840 and 6 from LOC_Os02g09650) and 10 SNPs were informative from RNA-Seq series (all from LOC_Os09g21120) (Table 1). These three genes were located at positions S2_4959308, S8_8258609 and S9_12748551, respectively (from Nippon bare reference genome, available at RGAP website - http://rice.plantbiology.msu.edu). These SNPs presented the same classification pattern for the 20 most tolerant and susceptible genotypes (Table 2), i.e., they showed the AA SNP pattern for the same tolerant accessions, and for that reason, a single SNP (S2_4959308, S2 acronym, here-after) was considered from now on.

Table 2: SNP profile of the 20 most productive (20+) and 20 least productive (20-) accessions (2012 drought experiment in the field) for the three selected markers.

Of the 175 accessions evaluated with the SNP S2, 60 were genotyped as AA and 115 as GG. In seven accessions, was not possible to obtain the grain yield data for all three years, and for this reason, they were removed from the yield joint analysis. Considering the 168 remaining accessions, in all experiments and joint analysis the average yield of 56 AA genotypes were significantly higher than the average for 112 GG genotypes (Figure 1).

The percentage of SNP S2 in selecting genotypes with an above-average yield for two or three years (selection for AA genotypes) was just 50%. On the other hand, considering the genotypes with below-average yield for two or three years (selection for GG genotypes), the percentage was 76 and 80%, respectively (Table 3).

Table 3: Success of SNP S2 (AA or GG patterns) to predict the yield of accessions above or below the average in three years of drought experiments

The setup of three drought tolerance experiments showed to be fundamental to discriminate accessions with variable performance for this trait. For example, the cultivars Primavera and Soberana are known to be susceptible to drought, but in one experiment, their yield was above average. When the phenotypic and genotypic results were compared in the pairwise analysis of the experiments, the percentage of SNP S2 to discriminate for drought tolerance ranged from 45% (2011 × 2012) to 62% (2010 × 2012), whereas for susceptibility, it varied from 61% (2010 × 2012) to 87% (2010 × 2011) (Table 4).

Table 4: Coincidence percentages between yield and the SNP S2 pattern (AA or GG) for 168 upland rice accessions evaluated in three years of drought experiments. A: above the yield average; B: Below the yield average.

The high percentage of susceptible accessions identified by AA allele of the SNP S2 indicates that selection to discard GG accessions from the total set of accessions is advantageous for breeding purposes. Instead of evaluating 168 accessions, the 56 accessions with the AA genotype identified by SNP S2 would be selected, that is, only 33.3% of the accessions, and these would then be evaluated in smaller experiments, with a higher number of repetitions in the environment for which the drought-tolerant cultivars under selection are intended. Of these 56 AA selected accessions, 38 (68%) are landraces (LRC), while 9 (16%) are improved varieties from Brazil (IVB) and 9 (16%) are improved varieties from other countries (IVC). For the 112 GG accessions, the accessions were better distributed among groups, with 42 (38%) LRC, 36 (32%) are IVC and 34 (30%) are IVB.

Adaptability and stability analysis for drought yield

The parameter used to determine genotype adaptability was RPGV; HMGV was used to determine stability, whereas HMRPGV evaluated both adaptability and stability. The classification of genotypes was based on adaptability since the aim is to identify a genotype that is locally adapted to drought. However, the results were quite similar between the different parameters, including the yield (empirical best linear unbiased predictor, Eblups) of the joint analysis (Table 5). Considering the 20 most adapted genotypes, 8 are LRC, 7 are IVB and 5 are IVC; among the 20 least adapted, 9 are LRC, 5 IVB and 6 IVC (Table 5). SNP S2 genotyping analysis showed that 19 of the 20 drought susceptible genotypes (or 95%) in the joint analysis have the GG pattern.

Table 5: Estimates of adaptability (RPGV), stability (HMGV) and harmonic mean for adaptability and stability (HMRPGV) of the 20 most productive (20+) and 20 least productive accessions (20-), obtained from the joint analysis (BLUPj). The values between parenthesis are related to the ranking of the accessions according RPGV*GA. GA, yield general average.

Discussion

One of the great challenges of molecular genomics studies is to transfer the results to the routine practice of breeding programmes. The use of markers for assisted selection aims to accelerate the process of development of lines by prior identification of those lines possessing the favourable allele of genes of interest or markers linked to those genes of interest. This would increase the efficiency of breeding programmes because materials with a low potential for success in the programme would not be maintained in the experiments. Specifically, for the drought-tolerance trait, an initial screening that results in the reduction in the number of evaluated genotypes represents a significant saving of resources. This type of experiment is difficult to conduct due to the need for constant monitoring of the water-deficit conditions so that the deficit is high enough to kill the plants or low enough to prevent discrimination between the tolerant and susceptible plants. The translation of these markers to routine breeding based on marker assisted selection is difficult, because these traits are controlled for many loci that individually explain a small percentage of the phenotypic variation and the markers might be specific to a genetic background distinct from those evaluated or in the selection for environments different from those used in the experiments.

This work began with a gene selection, which were used to sequence 175 upland rice accessions. This set of accessions was evaluated at the Drought Phenotyping Experimental Station in Porangatu for three years (2010 to 2012). Considering the 276 SNPs obtained from the 32 sequenced genes, three were consistently associated with the 20 most productive accessions (SNP pattern AA) and the 20 least productive accessions (SNP patterns GG for S2, TT for S8 and CC for S9) in the 2012 experiment. The fact that the three genes located on chromosomes 2, 8, and 9 have the AA pattern related to higher productivity in the drought could have occurred due to chance, or somehow, being in linkage disequilibrium. This should be investigated in the future. According to prioritization (RiceNet V2; Lee et al., 2015) and co-expression network analyses (RiceFrend; Sato et al., 2013), these three genes are not connected. The gene LOC_Os02g09650 has an AP2 protein domain, and the ontology gene analysis shows that it has transcription factor activity (GO: 0003700), signal transduction (GO: 0007165), response to stress (GO: 0006950) and response to abiotic stimulus (GO: 0009628). In its gene model, it has only one exon and one transcript. This gene has been described as having a direct role in the response to drought (Dietz et al., 2010; Kazan, 2015). The gene LOC_Os08g13840 is the transcription factor WRKY25 (GO: 0003700 and GO: 0006355), and its gene model has two transcripts and two exons. This transcription factor has previously been associated with the response to drought by Tripathi et al. (2014). The action of transcription factors in response to drought has been well documented (Singh and Laxmi, 2015; Samad et al., 2017). The gene LOC_Os09g21120 (OsPUB6) is associated with catalysis of the transfer of ubiquitin from one protein to another (GO: 0004842), interacting selectively and non-covalently with any protein or protein complex (GO: 0005515) and the process in which one or more ubiquitin groups are added to a protein (GO: 0016567); furthermore, its gene model has only one exon and one transcript. Ubiquitination is a major type of post- translational modification of proteins that occurs in eukaryotic cells, and a cell death assay, in which one of the U-box genes was silenced, caused cell death in rice protoplasts, suggesting that ubiquitination might be involved in rice cell defence responses (Zeng et al., 2008).

Of the three genes associated with drought tolerance or susceptibility, two of these were identified by GWAS (LOC_Os02g09650, SNP S2; and LOC_08g13840, SNP S8), whereas one was identified by RNA-Seq (LOC_Os09g21120, SNP S9). The advantage of the GWAS analysis was that it was performed with phenotypic data obtained in field experiments, which indicates that such genes would be expressed, at least theoretically, in conditions similar to those found in crop cultivation. On the other hand, the operational cost to obtain the phenotypic data in the field of a given study population is much higher than that required to establish an RNA-Seq experiment in a greenhouse because field trials must be conducted at a phenotyping site where rainfall is prevented. In addition, the genotyping of tens to hundreds of accessions is costly. The analysis by RNA-Seq allows few accessions (even a single accession) to be evaluated in a greenhouse experiment, which allows more precise control of the water deficit. However, to identify the SNPs in these genes, a second round of analysis is necessary, such as by Capture-Seq approach used in this work. The SNPs related to drought tolerance obtained by RNA-Seq are identified in the transcribed genes, whereas the SNPs derived from a structural genome sequencing approach (such as genotyping by sequencing, or GBS) can be located both in genes and in intergenic sequences. The SNP 7, for example, selected to compose a marker panel for assisted selection for drought in rice (Pantalião et al., 2016), lies in an intergenic region, allowing a more comprehensive exploration of the genetic variability of structural DNA than the functional genome of rice.

As the SNPs derived from the three genes showed the same genotyping pattern (AA for the most productive genotypes and GG, TT and CC for the least productive genotypes under drought conditions), any of these three SNPs could be used in marker-assisted selection, with S2 being the one selected for this study.

Phenotypic analysis × genotypic analysis for drought tolerance

The genotype interaction per year was quite significant for grain yield. Of the 168 accessions analysed in all three experiments, 20 had above-average yield, while 35 were below average; that is, the remaining 113 accessions, or 67%, presented different productive performance in the three years of experiments. The comparative analysis of the AA pattern with the GG pattern of S2 marker showed that the most efficient way to perform the assisted selection was to select for the GG pattern of genotypes because it corresponded to 19 of the 20 least adapted genotypes in the joint analysis.

No distinction in adaptability was observed between the LRCs and improved varieties. That is, the breeder may prioritize the drought-tolerant materials derived from the breeding programme, as those materials already have other desirable phenotypes, such as plant architecture, cycle and grain quality. In contrast, when choosing landraces, additional work would be required to ensure that the derived inbred lines would have the favourable phenotype, due to the linkage between drought tolerance and undesirable traits (Vikram et al., 2015).

As drought tolerance is a complex trait, the identification of markers robust enough for use in assisted selection that can be used in multiple environments over several years and several populations is a difficult task (Kamoshita et al., 2008). The development of a technological product with these characteristics is difficult because of the different responses of the genotypes to the different environmental conditions that occur each year and at each site. The proposal to use the SNP S2 for the breeding programme is different, as this marker is based on the identification of the accessions that have the SNP GG pattern, which was mostly associated with accessions with grain yield below the experiments average. This would allow the selection of accessions with genotype AA, which, in the present study, represented 33.3% of the accessions evaluated. Fewer accessions (but with potential higher yield average under drought) would allow experiments with a greater number of replications to be performed, in at least three years, for a more consistent evaluation of the genotypes. This is a key point in regard to studying drought tolerance: the experiment must be repeated as often as possible to adequately estimate the productive potential of each genotype when grown under water deficit conditions. When selecting for a complex trait, such as drought tolerance, the following is important to consider: achieving balance between the size of the trial, the number of repetitions within the experiment, the number of years and the number of experimental sites. This decision should also consider the cost of implementing these trials.

Based on the results obtained in this work, SNP 2 assisted selection should be performed in a wide set of accessions and the selected ones would be evaluated in experiments with more repetitions for each location and for more than one experimental season. The removal of the accessions with the GG SNP 2 pattern could eliminate some of the most productive genotypes under drought conditions. However, considering the availability of drought tolerant materials in rice gene banks, the impact of the elimination of some of these accessions is reduced due to the discard of the clear majority of consistently drought susceptible accessions.

Conclusion

1) Both the RNA-Seq and GWAS methodologies were efficient in identifying SNP markers associated with drought tolerance or susceptibility. RNA-Seq is faster in identifying genes associated with the trait but requires a subsequent step to identify SNPs in these genes. GWAS requires the genotyping of many accessions and field experiments, but directly identifies SNP markers of both genes and intergenic regions.

2) The SNP markers derived from candidate genes obtained in the literature evaluated in this study were not useful for marker-assisted selection because they did not discriminate between the drought-tolerant and drought-susceptible individuals.

3) Of the 32 genes used in the identification of drought-tolerant genotypes, three could identify the susceptible genotypes consistently.

4) This work showed that at least three experiments per target site are required for genetic improvement for drought tolerance, regardless of the experimental design used.

5) The use of only one SNP marker means a significant reduction in the cost of genotyping, allowing many genotypes to undergo this screening. The selected accessions can then be used to conduct smaller, and consequently cheaper and more precise, experiments.

Acknowledgments

The Coordination for the Improvement of Higher Education Personnel/Ministry of Education (CAPES/MEC) for the grants to AFV; National Council for Scientific and Technological Development (CNPq) for the grants to CB and RPV and the Brazilian Agricultural Research Corporation Rice and Beans (EMBRAPA Rice and Beans) for support for this research.