Preliminary analysis of two NAC transcription factor expression patterns in Larix olgensis

Qing Cao·Peiqi An·Sufang Zhang·Junhui Wang·

Hanguo Zhang1·Lei Zhang1

Abstract The NAC transcription factor family is plantspecific with various biological functions. However, there are few studies on the NAC gene involving coniferous species. Bioinformatics research and expression analysis of NAC genes in Larix olgensis can be used to analyse the function of the NAC gene in the future. Screening of excellent genetic materials and molecular breeding have been utilized to cultivate high-quality, stress-resistant larches. According to the transcriptome data for L. olgensis, the genes Unigene81490 and Unigene70699 with complete ORFs (open reading frames) were obtained by conserved domain analysis and named LoNAC1 and LoNAC2, respectively. The cDNAs of LoNAC1 and LoNAC2 were 1971 bp and 1095 bp in length, encoding 656 and 364 amino acids, respectively. The molecular weights of the proteins encoded by the two genes were predicted to be 72.61 kDa and 41.13 kDa, and subcellular localization analysis indicated that the proteins were concentrated in the nucleus. The results of real-time quantitative PCR analysis showed that at different growth stages and in different tissues of L. olgensis, the relative expression levels of the two NAC genes were highest in the stem, and the expression differences were more obvious in non-lignified tissues. After drought, salt and alkali stress and hormone treatment, expression was induced to different degrees. The expression levels of LoNAC1 and LoNAC2 in semi-lignified L. olgensis were higher than in the other two periods (non-lignified and lignified), and expression levels significantly increased under drought and salt stress. Relative expression levels changed under hormone treatment. It is speculated that these two genes may not only be related to drought and salt stress and secondary growth but may also be induced by hormones such as abscisic acid. Overall, LoNAC1 and LoNAC2 are genetic materials that can be used for molecular breeding of larch.

Keywords Larix olgensis·NAC transcription factor·Bioinformatics analysis·Expression patterns

Introduction

The NAC transcription factor family is plant-specific and widely distributed in terrestrial plants. It is also regarded as one of the gene families with the most transcription factors (Riechmann et al. 2000; Olsen et al. 2005). More than 8000 transcription factors of the NAC family have been found in plants, and 151 inArabidopsis thaliana(L.) Heynh. alone (de Oliveira et al. 2011). This study found that the N-terminus of the NAC protein has a highly conserved domain, while the C-terminus has a transcription activation domain and presents diversity, which is an important recognition-related feature of the NAC protein structure (Aida et al. 1997). There are approximately 150 amino acids in the N-terminal domain of the NAC protein, which can be divided into five subdomains (Chen et al. 2019). The C-terminus has a simple amino acid sequence with high repetition and contains a greater number of Ser, Thr, and Glu and some acidic amino acid residues than the N-terminal domain (Olsen et al. 2005). By aligning the NAC protein sequences ofArabidopsis, researchers have shown that some common sequences can be found even at the C-terminus (Ooka et al. 2003).

In the NAC gene family, thePetuniaNAM gene was the first one found to be related to plant morphogenesis (Souer et al. 1996). Researchers found that inArabidopsis, the key factor regulating secondary wall thickening of xylem fibre cells was the specific expression of the NAC familyNST1andNST3/SND1genes (Mitsuda et al. 2007). Huang et al. (2015) found that the genesNAC29andNAC31in rice are related to the regulation of cellulose synthesis. He et al. (2005) reported that theArabidopsistranscription factorAtNAC2was highly expressed in roots under high salt conditions, the lateral roots of plants overexpressing this gene were well-developed, and changes in auxin and ethylene were observed.Arabidopsis ANAC096can help plants survive dehydration-related osmotic stress and is related to ABA-induced genes (Xu et al. 2013). Liu et al. (2018) noted thatTsNAC1can target an important proton transporter to improve salt tolerance. Researchers have studied the function of the NAC gene family but most of these have focused on plants such asArabidopsis, tobacco, and rice. The function of NAC family of genes in coniferous species remains to be explored and verified.

Larix olgensisA. Henry, belonging to the larch genus of the pine family, is a fast-growing timber species as well as a species for soil and water conservation in China (Zhang 2012). With the advance of science and technology resulting in the genetic improvement of larch, molecular breeding has been combined with traditional breeding to improve the efficiency of genetic improvement and to accelerate the improvement process (Levee et al. 1997). In this study, two full-length NAC genes ofL. olgensiswere used, and their functions initially estimated from the expression levels of genes in different tissue parts and under drought and salinealkali conditions and hormone treatment. This will provide the basis for the verification of the NAC gene in subsequent experiments and for the screening of genetic material by genetic engineering-based breeding.

Materials and methods

Plant materials

Seedlings ofL. olgensisat different growth stages, non-lignified (approximately 60 days), semi-lignified (approximately 120 days) and lignified (180 days), were used as experimental material. Roots, stems and needles were separated, wrapped in aluminium foil, labelled, rapidly frozen in liquid nitrogen and stored at - 80 °C. The 180-dayL. olgensismaterial was treated with solutions of PEG6000 (25% w/v), NaCl (250 m mol/L), and NaHCO3(50 m mol/L). Samples were taken at 0 h, 12 h, 24 h, 48 h and 96 h (three samples for each treatment and time). Gibberellin A3 (GA3), abscisic acid (ABA), N-(phenylmethyl)-9H-purin-6-amine(6-BA), indole-3-acetic acid(IAA), 2,4-dichlorophenoxyacetic acid (2,4-D) and methyl jasmonate (MeJA), each at 50 mg/L, were used to treat samples (at approximately 90 days) for 2 h, 4 h, 8 h, 12 h, 24 h, 48 h, 72 h, and 96 h (Lv et al. 2014; Zhao 2019). UntreatedL. olgensisseedlings were used as controls.

Prediction of gene sequence structure and function

Through the National Center for Biotechnology Information (NCBI) online tool Blastx, more than 10 candidate sequences of NAC genes obtained in the laboratory were compared, and the structural domains of the sequences were predicted and analysed with the online tool CD-search (Marchler-Bauer et al. 2015).Unigene81490(LoNAC1) andUnigene70699(LoNAC2) were selected with the NAC family conserved structural domain and complete open reading frame. MEGA5.0 software (Tamura et al. 2011) was used to construct the neighbour-joining tree. Amino acid sequences corresponding to full-length genes in Blastx and to several full-length genes with the highest similarity in the evolutionary tree were combined, and these were used to perform multiple sequence comparisons by BioEdit software. Protparam (https:// web. expasy. org/ cgi- bin/ protp aram/ protp aram) was used to predict and analyse the physicochemical properties of the protein encoded by the full-length NAC gene. GOR4 (https:// npsa- prabi. ibcp. fr/ cgi- bin/ npsa_ autom at. pl? page= npsa_ gor4. html) was used to predict the secondary structure, and WOLF PSORT (https:// wolfp sort. hgc. jp/) and SwissModel (https:// www. swiss model. expasy. org/) to predict the subcellular localization and three-dimensional structure of the protein, respectively.

Analysis ofLoNAC1 and LoNAC2 expression inL. olgensis

Total RNA was extracted using Universal Plant Total RNA Extraction Kit (BIOTEKE, Beijing, China), and cDNA obtained by reverse transcription of total RNA was used as a template (ReverseScriptRTreagent Kit, TaKaRa). Primer Premier 5.0 software designed the qRT-PCR primers (Table 1), the primers specifically screened by gel electrophoresis detection.LoTublinwas selected as a reference gene and amplification was performed automatically accordingto the Real Master Mix (SYBR Green) kit instructions and an ABI 7500 real-time PCR instrument was used. Three replicates were set up in the quantitative PCR instrument and the results analysed by the 2-ΔΔCtmethod (Relative Expression = 2-((Ctgene-CtLoTublin)sample-(Ctgene-CtLoTublin)con-trol)= 2-(ΔCtsample-ΔCtcontrol)= 2-ΔΔCt) (Pfaffl 2001).

Table 1 List of primer sequences used in qRT-PCR

Results

Identification ofLoNAC1 andLoNAC2

The NCBI online tool blastx was used to compare the genes from theL. olgensistranscriptome database obtained in the laboratory. Two genes,Unigene81490andUnigene70699, with the NAC family special conserved domain (NAM), were obtained and namedLoNAC1 andLoNAC2, respectively. The cDNAs were 1974 bp and 1098 bp in length, encoding 657 and 365 amino acids, respectively (Fig. 1).

Fig. 1 Nucleotide and encoded amino acid sequences of LoNAC1 and LoNAC2

Predicted physicochemical properties ofLoNAC1 andLoNAC2 proteins

According to the prediction of physical and chemical properties, the theoretical molecular weights of theLoNAC1andLoNAC2proteins are 72.61 kDa and 41.13 kDa (1 kDa = 1000 Da = 1000 g/mol), respectively, and the predicted isoelectric points are 4.99 and 4.79, respectively. TheLoNAC1andLoNAC2proteins contain 63 and 40 positively charged amino acids, respectively, and 91 and 61 negatively charged amino acids, respectively. The instability coefficient of theLoNAC1gene is 38.13, the instability coefficient of theLoNAC2gene is 50.12, and the protein hydrophobicity of the two sequences is - 0.460 and - 0.667, respectively. Subcellular localization analysis predicted that the proteins were concentrated in the nucleus.

Using the GOR4 webpage to predict the secondary structure of the two genes (Fig. 2) showed that they also have certain similarities in secondary structure.LoNACland LoNAC2 proteins are mainly composed of random coils.LoNAC1andLoNAC2contain 21.31% and 26.58% α-helices, 18.57% and 16.44% extension chains, and 60.12% and 56.99% random coils, respectively, and lack β-turns. SwissModel homology modelling was used to predict the tertiary structure of the two proteins. As shown in Table 2 and Fig. 3, the N-terminus ofLoNAC1andLoNAC2proteins have a highly conserved domain (NAM), which is an important recognition-related feature of the NAC protein structure and the prediction results for the tertiary structure of the two proteins were very similar.

Fig. 2 Secondary structure of the a LoNAC1 and b LoNAC2 proteins; vertical blue lines represent a spiral; vertical redmiddle line a fold; vertical purple line a curl

Table 2 The result of CD-Search of LoNAC1 and LoNAC2

Fig. 3 LoNAC1 (a) and LoNAC2 (b) protein tertiary structure prediction

Fig. 4 Phylogenetic tree of the LoNAC1 and LoNAC2 proteins

Sequence alignment and evolution tree analysis ofLoNAC1 and LoNAC2 proteins

LoNAC1andLoNAC2were translated to obtain the amino acid sequence and the sequence was compared with the sequences of NAC family proteins ofArabidopsis thalianaby MEGA5.0 software to construct a phylogenetic tree (Fig. 4). Figure 4 shows thatLoNAC1andLoNAC2ofL. olgensisandAtNAC053,AtNAC078,AtNAC082, andAtNAC103 ofArabidopsis thalianaare clustered on the same branch. It is speculated that these genes are relatively similar in evolutionary kinship. Using BioEdit software,LoNAC1andLoNAC2were multiplexed and compared with theArabidopsisNAC membersAtNAC053,AtNAC078,AtNAC082, andAtNAC103 (Fig. 5). There is a highly conserved domain at the N-terminus of the amino acid, and this domain can be divided into multiple subdomains. The C-terminus shows diversity, a highly variable transcription activation region, which is consistent with the characteristics of NAC family transcription factors.

Analysis of the tissue expression patterns ofLoNAC1 andLoNAC2

Different seedling growth stages, i.e., non-lignified (approximately 60 days), semi-lignified (approximately 120 days) and lignified (180 days) (Fig. 6) stored at - 80 °C was used to extract plant RNA and reverse transcribe it to the corresponding cDNA. The gene expression in different tissues (roots, stems, needles) was measured at different growth stages by qRT-PCR (Fig. 7). In the needles ofL. olgensis, both genes had the highest relative expression levels during the semi-lignification period, and the expression levels at different stages were in the order: semi-lignified ＞ lignified ＞ non-lignified; among the stems, the relative expression of genes was also highest during the semi-lignification period, and the relative expression levels at different stages were in the order: semi-lignified ＞ non-lignified ＞ lignified. In the roots, the relative expression level was again highest in the semi-lignification stage, and the expression levels at different stages were in the order: semi-lignified roots ＞ lignified roots ＞ non-lignified roots, similar to the needles. In the three tissues, the two genes had the highest relative expression levels during the semi-lignification period. The relative expression levels in the roots were significantly different among the three growth stages (Fig. 7).

Fig. 5 Multiple sequence alignment analysis of LoNAC1 and LoNAC2 with AtNAC053, AtNAC078, AtNAC082 and AtNAC103; “*” represents a conserved amino acid; “:” a conservative replacement; “.” a non-conservative replacement

Fig. 6 Different growth stages of L. olgensis

Expression analysis ofLoNAC1 andLoNAC2 under abiotic stress

Under drought stress, the expression levels ofLoNAC1andLoNAC2were upregulated at five time periods (Fig. 8). The relative gene expression levels inL. olgensisseedlings for 96 h reached the highest levels, increasing 12.7-fold (LoNAC1) and 12.8-fold (LoNAC2), respectively. Under salt stress,LoNAC1was down regulated at 12 h and 24 h, whileLoNAC2was down regulated at 12 h. In addition, both expression levels were up regulated and the relative expression level reached the highest level at 96 h, approximately 10 times that before treatment. Under alkali stress, the relative expression levels ofLoNAC1at 12 h and 24 h andLoNAC2at 24 h and 48 h were lower than in the controls, and expression was inhibited. Expression was induced at the remaining time periods, and the relative expression level was highest 96 h after treatment, reaching 2-3 times the value of the controls.

Fig. 7 Expression analysis of LoNAC1 and LoNAC2 at different growth stages of L. olgensis

Fig. 8 Expression analysis of LoNAC1 and LoNAC2 under different stress levels: a PEG, b NaCl, c NaHCO3

Analysis of the expression patterns ofLoNAC1 andLoNAC2 under hormone induction

RT-PCR was used to measure gene expression inL. olgensistreated with six hormones (Fig. 9). At the same time, according to the results, the relative expression levels ofLoNAC1andLoNAC2changed under different hormone induction conditions. Under IAA treatment, bothLoNAC1andLoNAC2genes were up regulated at most time periods;LoNAC1 expression reached a maximum of approximately eightfold at 96 h, andLoNAC2 of approximately threefold at 96 h. Under MeJA treatment,LoNAC1 was only slightly up regulated at 96 h andLoNAC2 slightly up regulated at 4 h and 72 h. In addition, the expression of both genes was down regulated. Under 2,4-D treatment,LoNAC1andLoNAC2were both up regulated;LoNAC1expression reached a maximum of 7.5-fold at 72 h, andLoNAC2 5.7-fold at 12 h. Under ABA treatment, both genes were up regulated, except for a slight down regulation at 4 h and 48 h, and their relative expression levels reached a maximum at 96 h, increasing 7.0- and 5.4-fold, respectively. Under GA3treatment, both genes were up regulated, except for a slight down regulation at 24 h, and their relative expression levels reached a maximum at 2 h. Under 6-BA treatment,LoNAC1was down regulated at 2 h, 4 h, 8 h, 12 h and 24 h, but was most highly down regulated at 8 h.LoNAC2was up regulated, except for the down regulation observed at 8 h, and the relative expression reached a maximum at 48 h.

Fig. 9 Expression analysis of LoNAC1 and LoNAC2 under different hormone treatments: a IAA, b MeJA, c 2,4-D, d ABA, e GA3, f 6-BA

Discussion

Bioinformatics analysis ofLoNAC1 andLoNAC2

Bioinformatics analysis showed that the full-length sequences of both theLoNAC1 andLoNAC2 genes contain special NAM domains. The coefficient of instability ofLoNAC1 is less than 40, and that ofLoNAC2 greater than 40, i.e., theLoNAC1 protein is stable and theLoNAC2 protein unstable. It is speculated that theLoNAC1 protein may be present for a long duration inL. olgensis, while theLoNAC2 protein may appear at certain stages. Both proteins had negative hydrophobicity results and are presumed to be hydrophilic proteins. The evolutionary tree and homology analysis of the amino acids encoded byArabidopsisNAC family genes showed thatLoNAC1 andLoNAC2 are clustered on the same branch withAtNAC053,AtNAC078,AtNAC082, andAtNAC103, and it is speculated that these genes are close in evolutionary kinship. Their structures and functions may be similar.Under the condition of protein toxicity,AtNAC053 andAtNAC078 work together to activate the expression of many factors so that the plant produces sufficient protein homeostasis factors, such as the 26S proteasome to regulate the protein toxicity stress response. TheAtNAC053andAtNAC078proteins play a role as central regulators in this process (Gladman et al. 2016). However, the induction ofAtNAC103 expression inArabidopsisdepends onbZIP60, which participates in the growth and development ofArabidopsisand in the endoplasmic reticulum stress response (Sun et al. 2018). It is believed that the functions ofLoNAC1 andLoNAC2 are similar and in the phylogenetic tree analysis, the affinity ofLoNAC1 andLoNAC2 reached 97%. These two genes may cooperate in some regulatory mechanisms. Due to the differences in gene functions between different plants and the influence of distant and close relationships, this conclusion requires further study.

Analysis ofLoNAC1 andLoNAC2 qRT-PCR results

Abiotic stress restricts the growth and development of plants, leading to a reduction in the yield and quality of agricultural and forestry crops. Research indicates that environmental stresses have caused almost half crop losses globally. Among environmental stress factors, drought is considered the most important factor restricting the development of global agriculture (Boyer 1982). As soil salinity becomes an increasingly serious problem, saline-alkali stress has also emerged as one of the factors that restrict crop growth. Plant hormones are a class of substances that regulate growth and development and are related to processes of plant environmental adaptation. They not only independently but also cooperatively regulate seed maturation, dormancy, and germination, vegetative and reproductive growth, and plant adaptation to abiotic and biological stresses during growth. Previous studies have found that NAC genes participate in both the stress response of plants(Lu et al. 2007; Mao et al. 2014), secondary growth (Wang et al. 2015) and also regulate plant cells and tissue death effects (Tran et al. 2009; Ma et al. 2018), which are related to hormone synthesis and regulatory networks (Fujita et al. 2004; Gao et al. 2010; Mao et al. 2017). NAC is a family of transcription factors with various biological functions discovered in recent research.

According to the analysis of the qRT-PCR results, during the non-lignification period, the distribution ofLoNAC1andLoNAC2in different tissues varied greatly. There was little difference in the lignification period, and the two genes had the highest relative expression levels in the semi-lignification period, indicating that the genes may participate in the secondary growth ofL. olgensis. Under drought and salt stresses, the relative expression levels ofLoNAC1andLoNAC2were quite different, while under alkaline stress, the differences were relatively small. These results indicate that the genesLoNAC1 andLoNAC2 respond to drought and salt stress. Under hormone treatment, bothLoNAC1andLoNAC2were induced to different degrees, and the relative expression levels in the 2,4-D, ABA and GA3treatments were significantly different, suggesting that the expression of the two genes was highly correlated with these three hormones.

2,4-D is a representative, artificially synthesized plant hormone that is an auxin analogue and widely used as a growth regulator for some crops (Hu et al. 2019). Studies have shown that 2,4-D can delay senescence in citrus plants. The levels of many endogenous hormones changed, and defence-related genes and proteins were up regulated, which improved the ability of defence under adversity. Some NAC family genes were over expressed (Ma et al. 2014). ABA plays a vital role in stress responses and regulates various developmental processes such as seed maturation and dormancy, organ shedding, and leaf senescence (Erik and Stokstad 2010). There also exist ABA-dependent regulatory pathways that respond to various abiotic stresses such as drought, high salt stress, and cold stress (Yamaguchi-Shinozaki and Shinozaki 2005).

This study found that some NAC genes are related to the corresponding pathways of 2,4-D and ABA stress. For example, in 2,4-D-treated citrus plants, some NAC family genes were up regulated which improved the ability of defence under adversity (Hu et al. 2019).OsNAC52may respond to ABA to increase drought tolerance in transgenic plants (Gao et al. 2010). In this study, analysis of the qRT-PCR results show that the expression of two NAC genes was significantly up regulated under drought and salt stress, and these genes were induced by treatment with the hormones 2,4-D and ABA. It may also be related regulatory networks inL. olgensis. The promoters of the two genes are currently unknown but will be identified by RACE technology in the future, and promoter sequence analysis will be completed to establish the regulatory network of these two genes.

Conclusion

LoNAC1andLoNAC2have special conserved NAM domains and it is believed that their functions are similar to those of NAC family genes. Analysis of the qRT-PCR results show that the two genes were induced by 2,4-D, ABA, and GA3and participated in the secondary growth process ofL. olgensisand in the response to drought and salt stress. These genes can be used as excellent genetic material for molecular breeding.

Author ContributionsQC and LZ conceived and designed the study. QC, PA and SZ performed the experiments. QC wrote the paper. JW, HZ and LZ reviewed and edited the manuscript. All authors read and approved the manuscript.

Open AccessThis article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article’s Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article’s Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http:// creat iveco mmons. org/ licen ses/ by/4. 0/.