Abstract
Populations can adapt to stressful environments through changes in gene expression. However, the role of gene regulation in mediating stress response and adaptation remains largely unexplored. Here, we use an integrative field dataset obtained from 780 plants of Oryza sativa ssp. indica (rice) grown in a field experiment under normal or moderate salt stress conditions to examine selection and evolution of gene expression variation under salinity stress conditions. We find that salinity stress induces increased selective pressure on gene expression. Further, we show that trans-eQTLs rather than cis-eQTLs are primarily associated with rice’s gene expression under salinity stress, potentially via a few master-regulators. Importantly, and contrary to the expectations, we find that cis-trans reinforcement is more common than cis-trans compensation which may be reflective of rice diversification subsequent to domestication. We further identify genetic fixation as the likely mechanism underlying this compensation/reinforcement. Additionally, we show that cis- and trans-eQTLs are under different selection regimes, giving us insights into the evolutionary dynamics of gene expression variation. By examining genomic, transcriptomic, and phenotypic variation across a rice population, we gain insights into the molecular and genetic landscape underlying adaptive salinity stress responses, which is relevant for other crops and other stresses.
eLife assessment
Working with a diverse panel of field-grown rice accessions, this valuable study measures changes in transcript abundance, tests for patterns of selection on gene expression, and maps the genetic basis of variation in gene expression in normal and high salinity conditions. The authors provide solid evidence that salinity treatment increases the number of genes with mean expression levels away from the optimum, and that a relatively small number of genes are hotspots for genetic variants that affect genome-wide patterns of variation in gene expression under high salinity. The design, clarity, and interpretation of several statistical analyses can be improved, additional opportunities for integration among datasets and analyses could be realized, and genetic manipulation would be required to confirm the functional involvement of any specific genes in regulatory networks or organismal traits that confer adaptation to higher salinity conditions. The manuscript will not only be of interest to evolutionary biologists studying the genetics of complex traits, but it will also be a resource for plant biologists studying mechanisms of abiotic stress tolerance.
Introduction
Plants face numerous stresses that reduce their growth and fitness, and they have a variety of adaptations to help them deal with environmental challenges (Mareri et al., 2022). For crop plants, stresses such as drought, heat and salinity can greatly limit productivity and agricultural sustainability worldwide (Fahad et al., 2017; Kopecká et al., 2023); indeed, as a direct result of various abiotic stresses, an estimated ∼ 50-70% of crop yields are lost (Francini and Sebastiani, 2019). While there is a wealth of information on the physiological responses of plants, including crops, to stress (H. Zhang et al., 2022), we are still developing an understanding of the underlying genetic responses to stress and the complex interactions between stresses and environmental cues, networks of gene expression regulation, physiological responses and ultimately plant fitness (Cramer et al., 2011; Siddiqui et al., 2021).
Much of the research on the genetic basis of stress responses in crops has focused on identifying individual genes associated with specific responses and/or tolerance traits (Cui et al., 2020; Jiang et al., 2019; Kumar and Wigge, 2010; Laohavisit et al., 2013; Yuan et al., 2014). Stress responses, however, are complex traits that can be influenced by multiple genetic pathways with numerous interconnected genes. Expression levels of genes are controlled by regulatory elements and are often environmentally induced, suggesting a critical role of regulatory divergence in adaptive evolution (De Clercq et al., 2021; Wilkins et al., 2016). Furthermore, changes in gene expression can lead to changes in the transcript network correlational structure, reshaping regulatory pathways (Lea et al., 2019) and influencing the way in which gene expression variation can influence adaptation to stressful environments (Signor and Nuzhdin, 2018; Wagner and Lynch, 2008). As such, there is an increasing need to understand how evolutionary dynamics of gene expression and transcript abundance relate to the genetic underpinnings of stress responses and adaptation in crops (Groen et al., 2020; Ruffley et al., 2023; Siddiqui et al., 2021).
One of the most important stresses for many plants is salinity. Soil salinity causes osmotic imbalance between the plant and the soil, which impedes the uptake of water and other key nutrients (Hakim et al., 2014; Munns, 2002), possibly leading to acute ion toxicity (Liang et al., 2018). While some plants are considered halophytic and can thrive in saline environments, other plants are highly sensitive to salts and experience negative effects of salinity even at low concentrations (Carillo et al., 2011). Oryza sativa (Asian rice) is one such salt-sensitive crop, facing significant yield loss due to soil salinity (Hussain et al., 2017; Melino and Tester, 2023). Rice responds to salinity stress by adjusting physiological and biochemical processes involved in osmotic and ion homeostasis, nutritional balances and oxidative stress (Castillo et al., 2007; Miller et al., 2010; Munns and Tester, 2008; Qin and Huang, 2020; Wang et al., 2012). Although studies have identified multiple genes associated with these processes (for reviews, see Ponce et al., (2021) and Liu et al., (2022)), we still lack a systems-level understanding of how gene expression variation mediates responses to salinity stress and evolves at the molecular level.
Here, we use an integrative system genomics approach to comprehensively dissect the genome-wide molecular and phenotypic response to salinity stress in rice. This study builds on prior work by our group that examined selection on gene expression in rice in response to normal and dry conditions (Ćalić et al., 2022; Groen et al., 2021, 2020). Using genomic, transcriptomic, and phenotypic datasets obtained from 130 diverse accessions of rice subjected to moderate levels of salinity stress, we (i) explore the selection on gene expression variation under salinity stress, (ii) dissect the genetic architecture of gene expression variation under saline conditions, and (iii) identify genes, molecular pathways, and salt stress response traits as well as associated trade-offs in the saline environment. We demonstrate that salinity stress induces increased selective pressure on gene expression, and we identify variation in biological processes and physiological traits that is beneficial and detrimental to plants in a saline environment, providing novel insights into the molecular landscape underlying an adaptive response to excess salt. We integrate these datasets with genomic sequence information to elucidate the genetic architecture and regulatory networks governing rice’s response to salinity stress.
Results
Variation in transcript abundance
To understand the microevolutionary dynamics of gene expression variation under salinity stress, we first investigated variation in transcript abundance in 130 accessions of O. sativa ssp. indica. We conducted a field experiment in the dry season of 2017 (January-May) at the International Rice Research Institute in Los Baños, Laguna, Philippines. Three replicates of each accession were planted separately in a normal wet paddy field as well as a similar field in which plants were exposed to moderate salinity stress (salt levels maintained at 6 dSm-1) maintained until maturity. Average fecundity was significantly lower in the saline field than in normal conditions (Fig. 1a; two-tailed paired t-test P = 1.658 x 10-8) with most of genotypes having significantly lower fecundity in the saline field (n = 94; one-sample proportion test P = 3.639 x 10-7).
![](https://prod--epp.elifesciences.org/iiif/2/99352%2Fv1%2Fcontent%2F596807v1_fig1.tif/full/max/0/default.jpg)
The strength and pattern of selection on heritable gene expression.
a, The O. sativa ssp. indica populations showed higher average fitness in the normal (blue) and saline (pink) field (two-tailed paired t-test P = 1.658 x 10-8) and fitness further showed a significant effects of genotype (G) and environment (E); genotype × environment (G×E) was not significant. Analysis of variance (ANOVA) [G and E (P < 0.001), G×E (P = 0.49)]; n=130 accessions. b, Broad-sense heritability (H2) distribution of Oryza sativa spp. indica transcripts. Two-way ANOVA, genotype FDR-adjusted q < 0.001, n = 130 accessions. c-e, The strength of linear selection |S|, linear selection differentials (S), and quadratic selection differentials (C) for genome-wide gene expression in normal (blue) and saline (pink) conditions. f, Conditionally neutral (light gray), and antagonistically pleiotropic transcripts (blue and magenta represent beneficial expression in normal and saline conditions, respectively). Black represents transcripts experiencing selection in the same direction in both environments (expression is beneficial or detrimental in both environments).
To examine gene expression in the field, mRNA levels were measured in leaf blades from 780 plants (130 accessions in triplicates for each environment) via 3’-end-biased mRNA sequencing. Leaf blades were sampled 38 days after sowing (DAS), corresponding to 7 days of plants being exposed to stress in the saline field. Population variance in gene expression of 18,141 widely expressed transcripts was partitioned into genotype (G), environment (E), and genotype × environment (G×E) effects (Supplementary Table 1). At a conservative false discovery rate (FDR) of 0.001, all but 3 transcripts displayed a significant genotype effect, indicating that expression levels of most transcripts are heritable. This effect of genotype is reflected in the broad-sense heritability (H2) distribution of gene expression levels (Fig. 1b), which had a median value H2 = 0.53 (range of 0.012 - 0.987; Supplementary Table 1). In addition to the high levels of heritability, 16,371 transcripts had a significant G×E term, indicating that for many transcripts, genotypes showed heritably different levels of expression in different environments. The significant G×E indicates genetic variation for plasticity, indicating that this plasticity can evolve. Although we see a widespread heritable plastic response in gene expression, we found little evidence of genotype-dependent plasticity in fitness (G×E for fitness P = 0.49; Fig. 1a), indicating that the G×E of transcripts does not translate to the complex trait of fitness. This could be due to a combination of factors, like gene interactions (pleiotropy and epistasis) leading to little to no effect on fitness, or environment specific genotype-dependent gene regulation (environment specific eQTLs). Furthermore, only a relatively small number of transcripts (254) showed significant variation due to E, indicating genotype-independent plasticity in gene expression for only a small fraction of genes.
Selection on gene expression
To identify transcripts associated with high fitness in normal and saline environments, we measured the strength of selection on gene transcript levels. We did this using phenotypic selection analysis (Lande and Arnold, 1983), taking the total number of filled rice seeds/grains (total fecundity) as a proxy for fitness (Groen et al., 2020). We estimated the linear (S) and quadratic (C) selection differentials, which are estimates for directional (negative or positive S) selection, and stabilizing (negative C) or disruptive (positive C) selection. We calculated the raw (S and C), variance-standardized (Sσ and Cσ), and mean-standardized differentials (Sμ and Cμ). We found that both mean and variance of transcript expression vary significantly between conditions (Mann-Whitney U-test, mean and standard deviation: P < 2.2 x 10-16), which means that the interpretation of the strength of selection based on either of these differentials could be misleading (Supplementary Table 2). To overcome this, and given that we were interested in selection on gene expression, which were all measured in the same units, we used the non-standardized raw selection differentials for all downstream analyses. We also report the mean- and variance-standardized selection coefficients (Supplementary Table 3).
Most transcripts were (nearly) neutral (|S| < 0.1) in both the normal and saline environments (Fig. 1c). Additionally, within the normal environment, there was a general trend towards stronger positive compared to negative selection on gene expression, indicating that higher expression of transcribed genes is associated with greater fitness (Table 1, Mann-Whitney U-test, normal: P = 3.95 x 10-15). In comparison, though the strength of positive selection was higher than that of negative selection in the saline field, higher expression of most of the transcribed genes is associated with lower fitness (Table 1, Mann-Whitney U-test, saline: P = 0.0068). Although no transcripts cleared the Bonferroni correction threshold in saline conditions, under normal conditions seventeen transcripts cleared the threshold, thirteen of which were under positive directional selection (Supplementary Table 3).
![](https://prod--epp.elifesciences.org/iiif/2/99352%2Fv1%2Fcontent%2F596807v1_tbl1.tif/full/max/0/default.jpg)
Summary statistics of selection on gene expression.
We also found that a high proportion of transcripts experienced stabilizing selection (C < 0) in both normal and saline environments. This effect was pronounced under normal conditions, with the strength of stabilizing selection being stronger than that of disruptive selection (Table 1; Mann-Whitney U-test, normal: P < 2.2 x 10-16). In contrast, within the saline environment, there was no detectable difference between stabilizing and disruptive selection (Mann-Whitney U-test, normal: P = 0.514), due to a higher proportion of transcripts experiencing stronger disruptive selection. This indicates that in the saline environment, as compared to normal conditions, an increase in fitness is associated with extremes in transcript abundance.
Comparing the distribution of selection differentials between environmental conditions, we found that selection was stronger in the saline field compared to the normal wet paddy field (Table 1). This was true for the overall strength of directional selection (Mann-Whitney U-test, P=0.012), as well as for negative (P=8.27 x 10-60.009), stabilizing (P < 2.2 x 10-16) and disruptive (P < 2.2 x 10-16) selection, but not positive selection (P=0.735).
Gene expression trade-offs
Next we compared the proportions of genes showing an opposite direction of selection in across environments (antagonistic pleiotropy-AP) and those showing selection in only one environment (conditional neutrality-CN) (Anderson et al., 2011). To account for the inherent bias associated with the detection of CN, we used a more stringent P-value cutoff to define CN (P < 0.025) in comparison to AP (P < 0.05) (Groen et al., 2020; Siddiqui et al., 2021). We found that 10.80% of the transcripts showed selection patterns consistent with CN, while only 0.28% of transcripts showed AP (Fig. 1f; Supplementary Table 4). The proportion of transcripts showing CN was greater than expected by chance and greater than the proportion showing AP (two-tailed proportion z-test P < 2.2 x 10-16), even accounting for the fact that there is an inherent bias favoring detection of CN (Anderson et al., 2013). These results are consistent with a lack of trade-offs in which the expression of a gene is favored in one environment and disfavored in another.
Among the 51 AP transcripts, increased expression of 38 was beneficial only in normal conditions (higher expression associated with higher fitness in normal conditions, but lower fitness under saline conditions). Gene ontology (GO) term analyses of these 38 transcripts indicated that many of these genes were involved in photosynthesis and metabolic processes (Supplementary Fig. S1; Supplementary Table 4). This is consistent with an observed reduction of photosynthesis in rice under salinity stress (Radanielson et al., 2018; Tsai et al., 2019).
We then wondered whether any difference underlie gene regulation of these AP transcripts relative to the non-AP transcripts, potentially indicating the genetic basis for this trade-off. We identified single nucleotide polymorphisms (SNPs) associated with transcript expression levels in each environment separately using whole-genome polymorphism data (expression quantitative trait loci [eQTL] analyses; see below). We found SNPs regulating expression of two photosynthesis related AP transcripts (PSAN and CRR7) only in the normal environment (Supplementary Fig. S2). We further looked at the 13 AP transcripts that were beneficial only in the saline environment (Supplementary Table 4). Although there was no significant GO enrichment for these transcripts, among them we identified a cyclophilin-encoding transcript (OsCYP2), which has been shown to confer salt tolerance in rice (Lee et al., 2015; Roy et al., 2022; Ruan et al., 2011). However, no SNP was associated with the expression of this gene in either normal or saline conditions.
Biological processes under selection
To investigate the broader biological processes associated with differential selection (strong selection in only one environment), we ranked all GO biological processes by their median selection strength in each environment and identified the processes with significantly stronger selection relative to their respective environment. We identified 13 and 18 processes that were under strong differential selection in normal and saline conditions, respectively (Supplementary Table 5; Fig. 2a). Processes primarily involved in various aspects of growth and defense were under stronger selection in normal conditions, whereas processes associated with regulation of flowering, cell cycle control and reproduction showed stronger selection under saline conditions (Fig. 2a). This provides insight into specific biological processes involved in flowering time regulation and reduced yield associated with salinity stress, as has been observed in multiple species (Chang et al., 2019; Li et al., 2007; Van Zandt and Mopper, 2002; R. Zhang et al., 2022). Furthermore, studies have found that the osmotic stress induced by salinity causes a reduction in the cyclin-dependent kinases (CDKs) responsible for cell-cycle transitions (G1/S and G2/M) (Ma et al., 2015; Schuppler et al., 1998; West et al., 2004). Aligned with this, our study also supports the notion that salinity stress affects cell cycle regulation, and leads to reduced growth and reproduction.
![](https://prod--epp.elifesciences.org/iiif/2/99352%2Fv1%2Fcontent%2F596807v1_fig2.tif/full/max/0/default.jpg)
Biological processes and pathways with differential responses to selection under saline conditions.
a, GO biological processes under stronger selection in normal (blue) and saline conditions (pink). b, Linear selection gradients (β), along with direct (D), indirect (I) and total (T) responses to selection on suites of transcripts in normal (blue) and saline conditions (pink).
Gene expression usually operates within the context of robust gene interaction/regulatory networks (Amiri et al., 2018; Israel et al., 2016; Ko and Brandizzi, 2020). The fact that genes interact in networks means that selection on one gene can lead to indirect selection on interacting genes, which can constrain the response of a population to selection on gene expression (Agrawal and Whitlock, 2010; Groen et al., 2020; Kondrashov and Houle, 1997). To examine this phenomenon, we identified suites of correlated transcripts using principal component (PC) analysis on genome-wide gene expression levels. We then estimated the linear (β) and quadratic selection (³) gradients for PCs explaining over 0.5% of variance in each environment (Supplementary Table 6). Although quadratic selection was generally weak, linear selection on some PCs showed significant directional selection (Supplementary Table 7).
Using the breeder’s equation (Falconer and Mackay, 1996), we predicted the response to selection on these PCs to examine the constraints on a population’s microevolutionary response to selection on gene expression. We found that over half of PCs (7 of 11 and 7 of 12 PCs in normal and saline conditions, respectively) displayed opposite signs of direct and indirect responses to selection (Supplementary Table 7). This finding contrasts with prior work in rice under dry and wet conditions that found a lack of constraint in response to selection for most traits except seed size because direct and indirect responses to selection were largely similar (Ćalić et al., 2022), indicating that different types of stress may either constrain or facilitate the response to selection.
We found PCs enriched for metabolic pathways and biosynthesis of phenylpropanoids and secondary metabolites to be under selection in the normal environment (Fig. 2a; Supplementary Fig. S3; Supplementary Table 7). Different transcripts involved in these pathways were under positive and negative selection which may act to keep these pathways in a steady-state. Interestingly, circadian rhythm was found to be under positive selection, with an overall positive response to selection, in the saline environment (Fig. 2a; Supplementary Fig. S3; Supplementary Table 7). This is in alignment with the role of circadian clock genes in conferring salt tolerance (Kim et al., 2013; Wei et al., 2021; Xu et al., 2022), and indicates a tentative increase in expression of circadian clock genes with continuous exposure to soil salinity.
Salinity stress induces decoherence
Gene expression levels are generally correlated, but these correlations can be perturbed by environmental stresses, a phenomenon that is termed decoherence (Lea et al., 2019). Such decoherence has been demonstrated in humans and primates (Lea et al., 2019; Pu et al., 2022; Watowich et al., 2022), but little is known about how stress alters the correlation structure of transcript levels in plants. To examine decoherence in rice, we utilized the recently developed CILP (Correlation by Individual Level Product) method (Lea et al., 2019), which detects the systematic loss of correlation in gene expression among individuals. Since CILP calculates product correlations for all possible pairs of genes, we used transcripts with selection strengths greater than 0.1 (|S|> 0.1) in at least one environment with expression greater than 0 in at least 50% of individuals (2,051 transcripts; Supplementary Table 8) to reduce data dimensionality.
We found that the correlation structure for gene expression was broadly similar across the 2.10 M unique transcript pairs, but > 29,000 transcript pairs (involving 1,742 unique transcripts) showed significantly different correlations between normal and saline conditions (FDR < 5%) (Fig. 3a; Supplementary Table 8). This change in correlation structure indicates possible divergence in gene interactions between environments, which may presage a restructuring of gene networks. GO term enrichment analysis of the transcripts that show decoherence with significant pairs greater than the median (median significant pair per transcript = 12, n = 853) highlighted important pathways related to plant responses and potentially tolerance to excess salt (Fig. 3b). For instance, circadian rhythm genes like OsPRR37 and GIGANTA (GI), which have been shown to confer salt tolerance in rice and Arabidopsis (Kim et al., 2013; Wei et al., 2021), differed significantly in their interactions between the two environments. Similarly, to manage the energy requirements of salt stress, the relative abundance of metabolites involved in energy producing pathways like glycolysis, tricarboxylic acid (TCA) cycle and various other metabolic processes have been shown to be altered as an early response to salt stress (Che-Othman et al., 2017; Jacoby et al., 2011; Yang and Guo, 2018). Additionally, aligned with our results, multiple studies have reported a positive correlation between salt tolerance and levels of secondary metabolites and amino acids with osmoprotectant properties associated with lowering osmotic stress (Krishnamurthy and Bhagwat, 1990; Petrusa and Winicov, 1997; Tari et al., 2010). Our results suggest that salt exposure induces decoherence of gene expression in some transcripts, that this decoherence results from a restructuring of the gene expression network, and this restructuring can allow salt stress tolerance, providing a potential molecular mechanism underlying this tolerance.
![](https://prod--epp.elifesciences.org/iiif/2/99352%2Fv1%2Fcontent%2F596807v1_fig3.tif/full/max/0/default.jpg)
Salinity stress induces regulatory decoherence.
a, Pearson correlation coefficients between pairs of transcripts (|S| > 0.1 and expression greater than 0 in at least 50% individuals) in normal (x-axis) and saline conditions (y-axis). Pink and blue represent pairs with correlation stronger in saline and normal conditions, respectively; gray represents correlation that is not significantly different between conditions. b, Enrichment of transcripts with significant pairs greater than the median (median significant pair per transcript = 12, n = 853) involved in regulatory decoherence post salt exposure.
Selection on traits
In addition to gene expression data, we also collected phenotypic data for 13 organismal traits in normal and saline conditions, which provides us with the opportunity to examine how salt stress influences complex phenotypes and identify connections between variation in traits with fitness consequences and underlying patterns of gene expression. Since these traits were measured on different scales, we estimated variance-standardized selection gradients on these traits, focusing only on seven traits that were not strongly correlated to limit the contribution of indirect selection (Pearson correlation coefficient < 0.6; Supplementary Fig. S4). We identified three traits – leaf osmotic potential (LOP), chlorophyll a content (Chl_a), and flowering time (FT) – that displayed different selection patterns in normal versus saline environments (Fig. 4; Supplementary Table 8).
![](https://prod--epp.elifesciences.org/iiif/2/99352%2Fv1%2Fcontent%2F596807v1_fig4.tif/full/max/0/default.jpg)
Traits with different selection profiles under salt stress.
Linear (β) and quadratic (γ) selection gradients on the traits LOP (leaf osmotic potential), Chl_a (chlorophyll a content), and FT (flowering time). Error bars represent mean ± SE; dots and asterisks indicate significance of selection-gradient at two-sided unadjusted P < 0.1 and P < 0.05.
Leaf osmotic potential is a key trait associated with water transport in plants, and which decreases with increasing salinity, leading to low water uptake by the plant. This trait was under positive selection in the control environment but under stabilizing selection under saline conditions, which suggests that an optimal LOP is important under salt stress (Fig. 4). Furthermore, we found Chl_a content at the reproductive stage to be under negative selection in saline conditions. Although ion toxicity has previously been shown to reduce chlorophyll content (Ashraf and Bhatti, 2000; Taïbi et al., 2016), our results showed that this reduction is associated with increased survival and reproductive fitness, consistent with the general trend for reduced photosynthesis under salinity stress. We also found that FT was under positive selection (selection for later flowering) in normal conditions but under stabilizing selection in saline conditions (Fig. 4). Moreover, FT was significantly reduced under saline conditions (two-tailed paired t-test P = 0.0012; Supplementary Fig. S5), implying that salt stress selects for an optimal flowering time, with the trend shifted towards earlier flowering compared to that in the normal wet paddy.
Genetic architecture of gene expression variation
To dissect the genetic architecture of gene expression variation, we identified expression quantitative trait loci (eQTLs) and examined whether and how selection acts on these eQTLs. We identified cis- and trans-eQTLs (FDR < 0.001) regulating expression of 3,065 and 3,277 genes in normal and salinity stress conditions, respectively (Supplementary Table 9). We observed that 49.64% (29,623 of a total 59,669) of cis-eQTLs were common in both environments, as compared to 18.62% (25,528 of a total 137,100) of trans-eQTLs; this result was robust to FDR cutoff (FDR of 0.01 and 0.05). Moreover, this is consistent with previous observations from studies on other species that have found 48%-77% overlapping cis-eQTLs and 9%-60% common trans-eQTLs across environments (Smith and Kruglyak, 2008; Snoek et al., 2012; Sterken et al., 2023). It further indicates that trans-eQTLs might be more environment-specific than cis-eQTLs. We tested this explicitly by identifying loci showing gene-environment interaction (G × eQTL) and found that at FDR of 0.05, cis-eQTLs constituted merely 0.28% (142 of 50718) of the total identified G × eQTLs.
We identified eQTL hotspots, which are regions of the genome that are associated with expression variation of a large number of genes (Qu et al., 2018). These regions can occur either due to low amounts of recombination (high linkage disequilibrium-LD) or because they contain master-regulators that pleiotropically control expression of multiple functionally-associated genes (Hammond et al., 2011; Kliebenstein, 2009; West et al., 2007). To account for LD and identify regions likely to contain master-regulators, we chose to focus on the subset of genes regulated by the lead-SNPs (SNP with the most significant association) in a given 100-kb region. Through this approach, we identified 11 trans-eQTL hotspots (number of unique genes > 30) in saline conditions, but none in normal conditions (Supplementary Fig. S6; Supplementary Table 10). These results indicate that gene regulation under stress conditions may be dependent on a few master-regulators, as has been shown for drought stress in rice and maize (Kuroha et al., 2017; Liu et al., 2020). Interestingly, one of the hotspots, Chr7: 25.9-26.0Mb, influenced the expression of a disproportionately high number of genes (191 genes) and contains the gene OsFLP, a R2R3 myb-like transcription factor that regulates stomatal development (Wu et al., 2019). OsFLP has recently been associated with salt tolerance in rice (Qu et al., 2022).
It has been shown previously that when a gene is regulated via both cis and trans factors, their effects tend to drive target gene expression in opposite directions, canceling the combined effect on expression. This is commonly referred to as cis-trans compensation, which is thought to arise due to stabilizing selection that maintains similar gene expression over evolutionary timescales in the face of new mutations (Coolon et al., 2014; Goncalves et al., 2012; Landry et al., 2005; Lovell et al., 2018). Contrary to this expectation, we found that most genes appear to have a reinforcing (cis and trans effect in the same direction) rather than a compensatory pattern (two-tailed proportion z-test P = 1.65 x 10-11; Fig. 5a; Supplementary Table 11). Among the 524 cis-trans co-occurring genes in saline conditions, 317 were reinforcing (∼60.5%) and 207 were compensatory (39.5%). This pattern held for normal conditions as well (Supplementary Fig. S7a; Supplementary Table 11). Thus, expression patterns of these genes are expected to either increase or decrease expression over evolutionary timescales, indicating the presence of directional selection.
![](https://prod--epp.elifesciences.org/iiif/2/99352%2Fv1%2Fcontent%2F596807v1_fig5.tif/full/max/0/default.jpg)
Genetic architecture of gene expression variation during salt stress.
a, Effect sizes of genes with both cis and trans factors under saline conditions showing excess of reinforcing cis-trans (teal) in comparison to compensating cis-trans (salmon). b, Inter-varietal variation in gene expression for genes under compensating control is significantly lower than for those under reinforcing control; one-sided Mann-Whitney P = 0.0028. c, Frequency distribution of MAF (minor allele frequency) for cis-eQTLs (blue) and trans-eQTLs (yellow) in saline conditions against the genome-wide background (gray).
We hypothesized that this excess of reinforcing rather than compensatory effects could be driven by the extensive diversity among rice landraces that arose from crop diversification following domestication (Meyer and Purugganan, 2013). To test this hypothesis, we examined whether genes under cis-trans reinforcement showed evidence of higher inter-varietal variation in gene expression (population-wide variance between mean accession expression levels) as compared to genes under cis-trans compensation. Supporting our hypothesis, we found significantly higher inter-varietal expression variation for genes under cis-trans reinforcement in saline (Fig. 5b; Mann-Whitney U-test = 0.002; mean-inter-varietal expression variation compensating = 1.12, reinforcing = 1.27) and in normal (Supplementary Fig. S7b) conditions. Furthermore, it has been suggested that cis-trans compensation/reinforcement can arise due to the genetic fixation of compensating/reinforcing trans-regulatory variants (McManus et al., 2014; Signor and Nuzhdin, 2018), which could lead to elevated linkage disequilibrium (LD) between the cis- and trans-regulatory variants acting on a gene. This was indeed the case: estimated LD between all pairs of cis- and trans-regulatory variants for the identified compensating and reinforcing genes was significantly higher as compared to the background LD (mean r2 compensating: 0.33, reinforcing: 0.36, background = 0.056; two-tailed permutation test P = 0.0019).
To examine the pattern of past selection acting on cis- and trans-eQTLs (before, during and after rice domestication, including the crop diversification phase) we plotted the folded site-frequency spectrum (SFS) of the minor allele frequencies and inferred the relative strength of historic selection (Joly-Lopez et al., 2020) We found that trans-eQTLs for both environments had been under strong purifying selection, with the SFS being significantly left-shifted compared to background SNPs (two-tailed t-test P < 10-16; mean MAF background = 0.20, trans normal = 0.0936, trans saline = 0.0910). Not surprisingly, given their highly pleiotropic nature, this effect was more pronounced for trans-eQTLs regulating multiple genes (two-tailed t-test P < 10-16; normal mean MAF in unique trans-eQTLs = 0.099, multiple trans-eQTL = 0.073; saline mean MAF in unique trans-eQTLs = 0.097, multiple trans-eQTL = 0.071). In contrast to trans-eQTLs, we found that cis-eQTLs in both normal and saline conditions had potentially been under balancing selection (Fig. 5c, Supplementary Fig. S8) with the SFS being significantly right-shifted for cis-eQTLs relative to background SNPs (two-tailed t-test P < 10-16; mean MAF background = 0.20, cis normal = 0.238, cis saline = 0.237). Moreover, for both conditions, the estimated nucleotide diversity (π) for 100-kb regions flanking cis-eQTLs was also significantly higher as compared to genome-wide 100-kb blocks (two-tailed t-test P < 10-11; mean π background = 0.286, cis normal = 0.297, cis saline = 0.297).
Comparing patterns of past selection for cis-eQTLs specific to each condition revealed no significant difference between the normal and saline environments. However, we did observe trans-eQTLs in saline conditions to be under stronger purifying selection (two-tailed t-test P = 0.017). Taken together, this indicates that cis- and trans-eQTLs are under different selection regimes in rice, and that purifying selection is stronger on trans-eQTLs regulating gene expression under salt stress conditions, giving us insights into the macroevolutionary dynamics of gene-expression variation.
Discussion
Gene expression is a key link in the chain of adaptive organismal responses of organisms to environmental challenges. In this study, we used a systems genomics approach to examine genome-wide transcript levels in rice in both normal and moderately saline field conditions to examine the evolutionary response and genetic architecture of gene expression under salinity stress. We found selection on expression of a set of genes, and for some genes, selection on expression differed among environments, indicating that salinity can select for changes in gene expression. We saw that the genetic regulatory pathways can be modified by salinity, as indicated by decoherence, which provides a potential genetic mechanism underlying salinity tolerance. We found limited evidence for trade-offs given a lack of antagonistic pleiotropy in the fitness effects of expression compared across environments, and we also found that there did not appear to be strong cis-trans compensation in gene regulation. These results provide us with insights into the fitness consequences and genetic architecture of gene expression under salt stress in rice, as we discuss below.
One of the primary goals of this study was to characterize selection on gene expression and how this selection varies between saline and normal conditions. Before discussing our results, it is important to note that our conclusions are influenced by the way in which selection differentials are standardized. Selection differentials are estimated through regression coefficients in the relationship between trait values and fitness (Lande 1979). Because traits are often measured on different scales, selection differentials are usually standardized to allow for meaningful comparisons and interpretation. The most common method of standardization is variance standardization, in which traits are standardized to a mean of zero and a standard deviation of one (Lande and Arnold, 1983). With this approach, the selection differential, multiplied by the heritability, gives the expected response to selection in standard deviation units. So while this is a valid approach, there is concern that when variance-standardized selection differentials are interpreted as reflecting the “strength” of selection based on their magnitude, they can be misleading because trait variance can vary among traits and environments, such that the trait variance and the magnitude of selection are conflated (Hereford et al., 2004). Instead, the use of mean-standardized selection estimates have been suggested (Hereford et al., 2004; Matsumura et al., 2012). However, it has not been widely recognized that mean-standardized selection estimates can also be biased if mean trait values vary consistently across conditions, such that magnitude of selection is conflated with trait means. We show here that this was indeed the case: the mean gene expression values differed consistently between saline and normal conditions, leading to biases when comparing the strength of selection across treatments when mean-standardized selection differentials were used. Specifically, using mean-standardized selection coefficients, we found selection on gene expression to be stronger under normal conditions, whereas unstandardized and variance-standardized selection coefficients indicated stronger selection under saline conditions. This result indicates that it is important for researchers to carefully consider the advantages and potential drawbacks of both variance-standardized and mean-standardized, as well as unstandardized, selection differentials in making decisions about which to use in any particular system and depending of the goals of the study, but also to consider these ideas in interpreting selection differentials as reflecting the strength of selection.
Because we found that both variance-standardized and mean-standardized selection gradients were biased, as explained above, we draw our conclusions about patterns of selection using unstandardized selection gradients. We believe this approach to be appropriate given that gene expression values are measured on the same scale and normalized for every transcript, and our goal is to identify transcripts under selection and compare selection across treatments. Using unstandardized selection coefficients, our study shows that most genes appear to have nearly neutral levels (|S| < 0.1) of selection on transcript levels in both normal and saline environments. This is consistent with previous studies conducted in various plant and non-plant species, which have shown that most traits (including transcript abundance) are under very weak selection, and only a few traits experience strong selection at microevolutionary timescales (Ahmad et al., 2021; Hoekstra et al., 2001; Kingsolver et al., 2001). There does seem to be a tendency, however, towards stronger positive selection over negative selection on gene expression, meaning that higher expression of a majority of genes is generally associated with higher fitness. Looking at quadratic selection, we found evidence for both stabilizing and disruptive selection on gene expression. Our results are consistent with the expectation of stabilizing selection being more common than disruptive selection in normal conditions (Kingsolver et al., 2001), and reinforce the idea that disruptive selection is more widespread under stress conditions (Ahmad et al., 2021), potentially reflecting the prevalence of frequency- and density-dependent competition for resources under stress (Kingsolver and Pfennig, 2007).
Comparing the distribution of selection coefficients, we further found that selection was stronger in moderate stress saline field conditions compared to normal wet paddy field conditions, indicating that salinity stress induced increased selective pressure on gene expression at microevolutionary timescales. Although exposure to salinity stress increased the strength of selection on gene expression, this increase was relatively low compared to what has been reported for drought stress in rice (Agrawal and Whitlock, 2010; Groen et al., 2020; Kondrashov and Houle, 1997). This could potentially be attributed to the levels of stress experienced by the plants – the drought stress treatment was more severe in terms of fitness after stress exposure (Agrawal and Whitlock, 2010; Groen et al., 2020; Kondrashov and Houle, 1997) in comparison to the moderate salinity stress in this study – but further work is needed to examine whether there is indeed a relationship between stress intensity and selection strength. While this work provides insight into the direction of selection and the predicted degree of evolutionary change under the specific conditions of this study, more work would need to be done over multiple generations and under additional conditions to make more realistic predictions of microevolutionary change in the field.
Our work contributes to an ongoing debate regarding whether the intensity of selection on gene expression increases under stress. We found in this study of rice and salinity, and in a prior study on drought (Groen et al. 2020), that directional selection was greater under stress compared to normal conditions. This finding is consistent with prior studies on other systems like fruit flies (Drosophila melanogaster) (Groen et al., 2020; Jasnos et al., 2008; Kondrashov and Houle, 1997) and yeast (Groen et al., 2020; Jasnos et al., 2008; Kondrashov and Houle, 1997), but contrasts with a recent study in D. melanogaster reporting no change in selection due to increasing levels of nutritional stress (Arbuthnott and Whitlock, 2018). Our results also support prior findings in rice (Ćalić et al., 2022; Groen et al., 2020) and in Boechera stricta (Anderson et al., 2013) that antagonistic pleiotropy occurs and can be important for local adaptation, while conditional neutrality in gene expression may be much more common. Our finding indicates that it may be feasible to breed salinity tolerant rice varieties without a yield penalty under non-saline conditions, though further work is needed to verify this.
We examined the genetic architecture of gene expression variation as well, and found that trans-eQTLs rather than cis-eQTLs are primarily associated with rice gene expression under salinity stress, potentially via a few master-regulators. Trans-eQTLs may be important in environment-dependent gene expression changes, while cis-eQTLs may be more robust to environmental changes as has been observed in other species (Smith and Kruglyak, 2008; Snoek et al., 2012; Sterken et al., 2023). This can be attributed to trans-eQTLs’ larger mutational target along with a larger effect of drift than of positive selection in fixing them. These results further corroborates with other studies that have found a more dominant effect of trans-eQTLs on within-species variation than between-species variation (Emerson et al., 2010; Metzger et al., 2016; Wittkopp et al., 2008).
We found that cis- and trans-eQTLs show different patterns of selection, with cis-eQTLs showing evidence for balancing selection and trans-eQTLs showing purifying selection. We see this pattern for rice, a largely selfing species, in contrast to outcrossing species where both cis- and trans-eQTLs have been found to be under negative selection (Hernandez et al., 2019; Josephs et al., 2015). Future studies will be needed to investigate whether the mating system has an influence on the pattern of natural selection acting on variants regulating gene expression. Finally, and contrary to expectations, we show that under saline conditions, cis-trans reinforcement is more prevalent than cis-trans compensation. This result may be driven by rice domestication and subsequent population diversification. Additionally, we find significantly elevated levels of LD among the cis- and trans-regulatory variants for compensating and reinforcing genes, indicating the role of genetic fixation as an underlying mechanism for cis-trans compensation/reinforcement.
Systems genomic approaches provide insights into large-scale patterns across genome-wide data over multiple scales. This approach can also identify possible genes or genetic pathways that may prove critical in various processes. Our analysis, for example, shows that many of the potential antagonistically pleiotropic genes are involved in photosynthesis and metabolic processes. Importantly, we found that the circadian rhythm pathway was under positive selection (selection for increased expression) in saline conditions and that the genes involved in circadian rhythm showed significant decoherence between the two environments. These findings make sense in light of the fact that circadian rhythm has been shown to regulate salt tolerance and flowering time in rice (Liang et al., 2021; Wei et al., 2021), so responses to salt stress may alter the expression patterns of genes in the circadian rhythm network, leading to decoherence. In further support of this idea of a link between salt tolerance and circadian rhythm, our selection analyses on physiological traits shows that salt stress leads to selection for earlier flowering. Additionally, our decoherence analyses identified transcripts in important pathways related to plant responses and potentially tolerance to excess salts, including the carbon metabolic pathways (glycolysis and tricarboxylic acid cycle), and accumulation of secondary metabolites and sugar moieties with osmoprotectant properties (Jacoby et al., 2011; Krishnamurthy and Bhagwat, 1990; Petrusa and Winicov, 1997; Tari et al., 2010; Yang and Guo, 2018).
Other genes that are beneficial only in the saline environment include a cyclophilin-encoding transcript (OsCYP2), which has been shown to confer salt tolerance in rice (Ruan et al., 2011). Moreover, we also identified eQTL variants regulating expression of two photosynthesis related antagonistically pleiotropic transcripts (PSAN and CRR7) in the normal environment. Finally, our trans-eQTL analysis identified a hotspot on chromosome 7 that contains OsFLP, a R2R3 myb-like transcription factor that regulates stomatal development (Wu et al., 2019), and appears to be involved in salt tolerance in rice (Qu et al., 2022). Together, these loci are possible new targets for functional and translational studies of salinity stress in rice. Coupled with the insights gained by a systems genomic approach in inferring large scale patterns from genome-wide information, this integration of data across multiple scales across a rice population has allowed us to provide an integrated examination of the molecular and genetic landscape underlying adaptive plant salinity stress responses.
Materials and methods
Plant Material
Domesticated rice is primarily classified into two distinct genetic subgroups, O. sativa ssp. indica and O. sativa ssp. japonica. These subgroups are grown in sympatry and are often recognized as subspecies, given the reproductive barriers between them (Nadir et al., 2018). Further analyses have identified a widely accepted classification consisting of genetically distinct varietal groups namely, indica, aus/circum-aus, aromatic/circum-basmati, tropical japonica, and temperate japonica (Garris et al., 2005). For this study, a total of 130 O. sativa ssp. indica (including indica and circum-aus groups) and 65 O. sativa ssp. japonica (including circum-basmati, tropical japonica, and temperate japonica) accessions were selected, including traditional varieties/landraces and three additionally replicated salt-sensitive and -tolerant test varieties. We focused our analyses on O. sativa ssp. indica since it is the predominant global varietal group (Supplementary Table 12). Seeds were obtained from the International Rice Genebank Collection (IRGC) at the International Rice Research Institute (IRRI) in the Philippines, and from a 2016 bulk seed collection obtained from plants grown under normal (wet and non-saline) conditions at IRRI during the course of a previous study (Groen et al., 2020).
Field Experiment
The field experiment was conducted in the dry season of 2017 at IRRI in Los Baños, Laguna, Philippines. Seeds from each accession were sown on December 16, 2016, and seedlings were then transplanted into the experimental fields at 17 days after sowing (DAS), on January 5, 2017. The field experiment was conducted across two locations: site L4 (14°09’34.6″N 121°15’42.4″E) was prepared as the non-salinized “normal” environment and site L5 (14°09’35.2″N 121°15’42.5″E) as the salinized environment, following Groen et al., (2020). Within each field environment, there were three blocks and three replicates of each genotype (accession), with each genotype planted once per block in a random location. Each plant was planted in a single-row with 0.2-m × 0.2-m spacing between them for a total of one focal plant and seven neighboring plants (included in the experiment) per plot. Each experimental plot included the accessions NSIC Rc 222 and NSIC Rc 182 that served as border rows (Supplemental Fig. S9; Supplemental Table S13). The application of salt in site L5 started on January 19, 2017, when the plants were 31 days old. The salinity level was monitored by recording electrical conductivity (EC), using EC meters installed in each of the parcels at a depth of 30 cm. The EC levels were recorded twice per day until reaching an EC=6 dSm-1 and then recorded daily. The salinity levels were then maintained at 6 dSm-1 (considered mild to moderate salinity stress) until maturity. Management and maintenance of the fields included the application of basal fertilizer, spraying of insecticides against thrips and removal of plants potentially infected with rice tungro virus.
Tissue collection for transcriptome sequencing
Leaf sampling was performed as previously described (Groen et al., 2020). Briefly, leaf collection in the non-saline and saline field was done from 10:00h to 12:00h at 38 DAS (8 days after the beginning of the salt treatment) in the non-saline and saline field from 10:00h to 12:00h. Both fields were sampled simultaneously and with individuals within a block collected in the same order. For each sample, about 10 cm of leaf length were cut into small pieces and placed in chilled 5-mL tubes containing 4 mL of RNALater (Fisher Scientific) solution for RNA stabilization and storage. Leaf samples from each of the 5-ml tubes were then transferred into pairs of 2-mL tubes (one for processing and one for backup), then stored at − 80 °C.
Yield harvesting and panicle trait phenotyping
A total of 780 plants were harvested individually and labeled such that the yield of all plants used for each type of measurement (mRNA sequencing, phenotypic measurements) was known. Individual seeds collected were further categorized as filled, partially filled, and unfilled, using manual assessment and a seed counter (Hoffman Manufacturing). Panicle length measurements and panicle trait phenotyping (PTRAP) of 30 seeds were also performed.
Functional trait phenotyping
In addition to yield-related measurements, we collected data on a number of physiological, morphological, and phenological traits to assess differences between rice accessions in response to soil salinity. In both the non-saline and saline fields, we recorded leaf osmotic potential (LOP) in the vegetative stage, performed chlorophyll analysis based on 1-mg leaf samples and measured ion content (analysis of sodium, Na+, and potassium, K+, analysis) based on 20-mg of leaf samples in both the vegetative and reproductive stages. We also measured plant height for growth rate (measured once a week until maturity). Flowering time was recorded as the day on which 50% of plants in a plot flowered; these plants included the focal plant and its seven neighboring plants. Whole plants were both harvested at the vegetative stage and at maturity to measure wet and dry biomass.
Extraction of total RNA for library construction
Leaf samples stored at −80 °C were thawed at room temperature briefly and excess RNALater was removed. Tissue samples were then flash-frozen in liquid nitrogen and ground using a TissueLyser II (Qiagen). After this, total RNA was extracted using the RNeasy Plant Mini Kit according to manufacturer’s protocol (Qiagen) and eluted in nuclease-free water. The integrity of total RNA was assessed by agarose gel electrophoresis, and RNA from a random subset of samples was further assessed by Agilent TapeStation (Agilent Technologies). RNA concentration was quantified on a Qubit (Invitrogen). Samples were stored at −80 °C until library preparation.
RNA-seq library preparation and sequencing
Library preparation for 780 samples was performed as described in Groen et al. (2020) and followed a plate-based 3′-end mRNA sequencing (3′ mRNA-seq) protocol. Briefly, total RNA from each sample was transferred individually into 96-well plates and normalized to a concentration of 10 ng in 50 μL nuclease-free water. Then, mRNA samples were reverse-transcribed using Superscript II Reverse Transcriptase (Thermo Fisher Scientific) and cDNAs were amplified using the Smart-seq2 protocol (Picelli et al., 2013) with modifications (Groen et al., 2020). This resulted in multiplexed pools of 96 samples each that were used for library preparation with the Nextera XT DNA sample prep kit (Illumina), returning 3′-biased cDNA fragments, similar to the Drop-seq protocol (Macosko et al., 2015). The resulting cDNA libraries were then quantified on an Agilent BioAnalyzer and sequenced at the NYU Genomics Core on an Illumina NextSeq 500 with the configuration HighOutput 1 x 75 base pairs (bp) and the settings: Read 1 of 20 bp (bases 1–12, well barcode; bases 13–20, unique molecular identifier (UMI)) and Read 2 of 50 bp. Raw sequence reads have been submitted to the SRA (BioProject PRJNA1010833).
RNA-seq data processing and data normalization
3′ mRNA-seq read data were processed as previously described in Groen et al. (2020). Briefly, Drop-seq tools v1.12 (https://github.com/broadinstitute/Drop-seq) and Picard tools v2.9.0 (https://broadinstitute.github.io/picard/) were used to generate the metadata. The reference genome, Nipponbare IRGSP 1.0 (ftp://ftp.ncbi.nlm.nih.gov/genomes/all/GCF/001/433/935/GCF_001433935.1_IRGSP-1.0), and annotations were indexed with STAR v020201 (Dobin et al., 2013). Prior to generating read counts, raw reads were converted from FASTQ to unaligned BAM format using Picard tools FastqToSam before being processed using the unified script for a FASTQ starting format. After this, digital gene-expression matrices displaying either UMI or raw read counts with transcripts as rows and samples as columns containing counts from reads with 96 expected sample barcodes were produced using the DigitalExpression utility. Sample barcodes corresponding to beads never exposed to rice total RNA were filtered out based on low numbers of transcribed elements as described previously (Groen et al., 2020; Macosko et al., 2015). Rice individuals that ended up being discarded due to low numbers of transcribed elements, were sequenced again using another library.
UMI counts per sample were normalized through dividing by the total number of detected UMIs in that sample and multiplying by 1 × 106 to obtain transcripts per million. The resulting data matrices were then merged into one digital gene-expression super-matrix, containing transcripts-per-million expression data for all samples. Elements with very low transcription levels (transcript models with a sigma signal < 20) were discarded, after which a robust normalization was conducted using an invariant set normalization protocol within the DChip utility v2010.01 (Li and Wong, 2001). All downstream analyses were done in log-space, using normalized expression levels (log2[normalized transcripts-per-million value + 1]) of transcribed elements estimated using R v3.4.3 (Computing, 2013; Robinson et al., 2010). In a final step of filtering, transcripts that were not detected in at least 10% of individuals across our populations and did not derive from protein-coding genes on nuclear chromosomes were removed prior to performing subsequent analyses.
Quantitative genetics of fecundity and gene expression
The effect of genotype (G), environment (E) and genotype-by-environment (G×E) on fecundity (number of filled grains) was assessed using two-way analysis of variance (ANOVA) with E as a fixed effect, and G and G×E as random effects. The significance of each term was determined using the F-tests (fixed effect) and likelihood ratio test (random effects). Fecundity, averaged by genotype, was further compared between the environments using a two-tailed paired t-test. Variation in gene expression was partitioned similarly using the same model as above, and significance of each term was tested using F-tests via a mixed-model ANOVA (Howell, 1997). Multiple testing was controlled using a False-Discovery Rate (FDR) of 0.001. Broad-sense heritabilities were estimated as H2 = σ2G/(σ2G + σ2GE/e + σ2E/re), with σ2G, σ2E and σ2GEas the variance explained due to G, E, and G×E; e and r represent the number of environments and number of replicates per environment, respectively. Inter-varietal differences in gene expression was estimated for each environmental condition as the population-wide variance between accession mean expression levels.
Univariate and multivariate selection analyses
Univariate selection differentials consider each trait separately and represent the total strength of selection acting on a trait (Lande, 1979; Lande and Arnold, 1983). Fitness was estimated as fecundity, which was the total number of filled rice grains per individual. Unstandardized linear selection differentials (S) on gene expression were estimated as coefficients of linear regression with fecundity as the dependent variable in each regression and transcript abundance and block as the independent variables. Unstandardized quadratic selection differentials (C) were estimated similarly as twice the coefficients of quadratic regression for transcript abundance and fecundity (Conner and Hartl, 2004; Stinchcombe et al., 2008). Using these we then estimated the variance-standardized, and mean-standardized selection differentials (Hereford et al., 2004; Lande and Arnold, 1983). Data preparation included filtering out individuals with zero fecundity followed by normalizing fecundity fitness by mean fitness. Further, transcript abundances were standardized by subtracting population mean abundance of the transcript and dividing by the standard deviation (SD) of the same transcript. To satisfy normality and remove noise inherent in expression data, transcripts with expression values more than 3 standard deviations from the mean were removed, which affected fewer than 1% of individuals. Selection differentials were estimated for transcripts that were expressed in at least 20 individuals and were estimated separately in each of the two environments. Multiple testing was controlled using Bonferroni correction (Bland and Altman, 1995).
Multivariate selection gradients represent the strength of direct selection acting on each trait, after removing indirect selection caused by correlations with other traits (Lande and Arnold, 1983). Principal component analysis (PCA) was performed on transcript abundance using the prcomp function in R (Computing, 2013; Kassambara, 2017) and PCs explaining over 0.5% of variance in each environment were chosen for multivariate selection analyses. Linear (β) and quadratic (³) selection gradients were estimated as coefficients of multiple regression with the normalized fecundity fitness as the dependent variable and the PCs as the independent variables (Conner and Hartl, 2004). We then chose the top one percent of transcripts showing the highest loadings of the PCs (n=182) and counted the number of transcripts showing evidence for positive selection (same directionality between loading of the transcript on the PC and univariate selection differential (S) acting on the transcript estimated above) versus negative selection (opposite directionality between loading of the transcript on the PC and univariate selection differential (S) acting on the transcript estimated above) and based on the majority assigned a directionality to the selection gradients on PCs. Variance-standardized multivariate selection gradients were estimated for functional traits without strong correlation to avoid collinearity among traits (Lande and Arnold, 1983; (Presotto et al., 2019), using a Pearson correlation coefficient < 0.6 as threshold, which was estimated using the cor function in R (Computing, 2013).
Selection analyses on gene ontology biological processes
Gene Ontology (GO) term annotations for rice genes/transcripts were downloaded from Monocots PLAZA 5.0 (Van Bel et al., 2022). All fourth-level biological-process terms were downloaded using GO.db v3.15.0 (Carlson, 2022) and only these terms were considered for further analyses. Next, terms with fewer than 20 transcripts in our dataset were filtered out to minimize redundancy, leaving a total of 670 terms and 10,235 associated transcripts. The selection strength on a biological-process term was estimated as the median selection strength of all transcripts annotated with that term. A term was considered to be under significantly stronger selection compared to the transcriptome-wide median if the median strength of selection for a term was over the transcriptome-wide median selection strength by at least the 95% confidence interval for the selection strength of that term. GO enrichments were done using ShinyGO (Ge et al., 2020).
Regulatory decoherence analyses
To examine regulatory decoherence in rice, a recently developed method, CILP (Correlation by Individual Level Product), was used (Lea et al., 2019). Since CILP calculates product correlations for all possible pairs of genes, only transcripts with a selection strength greater than 0.1 (|S| > 0.1) in at least one environment with expression greater than 0 in at least 50% individuals were included to reduce the dimensionality (leaving 2,318 transcripts). Multiple testing was controlled using a False-Discovery Rate (FDR) of 0.05.
Genotype data and SNP calling
Raw FASTQ files were downloaded from the Sequence Read Archive (SRA) website under BioProject PRJEB6180 (Wang et al., 2018) and under bioprojects PRJNA422249 and PRJNA557122 for 92 accessions (Gutaker et al., 2020) (Supplementary Table 12). Further, genomes of 19 accessions were re-sequenced, and submitted to the SRA (BioProject PRJNA1012700), leading to a total of 125 accessions for which genomic data was available (Supplementary Table 12).
Raw reads were processed for quality control and adapter trimming using the bbduk program of BBTools version 37.66 (https://jgi.doe.gov/data-and-tools/bbtools/) using the options: minlen = 25 qtrim = rl trimq = 10 ktrim = r k = 25 mink = 11 hdist = 1 tpe tbo. The output from this program was mapped to the reference genome O. sativa Nipponbare IRGSP 1.0 genome that was downloaded from NCBI Genome (https://www.ncbi.nlm.nih.gov/genome/?term=txid4530[orgn]) using bwa-mem2 v2.1 (Vasimuddin et al., 2019). PCR duplicates were marked and removed using the Picard tools version 2.9.0. SNPs were called using GATK HaplotypeCaller v4.2.0.0 to obtain a multi-accession joint SNP file. Only SNPs that were above 5bp distance from an indel variant were taken. Next, SNPs were filtered using the recommended GATK hard filtering (Van der Auwera et al., 2013). Further, using vcftools v0.1.16 (Danecek et al., 2011), SNPs with at least 80% genotype calls and a minor allele frequency of 0.05 were retained (--max-missing 0.8 --maf 0.05). Since rice is a inbred species, we also removed any SNPs that displayed heterozygosity of over 5% identified using vcftools v0.1.16 –hardy (Danecek et al., 2011). Next, missing genotype calls were imputed and phased using Beagle v4.1 (Browning and Browning, 2016), and using vcftools v0.1.16 -m2 -M2 (Danecek et al., 2011) only biallelic SNPs were retained for further analyses. Finally, SNPs were randomly pruned such that one SNP per 1000bp was retained using vcftools v0.1.16 –thin (Danecek et al., 2011), leaving a SNP dataset of 246,714 markers.
G-matrix estimation and prediction of short-term phenotypic evolution
A G-matrix (G) representing the additive genetic variance and covariance was estimated for the principal component axes (PCs) by taking the eigengene. This was done by deploying GREML v1.94 (Yang et al., 2011). Although the principal components are by definition uncorrelated at the level of the individual replicate plants, they start showing genetic covariances when loading values of replicates from each genotype are averaged. Next, using the multivariate breeder’s equation (Δ z = G β), we predicted the response to short-term phenotypic selection (Δz) on the PCs.
Association mapping
Association mapping was performed between the SNP markers and gene expression values recorded in the normal and saline environments. For this, the linear model in Matrix eQTL was used (Shabalin, 2012). The normalized gene expression values were averaged over the replicates in each environment separately and these averages were subsequently used to test for associations. The first five principal components (PCs) of the kinship matrix were estimated using GAPIT v3 (Wang and Zhang, 2021) and added as covariates to control for population structure. Associations were considered significant at a false-discovery rate (FDR) < 0.001, and when significant were included in downstream analyses. Due to long stretches of homozygosity attributed to the highly inbred nature of rice, trans-eQTLs were defined as being on a different chromosome or at least 1 Mb away from a gene under its influence on the same chromosome; cis-eQTLs were defined as <100 kb away from an associated gene. To identify significant G × eQTLs, we ran Matrix eQTL on the difference of expression in the normal and saline conditions (normal – saline) at FDR 0.05 (Smith and Kruglyak, 2008).
For each gene regulated by a SNP in cis or trans, a lead SNP was identified as the SNP with the most significant association within a 100-kb region. Furthermore, trans-eQTL hotspots were identified through analyzing the number of unique genes regulated by lead SNPs in a given 100-kb region. To detect genes under cis-trans compensation or reinforcement, the effect sizes of all lead trans-eQTLs were averaged for each gene and compared to the lead cis-eQTL for the same gene. Next, for these genes we estimated the mean proportion of individuals with opposite and same direction cis-trans allelic configuration. Genes were defined as compensating and reinforcing if they had at least 60% of individuals with opposite and same cis-trans allelic configuration, respectively. To examine the LD structure, we estimated r2 using plink v1.9 (Purcell et al., 2007) between (1) all pairs of cis- and trans-variants for the identified compensating and reinforcing genes and (2) 1,000 datasets of randomly selected SNP pairs (with equal numbers of variant pairs and a similar distribution of distances between the cis- and trans-variants as in 1). Next, we compared r2 between (1) and (2) using a two-tailed permutation test.
To examine the patterns of past selection on eQTLs, we used the minor allele frequency (MAF) of the 246,714 SNP markers and compared these using the t.test function in R (Computing, 2013). Furthermore, we estimated the site-wise nucleotide diversity (π) and averaged it over 50-kb flanking regions around each cis-eQTLs (100-kb region total). We compared this π to the background nucleotide diversity, estimated as π averaged over 100-kb blocks throughout the genome minus the 100-kb cis-eQTL region above. The difference in mean of nucleotide diversity was tested using the t.test function in R (Computing, 2013).
Acknowledgements
We thank the New York University Center for Genomics and Systems Biology GenCore Facility for sequencing support, and the New York University High Performance Computing for supplying computational resources. We would also like to thank William Mauck and Nicholas Rogers for technical support, as well as Adrian E. Platts and Jae Young Choi, for their help with post-sequencing data processing. We are grateful to current members of the Purugganan laboratory (particularly J. Flowers, A. Kurbidaeva, O. Alam) and Elena Hamann at Fordham University (currently Heinrich Heine University Düsseldorf) for insightful discussions. S.C.G. is supported by a grant from the Gordon and Betty Moore Foundation/Life Sciences Research Foundation (award no. GBMF2550.06 to S.C.G.), startup funds from the University of California Riverside, and a grant from the National Institute of General Medical Sciences of the National Institutes of Health (award no. R35GM151194 to S.C.G.). This work was funded in part by grants from the US National Science Foundation Plant Genome Research Program (IOS 1546218 and 2204374) and the Zegar Family Foundation.
References
- Environmental duress and epistasis: how does stress affect the strength of selection on new mutations?Trends Ecol Evol 25:450–458https://doi.org/10.1016/j.tree.2010.05.003
- The strength and form of natural selection on transcript abundance in the wildMol Ecol 30:2724–2737https://doi.org/10.1111/mec.15743
- Transcriptome and epigenome landscape of human cortical development modeled in organoidsScience 362https://doi.org/10.1126/science.aat6720
- Genetic trade-offs and conditional neutrality contribute to local adaptationMol Ecol 22:699–708https://doi.org/10.1111/j.1365-294X.2012.05522.x
- Evolutionary genetics of plant adaptationTrends Genet 27:258–266https://doi.org/10.1016/j.tig.2011.04.001
- Environmental stress does not increase the mean strength of selectionJ Evol Biol 31:1599–1606https://doi.org/10.1111/jeb.13351
- Effect of salinity on growth and chlorophyll content in ricePakistan journal of scientific and industrial research 43:130–132
- Multiple significance tests: the Bonferroni methodBMJ 310https://doi.org/10.1136/bmj.310.6973.170
- Genotype Imputation with Millions of Reference SamplesAm J Hum Genet 98:116–126https://doi.org/10.1016/j.ajhg.2015.11.020
- The influence of genetic architecture on responses to selection under drought in riceEvol Appl 15:1670–1690https://doi.org/10.1111/eva.13419
- Salinity stress and salt tolerance In: Shanker A, Venkateswarlu B, editors. Abiotic Stress in Plants - Mechanisms and Adaptations. LondonEngland: InTech https://doi.org/10.5772/22331
- GO. db: A set of annotation maps describing the entire Gene Ontologyhttps://doi.org/10.18129/B9.bioc.GO.db
- Response to Salinity in Rice: Comparative Effects of Osmotic and Ionic StressesPlant Prod Sci 10:159–170https://doi.org/10.1626/pps.10.159
- Morphological and metabolic responses to salt stress of rice (Oryza sativa L.) cultivars which differ in salinity tolerancePlant Physiol Biochem 144:427–435https://doi.org/10.1016/j.plaphy.2019.10.017
- Connecting salt stress signalling pathways with salinity-induced changes in mitochondrial metabolic processes in C3 plantsPlant Cell Environ 40:2875–2905https://doi.org/10.1111/pce.13034
- R: A language and environment for statistical computingVienna: R Core Team
- A primer of ecological genetics
- Tempo and mode of regulatory evolution in DrosophilaGenome Res 24:797–808https://doi.org/10.1101/gr.163014.113
- Effects of abiotic stress on plants: a systems biology perspectiveBMC Plant Biol 11https://doi.org/10.1186/1471-2229-11-163
- CYCLIC NUCLEOTIDE-GATED ION CHANNELs 14 and 16 Promote Tolerance to Heat and Chilling in RicePlant Physiol 183:1794–1808https://doi.org/10.1104/pp.20.00591
- The variant call format and VCFtoolsBioinformatics 27:2156–2158https://doi.org/10.1093/bioinformatics/btr330
- Integrative inference of transcriptional networks in Arabidopsis yields novel ROS signalling regulatorsNat Plants 7:500–513https://doi.org/10.1038/s41477-021-00894-1
- STAR: ultrafast universal RNA-seq alignerBioinformatics 29:15–21https://doi.org/10.1093/bioinformatics/bts635
- Natural selection on cis and trans regulation in yeastsGenome Res 20:826–836https://doi.org/10.1101/gr.101576.109
- Crop Production under Drought and Heat Stress: Plant Responses and Management OptionsFront Plant Sci 8https://doi.org/10.3389/fpls.2017.01147
- Introduction to quantitative genetics
- Abiotic Stress Effects on Performance of Horticultural CropsHorticulturae 5https://doi.org/10.3390/horticulturae5040067
- Genetic structure and diversity in Oryza sativa LGenetics 169:1631–1638https://doi.org/10.1534/genetics.104.035642
- ShinyGO: a graphical gene-set enrichment tool for animals and plantsBioinformatics 36:2628–2629https://doi.org/10.1093/bioinformatics/btz931
- Extensive compensatory cis-trans regulation in the evolution of mouse gene expressionGenome Res 22:2376–2384https://doi.org/10.1101/gr.142281.112
- The strength and pattern of natural selection on gene expression in riceNature 578:572–576https://doi.org/10.1038/s41586-020-1997-2
- Fitness costs and benefits of gene expression plasticity in rice under droughtbioRxiv https://doi.org/10.1101/2021.03.16.435597
- Genomic history and ecology of the geographic spread of riceNat Plants 6:492–502https://doi.org/10.1038/s41477-020-0659-6
- Biochemical and anatomical changes and yield reduction in rice (Oryza sativa L.) under varied salinity regimesBiomed Res Int 2014https://doi.org/10.1155/2014/208584
- Regulatory hotspots are associated with plant gene expression under varying soil phosphorus supply in Brassica rapaPlant Physiol 156:1230–1241https://doi.org/10.1104/pp.111.175612
- Comparing strengths of directional selection: how strong is strong?Evolution 58:2133–2143https://doi.org/10.1111/j.0014-3820.2004.tb01592.x
- Ultrarare variants drive substantial cis heritability of human gene expressionNat Genet 51:1349–1355https://doi.org/10.1038/s41588-019-0487-7
- Strength and tempo of directional selection in the wildProc Natl Acad Sci U S A 98:9157–9160https://doi.org/10.1073/pnas.161281098
- Statistical Methods for PsychologyCengage Learning
- Effects of salt stress on rice growth, development characteristics, and the regulating ways: A reviewJ Integr Agric 16:2357–2374https://doi.org/10.1016/S2095-3119(16)61608-8
- Comparative Developmental Transcriptomics Reveals Rewiring of a Highly Conserved Gene Regulatory Network during a Major Life History Switch in the Sea Urchin Genus HeliocidarisPLoS Biol 14https://doi.org/10.1371/journal.pbio.1002391
- The role of mitochondrial respiration in salinity toleranceTrends Plant Sci 16:614–623https://doi.org/10.1016/j.tplants.2011.08.002
- Interactions between stressful environment and gene deletions alleviate the expected average loss of fitness in yeastGenetics 178:2105–2111https://doi.org/10.1534/genetics.107.084533
- Plant cell-surface GIPC sphingolipids sense salt to trigger Ca2+ influxNature 572:341–346https://doi.org/10.1038/s41586-019-1449-z
- An inferred fitness consequence map of the rice genomeNat Plants 6:119–130https://doi.org/10.1038/s41477-019-0589-3
- Association mapping reveals the role of purifying selection in the maintenance of genomic variation in gene expressionProceedings of the National Academy of Sciences 112:15390–15395https://doi.org/10.1073/pnas.1503027112
- Practical Guide To Principal Component Methods in R: PCA, M(CA), FAMD, MFA, HCPC, factoextra
- Release of SOS2 kinase from sequestration with GIGANTEA determines salt tolerance in ArabidopsisNat Commun 4https://doi.org/10.1038/ncomms2357
- The strength of phenotypic selection in natural populationsAm Nat 157:245–261https://doi.org/10.1086/319193
- Patterns and Power of Phenotypic Selection in NatureBioscience 57:561–572https://doi.org/10.1641/B570706
- Quantitative genomics: analyzing intraspecific variation using global gene expression polymorphisms or eQTLsAnnu Rev Plant Biol 60:93–114https://doi.org/10.1146/annurev.arplant.043008.092114
- Network-based approaches for understanding gene regulation and function in plantsPlant J 104:302–317https://doi.org/10.1111/tpj.14940
- Genotype—environment interactions and the estimation of the genomic mutation rate in Drosophila melanogasterProceedings of the Royal Society of London Series B: Biological Sciences 258:221–227https://doi.org/10.1098/rspb.1994.0166
- Abiotic Stress in Crop ProductionInt J Mol Sci 24https://doi.org/10.3390/ijms24076603
- Accumulation of choline and glycinebetaine in salt-stressed wheat seedlingsCurr Sci 59:111–112
- H2A.Z-containing nucleosomes mediate the thermosensory response in ArabidopsisCell 140:136–147https://doi.org/10.1016/j.cell.2009.11.006
- eQTLs Regulating Transcript Variations Associated with Rapid Internode Elongation in Deepwater RiceFront Plant Sci 8https://doi.org/10.3389/fpls.2017.01753
- Quantitative Genetic Analysis of Multivariate Evolution, Applied to Brain: Body Size AllometryEvolution 33:402–416https://doi.org/10.2307/2407630
- The measurement of selection correlated charactersEvolution 37:1210–1226https://doi.org/10.1111/j.1558-5646.1983.tb00236.x
- Compensatory cis-trans evolution and the dysregulation of gene expression in interspecific hybrids of DrosophilaGenetics 171:1813–1822https://doi.org/10.1534/genetics.105.047449
- Salinity-induced calcium signaling and root adaptation in Arabidopsis require the calcium regulatory protein annexin1Plant Physiol 163:253–262https://doi.org/10.1104/pp.113.217810
- Genetic and environmental perturbations lead to regulatory decoherenceElife 8https://doi.org/10.7554/eLife.40538
- OsCYP21-4, a novel Golgi-resident cyclophilin, increases oxidative stress tolerance in riceFront Plant Sci 6https://doi.org/10.3389/fpls.2015.00797
- The transcriptional repressor OsPRR73 links circadian clock and photoperiod pathway to control heading date in ricePlant Cell Environ 44:842–855https://doi.org/10.1111/pce.13987
- Plant salt-tolerance mechanism: A reviewBiochem Biophys Res Commun 495:286–291https://doi.org/10.1016/j.bbrc.2017.11.043
- Model-based analysis of oligonucleotide arrays: expression index computation and outlier detectionProc Natl Acad Sci U S A 98:31–36https://doi.org/10.1073/pnas.98.1.31
- GA signaling and CO/FT regulatory module mediate salt-induced late flowering in Arabidopsis thalianaPlant Growth Regul 53:195–206https://doi.org/10.1007/s10725-007-9218-7
- Salt tolerance in rice: Physiological responses and molecular mechanismsThe Crop Journal 10:13–25https://doi.org/10.1016/j.cj.2021.02.010
- Mapping regulatory variants controlling gene expression in drought response and tolerance in maizeGenome Biol 21https://doi.org/10.1186/s13059-020-02069-1
- The genomic landscape of molecular responses to natural drought stress in Panicum halliiNat Commun 9https://doi.org/10.1038/s41467-018-07669-x
- Highly Parallel Genome-wide Expression Profiling of Individual Cells Using Nanoliter DropletsCell 161:1202–1214https://doi.org/10.1016/j.cell.2015.05.002
- Environmental Stress and PlantsInt J Mol Sci 23https://doi.org/10.3390/ijms23105416
- Standardizing Selection Strengths to Study Selection in the Wild: A Critical Comparison and Suggestions for the FutureBioscience 62:1039–1054https://doi.org/10.1525/bio.2012.62.12.6
- CYCLIN-DEPENDENT KINASE G2 regulates salinity stress response and salt mediated flowering in Arabidopsis thalianaPlant Mol Biol 88:287–299https://doi.org/10.1007/s11103-015-0324-z
- Ribosome profiling reveals post-transcriptional buffering of divergent gene expression in yeastGenome Res 24:422–430https://doi.org/10.1101/gr.164996.113
- Salt-Tolerant Crops: Time to DeliverAnnu Rev Plant Biol 74:671–696https://doi.org/10.1146/annurev-arplant-061422-104322
- Contrasting Frequencies and Effects of cis- and trans-Regulatory Mutations Affecting Gene ExpressionMol Biol Evol 33:1131–1146https://doi.org/10.1093/molbev/msw011
- Evolution of crop species: Genetics of domestication and diversificationNature Reviews Genetics https://doi.org/10.1038/nrg3605
- Reactive oxygen species homeostasis and signalling during drought and salinity stressesPlant Cell Environ 33:453–467https://doi.org/10.1111/j.1365-3040.2009.02041.x
- Comparative physiology of salt and water stressPlant Cell Environ 25:239–250https://doi.org/10.1046/j.0016-8025.2001.00808.x
- Mechanisms of salinity toleranceAnnu Rev Plant Biol 59:651–681https://doi.org/10.1146/annurev.arplant.59.032607.092911
- An overview on reproductive isolation in Oryza sativa complexAoB Plants 10https://doi.org/10.1093/aobpla/ply060
- Proline status in salt-tolerant and salt-sensitive alfalfa cell lines and plants in response to NaClPlant Physiol Biochem
- Smart-seq2 for sensitive full-length transcriptome profiling in single cellsNat Methods 10:1096–1098https://doi.org/10.1038/nmeth.2639
- Advances in Sensing, Response and Regulation Mechanism of Salt Tolerance in RiceInt J Mol Sci 22https://doi.org/10.3390/ijms22052254
- Phenotypic selection under two contrasting environments in wild sunflower and its crop-wild hybridEvol Appl 12:1703–1717https://doi.org/10.1111/eva.12828
- A Novel Strategy to Identify Prognosis-Relevant Gene Sets in CancersGenes 13https://doi.org/10.3390/genes13050862
- PLINK: a tool set for whole-genome association and population-based linkage analysesAm J Hum Genet 81:559–575https://doi.org/10.1086/519795
- The phytohormonal regulation of Na+/K+ and reactive oxygen species homeostasis in rice salt responseMol Breed 40https://doi.org/10.1007/s11032-020-1100-6
- Lead Modulates trans- and cis-Expression Quantitative Trait Loci (eQTLs) in Drosophila melanogaster HeadsFront Genet 9https://doi.org/10.3389/fgene.2018.00395
- A Rice R2R3-Type MYB Transcription Factor OsFLP Positively Regulates Drought Stress Response via OsNACInt J Mol Sci 23https://doi.org/10.3390/ijms23115873
- Describing the physiological responses of different rice genotypes to salt stress using sigmoid and piecewise linear functionsField Crops Res 220:46–56https://doi.org/10.1016/j.fcr.2017.05.001
- edgeR: a Bioconductor package for differential expression analysis of digital gene expression dataBioinformatics 26:139–140https://doi.org/10.1093/bioinformatics/btp616
- OsCyp2-P, an auxin-responsive cyclophilin, regulates Ca2+ calmodulin interaction for an ion-mediated stress response in ricePhysiol Plant 174https://doi.org/10.1111/ppl.13631
- Proteomic identification of OsCYP2, a rice cyclophilin that confers salt tolerance in rice (Oryza sativa L.) seedlings when overexpressedBMC Plant Biol 11https://doi.org/10.1186/1471-2229-11-34
- Conflicts in natural selection constrain adaptation to climate change in Arabidopsis thalianabioRxiv https://doi.org/10.1101/2023.10.16.562583
- Effect of water stress on cell division and cell-division-cycle 2-like cell-cycle kinase activity in wheat leavesPlant Physiol 117:667–678https://doi.org/10.1104/pp.117.2.667
- Matrix eQTL: ultra fast eQTL analysis via large matrix operationsBioinformatics 28:1353–1358https://doi.org/10.1093/bioinformatics/bts163
- Genetics and genomics of root system variation in adaptation to drought stress in cereal cropsJ Exp Bot 72:1007–1019https://doi.org/10.1093/jxb/eraa487
- The Evolution of Gene Expression in cis and transTrends Genet 34:532–544https://doi.org/10.1016/j.tig.2018.03.007
- Gene-environment interaction in yeast gene expressionPLoS Biol 6https://doi.org/10.1371/journal.pbio.0060083
- Genetical Genomics Reveals Large Scale Genotype-By-Environment Interactions in Arabidopsis thalianaFront Genet 3https://doi.org/10.3389/fgene.2012.00317
- Plasticity of maternal environment-dependent expression-QTLs of tomato seedsTheor Appl Genet 136https://doi.org/10.1007/s00122-023-04322-0
- Estimating nonlinear selection gradients using quadratic regression coefficients: double or nothing?Evolution 62:2435–2440https://doi.org/10.1111/j.1558-5646.2008.00449.x
- Effect of salt stress on growth, chlorophyll content, lipid peroxidation and antioxidant defence systems in Phaseolus vulgaris LS Afr J Bot 105:306–312https://doi.org/10.1016/j.sajb.2016.03.011
- Salicylic acid increased aldose reductase activity and sorbitol accumulation in tomato plants under salt stressBiol Plant 54:677–683https://doi.org/10.1007/s10535-010-0120-1
- Chlorophyll fluorescence analysis in diverse rice varieties reveals the positive correlation between the seedlings salt tolerance and photosynthetic efficiencyBMC Plant Biol 19https://doi.org/10.1186/s12870-019-1983-8
- PLAZA 5.0: extending the scope and power of comparative and functional genomics in plantsNucleic Acids Res 50:D1468–D1474https://doi.org/10.1093/nar/gkab1024
- From FastQ data to high confidence variant calls: the Genome Analysis Toolkit best practices pipelineCurr Protoc Bioinformatics 43:11–11https://doi.org/10.1002/0471250953.bi1110s43
- Delayed and carryover effects of salinity on flowering in Iris hexagona (Iridaceae)Am J Bot 89:1847–1851https://doi.org/10.3732/ajb.89.11.1847
- Efficient Architecture-Aware Acceleration of BWA-MEM for Multicore Systems2019 IEEE International Parallel and Distributed Processing Symposium (IPDPS) :314–324https://doi.org/10.1109/IPDPS.2019.00041
- The gene regulatory logic of transcription factor evolutionTrends Ecol Evol 23:377–385https://doi.org/10.1016/j.tree.2008.03.006
- Effects of salt stress on ion balance and nitrogen metabolism of old and young leaves in rice (Oryza sativa L.)BMC Plant Biol 12https://doi.org/10.1186/1471-2229-12-194
- GAPIT Version 3: Boosting Power and Accuracy for Genomic Association and PredictionGenomics Proteomics Bioinformatics 19:629–640https://doi.org/10.1016/j.gpb.2021.08.005
- Genomic variation in 3,010 diverse accessions of Asian cultivated riceNature 557:43–49https://doi.org/10.1038/s41586-018-0063-9
- Natural disaster and immunological aging in a nonhuman primateProc Natl Acad Sci U S A 119https://doi.org/10.1073/pnas.2121663119
- Clock component OsPRR73 positively regulates rice salt tolerance by modulating OsHKT2;1-mediated sodium homeostasisEMBO J 40https://doi.org/10.15252/embj.2020105086
- Cell cycle modulation in the response of the primary root of Arabidopsis to salt stressPlant Physiol 135:1050–1058https://doi.org/10.1104/pp.104.040022
- Global eQTL mapping reveals the complex genetic architecture of transcript-level variation in ArabidopsisGenetics 175:1441–1450https://doi.org/10.1534/genetics.106.064972
- EGRINs (Environmental Gene Regulatory Influence Networks) in Rice That Function in the Response to Water Deficit, High Temperature, and Agricultural EnvironmentsPlant Cell 28:2365–2384https://doi.org/10.1105/tpc.16.00158
- Regulatory changes underlying expression differences within and between Drosophila speciesNat Genet 40:346–350https://doi.org/10.1038/ng.77
- Multiple transcriptional factors control stomata development in riceNew Phytol 223:220–232https://doi.org/10.1111/nph.15766
- The circadian clock ticks in plant stress responsesStress Biol 2https://doi.org/10.1007/s44154-022-00040-7
- GCTA: a tool for genome-wide complex trait analysisAm J Hum Genet 88:76–82https://doi.org/10.1016/j.ajhg.2010.11.011
- Unraveling salt stress signaling in plantsJ Integr Plant Biol 60:796–804https://doi.org/10.1111/jipb.12689
- OSCA1 mediates osmotic-stress-evoked Ca2+ increases vital for osmosensing in ArabidopsisNature 514:367–371https://doi.org/10.1038/nature13593
- Abiotic stress responses in plantsNat Rev Genet 23:104–119https://doi.org/10.1038/s41576-021-00413-0
- Study on the Effect of Salt Stress on Yield and Grain Quality Among Different Rice VarietiesFront Plant Sci 13https://doi.org/10.3389/fpls.2022.918460
Article and author information
Version history
- Sent for peer review:
- Preprint posted:
- Reviewed Preprint version 1:
Copyright
© 2024, Gupta et al.
This article is distributed under the terms of the Creative Commons Attribution License, which permits unrestricted use and redistribution provided that the original author and source are credited.