Genome-wide analyses of direct target genes of four rice NAC-domain transcription factors involved in drought tolerance
BMC Genomics volume 19, Article number: 40 (2018)
Plant stress responses and mechanisms determining tolerance are controlled by diverse sets of genes. Transcription factors (TFs) have been implicated in conferring drought tolerance under drought stress conditions, and the identification of their target genes can elucidate molecular regulatory networks that orchestrate tolerance mechanisms.
We generated transgenic rice plants overexpressing the 4 rice TFs, OsNAC5, 6, 9, and 10, under the control of the root-specific RCc3 promoter. We showed that they were tolerant to drought stress with reduced loss of grain yield under drought conditions compared with wild type plants. To understand the molecular mechanisms underlying this tolerance, we here performed chromatin immunoprecipitation (ChIP)-Seq and RNA-Seq analyses to identify the direct target genes of the OsNAC proteins using the RCc3:6MYC-OsNAC expressing roots. A total of 475 binding loci for the 4 OsNAC proteins were identified by cross-referencing their binding to promoter regions and the expression levels of the corresponding genes. The binding loci were distributed among the promoter regions of 391 target genes that were directly up-regulated by one of the OsNAC proteins in four RCc3:6MYC-OsNAC transgenic lines. Based on gene ontology (GO) analysis, the direct target genes were related to transmembrane/transporter activity, vesicle, plant hormones, carbohydrate metabolism, and TFs. The direct targets of each OsNAC range from 4.0–8.7% of the total number of up-regulated genes found in the RNA-Seq data sets. Thus, each OsNAC up-regulates a set of direct target genes that alter root system architecture in the RCc3:OsNAC plants to confer drought tolerance. Our results provide a valuable resource for functional dissection of the molecular mechanisms of drought tolerance.
Many of the target genes, including transmembrane/transporter, vesicle related, auxin/hormone related, carbohydrate metabolic processes, and transcription factor genes, that are up-regulated by OsNACs act as the cellular components which would alter the root architectures of RCc3:OsNACs for drought tolerance.
The global population is approximately 7.2 billion, but is projected to reach 9.7 billion by 2050 . In order to feed this growing population, it is estimated that current the annual crop production of 2 billion tons will need to double within this time frame. Plants are constantly exposed to environmental stresses, such as drought, salt, flooding, and high temperatures, which cause severe losses of crop yield. Drought stress, in particular, adversely affects plant growth, membrane integrity, osmotic adjustments, water relations, photosynthetic activity and grain yield . New innovative approaches are need to develop genetically engineered crops that are tolerant to drought stress conditions. Accordingly, an understanding of the molecular and cellular responses to drought stress is a prerequisite for the identification of key regulators that can be used to produce such crops [3,4,5].
Molecular and genomic analyses have revealed many drought-inducible transcription factors (TFs) that regulate the expression of stress-inducible downstream genes [6, 7]. These include TFs belonging to the AP2/ERF, MYB, bZIP and NAC families. TFs are regulators that bind to cis-regulatory sequences in the promoter and/or enhancer elements of target genes to activate or inactivate their expression . NAC [no apical meristem (NAM), Arabidopsis thaliana transcription activation factor (ATAF1/2) and cup shaped cotyledon (CUC2)] domain containing genes constitutes one of the largest plant-specific TF families in both A. thaliana and rice [9, 10]. NAC TFs have been shown to play a role in diverse processes involving shoot apical meristem development, embryogenesis, lateral root development [11,12,13], plant defense mechanisms, senescence , and abiotic stress responses . It is estimated that 2408 genes encode TF or TF-like proteins in rice (Oryza sativa), based on the Plant Transcription Factor Database v3.0 (http://planttfdb.cbi.pku.edu.cn). Among these, rice has 170 putative OsNAC genes, based on sequence analysis (Plant Transcription Factor Database, http://planttfdb.cbi.pku.edu.cn/), of which 41 and 29 are induced and repressed, respectively, in leaves of reproductive stage under drought stress conditions . However, only a few of these genes have been functionally characterized in the context of drought tolerance, although it is known that they are often functionally redundant [17, 18]. Despite their important function in modulating drought tolerance, little is known about the molecular mechanisms associated with OsNAC TFs and their target genes .
In A. thaliana, some of the NAC TF direct targets have been identified . Ohashi-Ito et al. (2010) reported that VND6 (VASCULAR-RELATED NAC-DOMAIN6) is a direct regulator of XCP1 (XYLEM CYSTEINE PROTEASE1), and that IRX5/CESA4 (IRREGULARXYLEM5/CELLULOSE SYNTHASE) is involved in plant cell death and secondary wall formation . In addition, the A. thaliana NAC TF, VNI2 (VND INTERACTING2) regulates leaf longevity under environmental stress conditions through direct binding to the promoters of the COR (COLD-REGULATED) and RD (RESPONSIVE TO DEHYDRATION) genes .
Chromatin immunoprecipitation (ChIP) can be used as a genome-wide approach to verify in vivo interactions between TF proteins and genomic DNA under physiological conditions, and to estimate the density of TF binding at specific loci . Recently, both RNA-Seq and ChIP-Seq technologies have been widely used to measure the transcript levels of TFs and to obtain information about genome-wide binding sites, respectively . Such approaches can provide insights into changes in protein association with regulatory regions of genes in response to environmental clues, and can help connect those changes with the transcription patterns of target genes.
Previously, we demonstrated that root-specific overexpression of each of four TF genes, OsNAC5, 6, 9, and 10, significantly enhanced drought tolerance, leading to high grain yield under field-drought conditions [25, 29,30,31]. The common phenotypic characteristics of each OsNAC overexpressing transgenic rice showed a change in root architecture. In particular, enlarged stele and aerenchyma (OsNAC5, OsNAC6 and OsNAC9), and enlarged stele, cortex, and epidermis (OsNAC10) can explain the altered root architecture, including root number and root diameter, and ultimately conferring increased drought tolerance phenotype. In this study, we used genome-wide ChIP-Seq and RNA-Seq analyses to identify the direct target genes for four OsNAC TFs that have been reported to be involved in influencing root growth, leading to drought tolerance. Our results reveal the regulatory networks governed by OsNAC TFs that contribute to drought-tolerance mechanisms, providing a foundation for the development of drought tolerant crops.
Plasmid construction and rice transformation
To generate RCc3:6MYC-OsNAC rice lines, the coding sequence of OsNAC genes were amplified using PrimeSTAR DNA polymerase (TaKaRa, Tokyo, Japan) with the respective gene-specific forward and reverse primers. After enzyme digestion of the PCR products with BamHI and NotI, the DNA fragment containing the coding sequence was then ligated into the multiple cloning site of the pE3n vector, which is flanked with a 6MYC tag coding sequence . Finally, the 6MYC-OsNAC sequences from the pE3n-OsNAC was cloned into the p700-RCc3 vector carrying a 1.3 kb RCc3 (Os02g0662000) promoter sequence using the Gateway system (Invitrogen, Carlsbad, CA) , and the final vectors were named RCc3:6MYC-OsNAC. Transgenic plants were obtained by Agrobacterium tumefaciens (strain LBA4404)-mediated transformation of embryogenic rice callus (Oryza sativa L. japonica cv. Illmi), as previously described .
Plant growth and sampling
RCc3:6MYC-OsNAC transgenic and non-transgenic (NT; Oryza sativa L. japonica cv. Illmi) seeds were germinated on MS (Murashige and Skoog) medium in the dark at 28 °C for 3 d (days) and transferred into light conditions for a half day. Seedlings were transplanted into the soil pots (4 × 4 × 6 cm; 3 plants per pot) and grown in a greenhouse (30 °C for 16 h for day/ 25 °C for 8 h for night cycle) under 60–80% humidity and 200-500 umol m-2 s-1 light intensity (37°32΄51.3”N 128°26΄26.6″E). Two weeks after transplanting to soil, the seedlings were completely submerged in water in the greenhouse for 3 d, the roots were harvested, immediately frozen in liquid nitrogen and stored at −80 °C. The root tissues for RNA-Seq and ChIP-Seq were simultaneously prepared as independent biological replicates.
Confirmation of expression of the OsNAC fused 6MYC protein using anti-MYC antibody
Proteins were extracted from the root samples with extraction buffer [50 mM Tris·HCl (pH 8.0), 150 mM NaCl, 5 mM EDTA, and 0.2% Triton X-100], and protein concentrations were determined using the Bradford method (Bio-Rad Laboratories, Inc., http://www.bio-rad.com). Protein extracts were then separated on 10% SDS-polyacrylamide gels and blotted onto a polyvinylidene difluoride (PVDF) membrane (Immobilon-P; Millipore Corporation, http://www.millipore.com). The immunoreactive proteins were detected using primary antibodies against c-MYC (#sc-789; Santa Cruz Biotechnology), and detected by chemiluminescence with Pierce Super Signal Substrate (Pierce, http://www.piercenet.com) in accordance with the manufacturer’s protocol.
Total RNA was prepared from the root tissue of the RCc3:6MYC-OsNAC transgenic and NT rice plants using Trizol reagent (Invitrogen, Carlsbad, CA) and then purified with an RNeasy Mini Kit (Qiagen, Valencia, CA). Contaminating genomic DNA was removed by treating with DNase I (Invitrogen), according to the manufacturer’s instructions. Handling of the RNA-Seq library construction, next-generation sequencing (NGS) and DEG (Differential Expressed Gene) analysis was completed by MACROGEN Inc., Seoul, Korea (http://macrogen.com). Raw sequence reads were trimmed to remove adaptor sequences, contaminant DNA, and PCR duplicates, and those with a quality lower than Q20 were removed. All reads were assembled with TopHat software, using annotated genes from the rapdb database [http://rapdb.dna.affrc.go.jp; IRGSP (v 1.0)]. DEGs were defined by an expression change ≥2-fold with a P value <0.05. The entire sets of original RNA-Seq data have been deposited in the NCBI Gene Expression Omnibus under GEO Series number GSE102919.
Quantitative real-time PCR qRT-PCR
To validate the results of the RNA-Seq analysis, total 23 genes, representing up-regulated, unchanged, and down-regulated, were selected from the DEGs (differentially expressed genes) in the RNA-Seq analysis (P value < 0.05). RNA isolated from the roots of NT, OsNAC5, 6, 9, and 10 transgenic rice lines was used for cDNA synthesis and qRT-PCR, using gene specific primers. The primers were designed to amplify 120–200 bp amplicons with a Tm at 60 °C. First-strand cDNA was synthesized from 1.0 μg of total RNA primed with oligo-dT using Superscript III (Invitrogen). qRT-PCR was performed using a Stratagene Mx3000P DNA analyzer (Agilent Technologies) with the Solg™ 2x real-time PCR smart mix with evagreen (Solgent, Seoul, Korea) for signal detection. To normalize the total amount of cDNA present in each reaction, a housekeeping gene encoding ubiquitin (Os06g0681400) was co-amplified. Primer sequences are listed in Additional file 1: Table S8. Values are the means ± SD (standard deviation) of three independnet experiments.
Chromatin immunoprecipitation (ChIP) and ChIP-Seq
Nuclei-enriched fractions were partially purified from rice roots as previously described with slight modifications . RCc3:OsNACs and NT rice seedlings grown in soil were hydroponically adapted in water for 3 d and cross-linked with 1% formaldehyde by vacuum infiltration for 15 min. Cross-linking was stopped by the addition of a final concentration of 125 mM glycine and vacuum infiltrating for 10 min. After washing the plants in cold water, the roots were excised, frozen in liquid nitrogen and stored at −80 °C.
Chromatin from the nuclei-enriched fractions was fragmented to an average size of approximately 100–200 bp by 15 cycles of sonication (30 s each) in 15 mL falcon tubes using a Bioruptor UCD-200 sonicator (Diagenode, USA). To determine the concentration of sheared genomic DNA, an aliquot (50 μL) of chromatin solution with 200 μL of elution buffer (1% SDS and 0.1 mM NaHCO3) was used to reverse the cross-linking of the chromatin fractions, and the eluent was mixed with NaCl to a final concentration of 0.3 M, then incubated at 65 °C for 6 h. After reversal of the cross-links, the resulting sheared DNA was purified using the QIAquick PCR purification kit (Qiagen 28,106) and was quantified using NANO drop.
For the ChIP assays with chromatin enriched fraction, containing 5 μg of sheared DNA, immunoprecipitation was performed overnight with polyclonal anti-MYC (sc-789, Santa Cruz). Protein A agarose beads bearing immunoprecipitates were then subjected to sequential washes, and DNA was purified using the QIAquick PCR purification kit (Qiagen 28,106). Purified DNA was used to prepare libraries with the TruSeq ChIP Library preparation Kit (Illumina, Inc), and purified using the Agencourt AMPure XP system (Beckman Coulter, Inc). Sequencing was performed on an Illumina Hi-Seq 2500 with 101 bp single end reads (MACROGEN, Seoul). The ChIP-Seq data have been deposited in the NCBI Gene Expression Omnibus under GEO Series number GSE102920.
Protoplast isolation and transactivation assay
Rice protoplasts were isolated according to the method described for A. thaliana by Yoo et al. (2007) and the method described for rice by Chen et al. (2006) and Zhang et al. (2011) with some modifications [27, 28, 34]. Rice seeds were grown on MS medium at 28 °C for 10 days, then leaves and stems were cut into 0.5 mm pieces using razors and plasmolyzed in a 0.6 M mannitol solution for 20 min. The plasmolyzed tissue was transferred to enzyme solution [1.5% Cellulase R-10 (Yakult Pharmaceutical Ind. Co., Japan), 0.75% Macerozyme R-10 (Yakult Pharmaceutical Ind. Co., Japan), 0.6 M mannitol, 0.1% bovine serum albumin (BSA), 1 mM CaCl2, 5 mM β-mercaptoethanol and 10 mM MES (pH 5.7)]. After incubation at 25 °C for 4 h with gentle shaking in the dark, an equal volume of W5 solution [154 mM NaCl, 125 mM CaCl2, 5 mM KCl, 5 mM D-glucose and 2 mM MES (pH 5.7)] was added and the protoplasts were collected by centrifugation at 300×g for 8 min. Protoplast-enriched pellets were resuspended in MMG solution [0.6 M mannitol, 15 mM MgCl2 and 4 mM MES (pH 5.7)], and 15 μl DNA (including 3 μg effector (OsNAC TFs) and 1.5 μg reporter (promoter:fLuc and internal control 35S:rLuc) was used to transfect the protoplast solution harboring up to 2–3 × 106 cells. Dual luciferase activity was analyzed using a Dual-Luciferase Reporter Assay System (Promega, USA) and measured with an Infinite M200 System (Tecan, USA). Three independent transfections of each sample were performed, and the relative luciferase activity was calculated as the ratio between fLuc and rLuc.
Root-specific overexpression of four 6MYC-tagged OsNAC genes
To perform functional analysis of these OsNAC TFs, we generated four root-specific 6MYC-tagged OsNAC overexpressing plants (RCc3:6MYC-OsNAC5, 6, 9 and 10) with the coding sequences of OsNAC5 (Os11g0184900), OsNAC6 (Os01g0884300), OsNAC9 (Os03g0815100), and OsNAC10 (Os11g0126900) fused to the 6MYC tag . Five independent lines per construct were propagated to obtain T3 seeds. Proteins extracts from the roots of the RCc3:6MYC-OsNAC transgenic lines were subjected to Western blot analysis using an α-MYC antibody (Fig. 1). Transgenic lines with the highest levels of 6MYC-tagged OsNAC protein were chosen for further experiments.
Identification of genes differentially regulated in RCc3: 6MYC-OsNAC roots using RNA-Seq
To investigate OsNAC regulatory networks, RNA-Seq analyses were performed of roots from 20 d-old non-transformed (NT) and RCc3:6MYC-OsNAC expressing plants grown under normal growth conditions. Each data set was obtained from two biological replications. A total of 107 million raw reads were generated to identify the genes with different expression profiles, and subsequently trimmed to 102 million reads. Approximately 45 million reads were aligned to the reference rice genome (IRGSP v.1.0, http://rapdb.dna.affrc.go.jp) and 88.4%–91.6% of the reads were well matched (Additional file 2: Table S1). After read mapping, gene expression levels were calculated as FPKM (Fragments per kilobase of transcript per million mapped reads). Of a total of 41,635 transcripts, 8457 were excluded because of a zero FPKM value in one of the four RCc3:6MYC-OsNAC plants, leaving 33,178 transcripts for further analysis (Fig. 2a). From the 10 independent libraries from duplicate samples of RCc3:6MYC-OsNAC5, 6, 9, 10 expressing and NT control roots, a total of 4766 DEGs were identified using the following criteria |fc (fold change) | ≥ 1.5 and p value <0.05 (Fig. 2b). We identified 2147 and 2619 genes that were up- and down-regulated in one or more of RCc3:6MYC-OsNAC sample, respectively (Fig. 2c). Specifically, 414, 601, 748, and 1172 genes were up-regulated by more than 1.5-fold in RCc3:6MYC-OsNAC5, 6, 9 and 10 transgenic roots, respectively, and 685, 1000, 1086, and 1687 genes were down-regulated in RCc3:6MYC-OsNAC5, 6, 9 and 10 transgenic roots, respectively. Venn diagrams representing the distribution of single gene shared amongst the four RCc3:6MYC-OsNAC transgenic roots are shown in Fig. 2. Of these genes, 45 and 213 genes were up- and down-regulated, respectively, in all of the RCc3:6MYC-OsNAC transgenic roots (Fig. 2d and e, Additional file 3: Table S2). We validated the expression patterns of 23 genes that were up-regulated in all four (Fig. 3a), three (Fig. 3b), and one or two of the RCc3:6MYC-OsNAC expressing root samples (Fig. 3c), and that were down-regulated in all four RCc3:6MYC-OsNAC root samples (Fig. 3d), by qRT-PCR. These results correlated well with those of the RNA-Seq analysis (Fig. 3).
To gain insight into the possible functions of the target genes of the OsNAC TFs, Gene Ontology (GO) term enrichment analysis was performed using the web-based tool, AgriGO (Additional file 4: Table S3). The ‘membrane-bounded vesicle’, ‘apoplast’, ‘cell wall’, ‘external encapsulating structure’, and ‘extracellular region’ terms of the ‘cellular component’ category were significantly enriched for the genes up-regulated in the RCc3:6MYC-OsNAC roots (Additional file 4: Table S3 and Additional file 5: Table S4), suggesting an association between these subcellular areas/structures and the mechanisms of drought tolerance in the OsNAC overexpressing plants.
Identification of OsNAC5, 6, 9, and 10 binding loci by ChIP-Seq
The RNA-Seq analyses using the four RCc3:6MYC-OsNAC transgenic root samples indicated both direct and indirect targets of the OsNAC proteins. To pinpoint the direct targets and to further understand how the OsNAC TFs regulate gene expression related to drought tolerance, we performed a ChIP sequencing (ChIP-Seq) analysis using DNA immunoprecipitated from roots of the four RCc3:6MYC-OsNAC transgenic and NT control plants. A MYC-specific antibody was used together with OsNAC bound genomic DNA fragments and ChIP DNA libraries were constructed to determine OsNAC binding sites by deep sequencing. The sequencing depth of the immunoprecipitated DNA fragments was >1 million mapped reads. Background ChIP-Seq peaks from the RCc3:6MYC-OsNAC roots were filtered by comparisons with those of NT control roots. The total number of read bases were calculated and are listed in Additional file 6: Table S5. Data sets were normalized to RPM (read per million) values. A total of 44,550 genomic regions were identified and annotated by comparison with the reference rice genome (IRGSP a.1.0). Figure 4a shows the distribution of genomic loci enriched by individual OsNAC proteins, such as intergenic (beyond the 10 kb region upstream of TSS, transcription start site), promoter (2 kb upstream of TSS), 5’UTR (untranslated region), exon, intron, 3’UTR, and TTS (transcription terminate site) with the numbers of enriched peaks. Interestingly, we found similar distribution patterns of genomic loci enriched by each of the four of OsNAC proteins (Fig. 4a). Overall, 72% of the binding loci were within the genic regions, including the promoter, and 28% were in intergenic regions. (Fig. 4b). Most of the binding loci were shared by four, three, and two of the OsNAC TFs (Fig. 4c and d), 1199 of which were shared by all four, and were distributed among the intergenic and promoter regions (Fig. 4d).
Identification of direct OsNAC targets by cross-referencing ChIP-Seq and RNA-Seq data
We identified 18,083, 19,637, 20,467, and 16,973 enriched regions for OsNAC5, 6, 9, and 10, respectively, using ChIP-Seq (Figs. 4e and 5a). To identify the direct OsNAC target genes, we combined the ChIP-Seq and RNA-Seq data. Initially, we removed the intergenic regions from the ChIP enriched regions and obtained 13,134, 13,994, 14,557, and 12,489 enriched regions for OsNAC5, 6, 9, and 10, respectively (Fig. 5b). These were sorted for genomic loci enriched by >2.0-fold (log2 ≥ 1.0) compared with the NT control. We then combined the RNA-Seq data from individual RCc3:6MYC-OsNAC transgenic roots, considering genes whose expression was up-regulated by >1.5-fold (p value <0.05) compared with the NT roots, with the selected genomic loci enriched by individual OsNAC proteins (Fig. 5c). A total of 975 genomic loci, distributed among 460 promoters, 11 5’ UTRs, 211 exons, 116 introns, 148 TTSs, and 28 3’ UTRs, were enriched at least one of the OsNAC samples (Fig. 5c; Additional file 7: Table S6). We focused on 471 binding loci (460 promoters and 11 5’ UTRs) involving 391 genes after redundant binding loci were removed (Fig. 6).
To elucidate the functions of these 391 direct OsNAC target genes, we performed a Gene Ontology (GO) based functional enrichment analysis using AgriGO (Go analysis toolkit and database for agriculture community, http://bioinfo.cau.edu.cn/agriGO/index.php) (Additional file 8: Fig. S1). Singular enrichment analysis (SEA) was used to GO enrichment analysis by comparing the list of differentially expressed genes to all expressed genes. The ‘membrane and vesicle related cellular component including vesicle’ (GO:0031982), membrane-bound vesicle (GO: 0031988), ‘cytoplasmic vesicle’ (GO:0031410), and ‘cytoplasmic membrane-bounded vesicle’ (GO:001623) terms in the ‘cellular component’ category were significantly enriched. Additional file 9: Table S7 shows the 391 target genes affiliated with 8 GO terms in ‘cellular components’. More specifically, the candidate genes were enriched in the ‘membrane and vesicle related cellular component’ category, including 192 genes in the ‘organelle’ related GO terms, 87 genes in the ‘vesicle’ related GO terms, and 224 genes in the ‘cell related’ GO terms. Some of these were classified into the ‘transmembrane/transporter activity containing proteins’, ‘vesicle related proteins’, ‘auxin/hormone related proteins’, ‘carbohydrate metabolism related proteins’ and ‘transcription factors’ terms. Taken together, a total of 121 target genes were selected by combining the RNA-Seq and ChIP-Seq analyses as summarized in Table 1. For the individual OsNAC proteins, we found 72 genes that were up-regulated in the RNA-Seq analysis and whose promoter regions were enriched in the ChIP-Seq analysis. The data suggested that 12, 13, 21 and 26 genes are direct target genes of the individual OsNAC5, 6, 9, and 10 genes, respectively.
To assess whether the OsNAC proteins can bind to the promoter region and activate the transcription of the target gene, we performed a transient expression assay in rice protoplasts using a dual-luciferase reporter system. We have chosen one target gene, OsABI5 (Os01g0859300) that was found to be regulated by OsNAC5, 6 and 9, but not by OsNAC10. A reporter plasmid with the 2-kb promoter region of ABI5 (Os01g0859300) fused to a firefly luciferase gene was constructed (pABI5:fLuc, Fig. 7a). Effector plasmids with OsNAC5, 6, 9, and 10, and a CaMV35S:rLuc (renilla luciferase) expression cassette as an internal control, were constructed (Fig. 7a). These plasmids were transiently co-expressed in rice protoplasts in the combinations indicated in Fig. 7b. Expression of the reporter gene was strongly activated by OsNAC5, 6, and 9, but not by OsNAC10 (Fig. 7b), validating the results of the RNA-Seq and ChIP-Seq analyses.
Resolving TF networks, including the identification of direct target genes, is a common goal in functional genomic and molecular biology studies. Direct target genes have previously been identified using in vitro methods, such as the electrophoretic mobility shift assay (EMSA) and yeast one-hybrid analysis. However, the ChIP method has emerged as a powerful tool to detect interactions between DNA-binding proteins and genomic DNA. Our group has characterized the drought-inducible TFs, OsNAC5, 6, 9 and 10 [25, 29,30,31], overexpression of these genes specifically in transgenic rice roots (RCc3:OsNAC5, 6, 9 or 10) caused enlarged root diameter and drought tolerance. Putative target genes involved in root growth, development and abiotic stress response were identified in the roots of these plants by microarray [25, 30, 31] or RNA-Seq  analyses. However, it is difficult to identify a regulatory network or a specific pathway using such data due to a lack of information about binding interactions between TFs and target genes [35, 36].
In this study, we performed RNA-Seq and ChIP-Seq analysis of roots of RCc3:6MYC-OsNAC5, 6, 9 and 10 plants and identified direct OsNAC target genes by cross-referencing transcriptome (RNA-Seq) information and binding loci (ChIP-Seq) information. The RNA-Seq experiments led us to identify 45 up- and 213 down-regulated genes in all four RCc3:6MYC-OsNAC transgenic roots (Additional file 3: Table S2). These included genes involved in the ‘cellular components’, ‘vesicles’, ‘cytoplasmic membrane-bounded vesicle’, and ‘cytoplasmic vesicle’ GO terms (Additional file 4: Table S3). We previously demonstrated that the RCc3:OsNAC transgenic roots had increased root diameter due to an enlarged xylem and augmented cortical cells. These changes in root phenotype appeared to be correlated with the up-regulation of genes involved in signal transduction (abscisic acid and Ca2+), barrier formation (lignin and suberin biosynthesis), and cell development and morphogenesis [25, 30, 31]. However, RNA-Seq data alone provide no information about whether the target genes are regulated directly by the TF binding to cis-elements of the target gene promoter, or indirectly through regulation of other TFs or signaling molecules. To determine the direct targets of the four OsNAC TFs, a genome-wide analysis of the interactions of OsNAC5, 6, 9, and 10 with the promoters of regulated genes was performed. Cross-referencing of the ChIP-Seq and RNA-Seq data sets from the four RCc3:OsNAC root samples allowed us to identify 391 direct OsNAC5, 6, 9, and 10 target genes (Additional file 9: Table S7). Our ChIP-Seq data revealed that 158, 202, 198, and 239 loci were bound by OsNAC5, 6, 9, and 10, respectively, and that 65, 91, 115, 186 genes, respectively, were up-regulated (Fig. 6). Of the up-regulated genes in the RNA-Seq experiments, 8.7, 6.5, 7.4, and 4.0% were direct targets in the RCc3:6MYC-OsNAC5, 6, 9 and 10 roots, respectively. The direct target genes could be divided into 5 subgroups according to GO analysis: ‘transmembrane/transporter activity’, ‘vesicle related’, ‘hormone related’, ‘carbohydrate metabolic process’, and ‘TFs’ (Table 1). Of the genes encoding transmembrane proteins that were directly regulated by the OsNAC TFs, PIPs (plasma membrane intrinsic proteins) and TIPs (tonoplast intrinsic proteins) are known to be more abundantly expressed in roots than in leaves, and their expression responds to drought stress, contributing to the drought tolerance of the RCc3:6MYC-OsNAC5, 6, 9 and 10 plants. Moreover, their localization is associated with the presence of root hydrophobic barriers . This idea is supported by the observation that transgenic A. thaliana lines overexpressing the PIP1 gene from Panax ginseng , TIP1;2 from the highly drought-resistant Thellungiella salsuginea  and TIP2 from Solanum lycopersicum  exhibited higher survival rates, reduced water loss, and tolerance of water deprivation under drought conditions. The TIPs might help promote increases in turgor pressure and cell elongation, and the mediation of water transport across cell membranes, including the tonoplast, in response to drought [41, 42]. Together, these proteins are involved in controlling water movement during drought stress. NODULIN genes, which we found to be directly regulated by the OsNAC TFs, are from the MtN3/SALIVA/SWEET family (Os01g0881300, Os11g0508600, Os02g0513100), MtN21/EamA-like TRANSPORTER/UMAMIT family (Os05g0409500), and EARLY NODULIN-LIKE (ENODL) family (Os03g0115000 and Os06g0681200). Furthermore, in A. thaliana, the sugar transporter SAG29/SWEET15 (At5g13170) is a direct downstream target of an AtNAC TF, ANAC092/AtNAC2/ORE1 . We also identified an ENODL gene (Os06g0681200), which is associated with root architecture by trait ontology in the NBRP database (http://shigen.nig.ac.jp/rice/oryzabase/), as being involved in drought tolerance (Table 1). Among the direct targets were also four genes (Os03g0216700, Os03g0571900, Os10g0195000 and Os06g0494400) encoding a multi antimicrobial extrusion proteins and four genes (Os02g0510600, Os03g0751600, Os08g0403300, and Os10g0532300) encoding heavy metal transport/detoxification proteins (Table 1). It has previously reported that the metal stress signaling pathway is associated with controlling drought stress signaling , and that the genes involved encode transmembrane/transporters that control physiological processes essential for optimal plant growth, as well as for abiotic stress tolerance.
Other subgroups of direct OsNAC targets that have previously been shown to be involved in drought tolerance include the drought-inducible germin-like proteins, of which there are ten in rice . Germin-like proteins are divided into two groups based on oxalate oxidase and superoxide dismutase activities, which catalyzes the oxidative breakdown of oxalate to CO2 and H2O2. They are important in cell wall fortification, mediating auxin responsive physiological processes [41, 46], and germin proteins have also been suggested to contribute to the preparation of cell walls  by producing H2O2 required by peroxidase to catalyze the crosslinking of several cell wall components . Thus, germin and germin-like proteins play a role in strengthening plant cell walls to resist abiotic stress [34, 39, 40]. Zhong et al. (2010) showed that SND1 and VND7, two A. thaliana NAC-domain TFs, directly regulate the expression of not only downstream TFs, but also genes involved in the secondary wall biosynthesis and polysaccharide biosynthesis . Similarly, we found that the OsNAC TFs targeted in this study directly regulate genes involved in cell wall polysaccharide biosynthesis (Table 1). Two types of TFs, AP2/EREBP and bHLH, were among the direct targets (Table 1). Six AP2/EREBP- and six bHLH-domain genes were directly regulated by one of the four OsNAC TFs, such as OsERF1 (Os04g0429050) in RCc3:OsNAC5, OsERF54 (Os01g0657400) in RCc3:OsNAC6, OsERF57 (Os07g0227600) and OsERF107 (Os02g0521100) in RCc3:OsNAC9, and OsERF120 (Os06g0222400) and OsERF9 (Os03g0263000) in RCc3:OsNAC10. A similar number of the bHLH TF genes were regulated in the transgenic rice roots (Table 1). Thus, drought inducible OsNAC TFs bind to their target gene promoters, up-regulating different classes of TFs including the AP2/EREBP- and bHLH-domain genes, indicating a role in the crosstalk inherent in abiotic stress responses.
Previous studies have predicted the functional redundancy of OsNAC5, 6, 9, and 10 by observing similar root phenotypes such as an enlarged root and drought tolerance in the OsNAC overexpressors [25, 26,27,28]. In the current manuscript, we suggest the functional redundancy of OsNACs (Fig. 6) by the target gene expression profiles in the OsNAC overexpressors. We present a strategy to identify direct OsNAC target genes involving ChIP-Seq and RNA-Seq analyses of RCc3:OsNAC5, 6, 9, and 10 expressing roots, which are drought tolerant. These results were validated using transient expression assays in rice protoplasts using a dual-luciferase reporter system (Fig. 7). Many of the direct target genes have previously been shown to be involved in drought tolerance, highlighting the important roles of the four OsNAC TFs in the associated gene networks.
Abiotic stress triggers a wide range of plant responses, from the alteration of gene expression and cellular metabolism to change in plant growth and development. Thus, understanding the complex mechanism of drought tolerance is crucial for crop production. Uncovering the molecular networks of the TFs that are involved in drought tolerance provides information about the individual processes of how the TFs interact in a drought signaling pathway. We obtained 391 direct targets of OsNACs by comparing genes from the results of RNA-Seq with those of ChIP-Seq. Among the genes selected from RNA-Seq results, the ratio of becoming a direct target of OsNACs ranged from 4.0% to 8.7% depending on RCc3:OsNACs plants. Interestingly, OsNACs target several cellular components, including transmembrane/transporter, vesicle, auxin/hormone, carbohydrate metabolic processes and transcription factors. Many of the target genes that are up-regulated by OsNACs act as the cellular components which would alter the root architectures of RCc3:OsNACs for drought tolerance. This would be a valuable resource for the functional dissection of the molecular mechanisms in plant drought response.
chromatin immunoprecipitation sequencing
quantitative reverse transcribed polymerase chain reaction
Ge Y, Bai G, Stoerger V, Schnable JC. Temporal dynamics of maize plant growth, water use, and leaf water content using automated high throughput RGB and hyperspectral imaging. Comput Electron Agric. 2016;127:625–32.
Bartels D, Sunkar R. Drought and salt tolerance in plants. Crit Rev Plant Sci. 2005;24:23–58.
Wang H, Wang H, Shao H, Tang X. Recent advances in utilizing transcription factors to improve plant abiotic stress tolerance by transgenic technology. Front Plant Sci. 2016;7:67.
Agarwal PK, Agarwal P, Reddy MK, Sopory SK. Role of DREB transcription factors in abiotic and biotic stress tolerance in plants. Plant Cell Rep. 2006;25(12):1263–74.
Zhu JK. Abiotic stress signaling and responses in plants. Cell. 2016;167(2):313–24.
Singh D, Laxmi A. Transcriptional regulation of drought response: a tortuous network of transcription factors. Front Plant Sci. 2015;6:895.
Nakashima K, Yamaguchi-Shinozaki K, Shinozaki K. The transcriptional regulatory network in the drought response and its crosstalk in abiotic stress responses including drought, cold, and heat. Front Plant Sci. 2014;5:170.
Franco-Zorrilla JM, López-Vidriero I, Carrasco JL, Godoy M, Vera P, Solano R. DNA-binding specificities of plant transcription factors and their potential to define target genes. Proc Natl Acad Sci U S A. 2014;111(6):2367–72.
Nuruzzaman M, Manimekalai R, Sharoni AM, Satoh K, Kondoh H, Ooka H, Kikuchi S. Genome-wide analysis of NAC transcription factor family in rice. Gene. 2010;465(1–2):30–44.
Jensen MK, Kjaersgaard T, Nielsen MM, Galberg P, Petersen K, O'Shea C, Skriver K. The Arabidopsis Thaliana NAC transcription factor family: structure-function relationships and determinants of ANAC019 stress signaling. Biochem J. 2010;426(2):183–96.
Nikovics K, Blein T, Peaucelle A, Ishida T, Morin H, Aida M, Laufs P. The balance between the MIR164A and CUC2 genes controls leaf margin serration in Arabidopsis. Plant Cell. 2006;18(11):2929–45.
Duval M, Hsieh TF, Kim SY, Thomas TL. Molecular characterization of AtNAM: a member of the Arabidopsis NAC domain superfamily. Plant Mol Biol. 2002;50(2):237–48.
He XJ, RL M, Cao WH, Zhang ZG, Zhang JS, Chen SY. AtNAC2, a transcription factor downstream of ethylene and auxin signaling pathways, is involved in salt stress response and lateral root development. Plant J. 2005;44:903–16.
Podzimska-Sroka D, O’Shea C, Gregersen PL, Skriver KNAC. Transcription factors in senescence: from molecular structure to function in crop. Plants. 2015;4(3):412–48.
Shao H, Wang H, Tang XNAC. Transcription factors in plant multiple abiotic stress responses: progress and prospects. Front Plant Sci. 2015;6:902.
Chung PJ, Jung H, Jeong DH, Ha SH, Choi YD, Kim JK. Transcriptome profiling of drought responsive noncoding RNAs and their target genes in rice. BMC Genomics. 2016;17:563.
Eulgem T, Rushton PJ, Robatzek S, Somssich IE. The WRKY superfamily of plant transcription factors. Trends Plant Sci. 2000;5(5):199–206.
Rijpkema AS, Gerats T, Vandenbussche M. Evolutionary complexity of MADS complexes. Curr Opin Plant Biolo. 2007;10(1):32–8.
Zong W, Tang N, Yang J, Peng L, Ma S, Xu Y, Li G, Xiong L. Feedback regulation of ABA signaling and biosynthesis by a bZIP transcription factor targets drought-resistance-related genes. Plant Physiol. 2016;171(4):2810–25.
Puranik S, Sahy PP, Srivastava PS, Prasad MNAC. Proteins: regulation and role in stress tolerance. Trends Plant Sci. 2012;17(6):369–81.
Ohashi-Ito K, Oda Y, Fukuda H. Arabidopsis VASCULARRELATEDNAC-DOMAIN6 directly regulates the genes that govern programmed cell death and secondary wall formation during xylem differentiation. Plant Cell. 2010;22(10):3461–73.
Nakashima K, Takasaki H, Mizoi J, Shinozaki K, Yamaguchi-Shinozaki KNAC. Transcription factors in plant abiotic stress responses. Biochim Biophys Acta. 2012;1819(2):97–103.
Kuo MH, Allis CD. In vivo cross-linking and immunoprecipitation for studying dynamic protein:DNA associations in a chromatin environment. Methods. 1999;19(3):425–33.
Angelini C, Costa V. Understanding gene regulatory mechanisms by integrating ChIP-seq and RNA-seq data: statistical solutions to biological problems. Front Cell Dev Biol. 2014;2:51.
Jeong JS, Kim YS, Redillas MCFR, Jang G, Jung H, Bang SW, Choi YD, Ha SH, Reuzeau C, Kim JK. OsNAC5 overexpression enlarges root diameter in rice plants leading to enhanced drought tolerance and increased grain yield in the field. Plant Biotechnol J. 2013;11(1):101–14.
Bowler C, Benvenuto G, Laflamme P, Molino D, Probst AV, Tariq M, Paszkowski J. Chromatin techniques for plant cells. Plant J. 2004;39(5):776–89.
Yoo SD, Cho YH, Sheen J. Arabidopsis mesophyll protoplasts: a versatile cell system for transient gene expression analysis. Nature Protoc. 2007;2(7):1565–72.
Chen S, Tao L, Zeng L, Vega-Sanchez ME, Umemura K, Wang GLA. Highly efficient transient protoplasts system for analyzing defence gene expression and protein-protein interactions in rice. Mol Plant Pathol. 2006;7(5):417–27.
Lee DK, Chung PJ, Jeong JS, Jang G, Bang SW, Jung H, Kim YS, Ha SH, Choi YD, Kim JK. The rice OsNAC6 transcription factor orchestrates multiple molecular mechanisms involving root structural adaptions and nicotianamine biosynthesis for drought tolerance. Plant Biotechnol J. 2016;15(6):754–64.
Redillas MCFR, Jeong JS, Kim YS, Jung H, Bang SW, Choi YD, Ha SH, Reuzeau C, Kim JK. The overexpression of OsNAC9 alters the root architecture of rice plants enhancing drought resistance and grain yield under field conditions. Plant Biotechnol J. 2012;10(7):792–805.
Jeong JS, Kim YS, Baek KH, Jung H, Ha SH, Choi YD, Kim M, Reuzeau C, Kim JK. Root-specific expression of OsNAC10 improves drought tolerance and grain yield in Rice under field drought conditions. Plant Physiol. 2010;153(1):185–97.
Dubin MJ, Bowler C, Benvenuto GA. Modified gateway cloning strategy for overexpression tagged proteins in plants. Plant Methods. 2008;4:3.
Jang IC, Nahm BH, Kim JK. Subcellular targeting of green fluorescentprotein to plastids in transgenic rice plants provides a high-level expression system. Mol Breed. 1999;5:453–61.
Zhang Y, Su J, Duan S, Ao Y, Dai J, Liu J, Wang P, Li Y, Liu B, Feng D, Wang J, Wang H. 2011. A highly efficient rice green tissue protoplast system for transient gene expression and studying light/chloroplast-related processes. Plant methods. 2011;7(1): 30.
Haynes BC, Maier EJ, Kramer MH, Wang PI, Brown H, Brent MR. Mapping functional transcription factor networks from gene expression data. Genome Res. 2013;23(8):1319–28.
Kogenaru S, Qing Y, Guo Y, Wang N. RNA-seq and microarray complement each other in transcriptome profiling. BMC Genomics. 2012;13:629.
Jarzyniak KM, Jasiński M. Membrane transporters and drought resistance - a complex issue. Front Plant Sci. 2014;5:687.
Peng Y, Lin W, Cai W, Arora R. Overexpression of a Panax Ginseng tonoplast aquaporin alters salt tolerance, drought tolerance and cold acclimation ability in transgenic Arabidopsis plants. Planta. 2007;226(3):729–40.
Sade N, Vinocur BJ, Diber A, Shatil A, Ronen G, Nissan H, Wallach R, Karchi H, Moshelion M. Improving plant stress tolerance and yield production: is the tonoplast aquaporin SlTIP2; 2 a key to isohydric to anisohydric conversion? New Phytol. 2009;181(3):651–61.
Wang LL, Chen AP, Zhong NQ, Liu N, XM W, Wang F, Yang CL, Romero MF, Xia GX. The Thellungiella Salsuginea tonoplast aquaporin TsTIP1;2 functions in protection against multiple abiotic stresses. Plant Cell Physiol. 2014;55(1):148–61.
Rietz S, Bernsdorff FE, Cai D. Members of the germin-like protein family in Brassica Napus are candidates for the initiation of an oxidative burst that impedes pathogenesis of Sclerotinia sclerotiorum. J Exp Bot. 2012;63(15):5507–19.
Wudick MM, Luu DT, Maurel CA. Look inside: localization patterns and functions of intracellular plant aquaporins. New Phytol. 2009;184(2):289–302.
Matallana-Ramirez LP, Rauf M, Farage-barhom S, Dortay H, Xue G. 2013. NAC transcription factor ore1 and senescence-constitute a regulatory cascade in Arabidopsis. Mol Plant. 2013;6(5):1438–52.
Oono Y, Yazawa T, Kawahara Y, Kanamori H, Kobayashi F, Sasaki H, Mori S, Wu J, Handa H, Itoh T, Matsumoto T. Genome-wide transcriptome analysis reveals that cadmium stress signaling controls the expression of genes in drought stress signal pathways in rice. PLoS One. 2014;9(5):e96946.
Ke Y, Han G, He H, Li J. Differential regulation of proteins and phosphoproteins in rice under drought stress. Biochem Biophys Res Commun. 2009;379(1):133–8.
Banerjee J, Maiti MK. Functional role of rice germin-like protein1 in regulation of plant height and disease resistance. Biochem Biophys Res Commun. 2010;394(1):178–83.
Zhong R, Lee C, Ye ZH. Global analysis of direct targets of secondary wall NAC master switches in Arabidopsis. Mol Plant. 2010;3(6):1087–103.
Xu Y, Buchholz WG, DeRose RT, Hall TC. Characterization of a rice gene family encoding root-specific proteins. Plant Mol Biol. 1995;27(2):237–48.
We thanks MACROGEN Inc., Seoul, Korea (http://macrogen.com) for the help with the ChIP-Seq, RNA-Seq and bioinformatics analyses.
This work was supported by the Basic Science Research Program (NRF-2013R1A6A3A04060627 to PJC) and by the Mid-career Researcher Program (NRF-2017R1A2B4007457 to J-KK) through the National Research Foundation of Korea, Ministry of Education. The funding body was not involved in the design of the study, analysis or interpretation of data, or writing the manuscript.
Availability of data and materials
The data supporting the conclusions of this article are included within the article and its additional files.
Ethics approval and consent to participate
The transgenic seeds, generated from Oryza sativa L. japonica cv. Illmi, used in this study were grown and maintained in the greenhouse facility of Seoul National University, with permission from Ministry of Science, ICT and Future Planning, Republic of Korea (LML14–491).
Consent for publication
The authors declare that they have no competing interests.
Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.
List of primers used in this study. (XLSX 12 kb)
Information about RNA-Seq data obtained by Illumina Hi-Seq 2500 sequencing. 1Raw seq. Reads no.: Total number of reads used in sequence analysis. 2Raw seq. Reads length: Total length of reads used in sequence analysis. 3GC (%): Percentage of GC found in sequence reads. 4Q20 (%): Percentage of bases with quality score of 20 or above. A base with a quality score of 20 has a probability of being an incorrect call of 1 in 100. (XLSX 14 kb)
List of up- and down-regulated differentially expressed genes (DEGs) (|fold change| ≥ 1.5, and p value <0.05,) in RCc3:OsNAC5, 6, 9, and 10 rice plants. (XLSX 54 kb)
Gene ontology analysis of genes that are up-regulated in RCc3: OsNAC5, 6, 9, and 10 transgenic roots. Genes were associated with gene ontology (GO) terms, including Biological Process, Molecular Function, and Cellular components. FDR, False Discovery Rate, p-value <0.05 using AgriGO and REVIGO. Main enriched processes are presented for each of RCc3:OsNAC5, 6, 9, and 10. (XLSX 15 kb)
List of genes associated with gene ontologies, including Biological Process, Molecular Function, and Cellular components. (XLSX 152 kb)
ChIP-Seq data obtained by Illumina Hi-Seq 2500 sequencing. (XLSX 16 kb)
Direct target genes of OsNAC5, 6, 9, and 10 identified by merging of chromatin immunoprecipitation (ChIP)-Seq and RNA-Seq data sets. Each color showed that Log2 ratio ≥ 1.0 in ChIP, and f.c. (fold change) ≥ 1.5, p-value <0.05 in RNA-seq, and red, sky blue, light green, and purple mean OsNAC5, OsNAC6, OsNAC9, and OsNAC10, respectively. 975 sheet means a total of 975 genomic loci, including 460 promoters, 11 5’ UTR, 211 exons, 116 introns, 148 TTs, and 28 3’ UTRs, and 460 sheet means 460 promoters that were enriched in at least one of the OsNACs. Distance to TSS, the distance between the nearest TSS and the center of the peak/region. (XLSX 1228 kb)
Singular enrichment analysis was performed using AgriGO (http://bioinfo.cau.edu.cn/agriGO/) to identify enriched gene ontologies up-regulated in OsNAC transgenic rice plants. a. OsNAC5, b. OsNAC6, c. OsNAC9, d. OsNAC10. (TIFF 4245 kb)
Integration of RNA-Seq and chromatin immunoprecipitation (ChIP)-Seq data sets. (XLSX 36 kb)
About this article
Cite this article
Chung, P.J., Jung, H., Choi, Y.D. et al. Genome-wide analyses of direct target genes of four rice NAC-domain transcription factors involved in drought tolerance. BMC Genomics 19, 40 (2018). https://doi.org/10.1186/s12864-017-4367-1