Genome-wide identification, evolution, and expression analysis of the NPR1-like gene family in pears

The NONEXPRESSOR OF PATHOGENESIS-RELATED GENES 1 (NPR1) plays a master regulatory role in the salicylic acid (SA) signal transduction pathway and plant systemic acquired resistance (SAR). Members of the NPR1-like gene family have been reported to the associated with biotic/abiotic stress in many plants, however the genome-wide characterization of NPR1-like genes has not been carried out in Chinese pear (Pyrus bretschneideri Reld). In this study, a systematic analysis was conducted on the characteristics of the NPR1-like genes in P. bretschneideri Reld at the whole-genome level. A total nine NPR1-like genes were detected which eight genes were located on six chromosomes and one gene was mapped to scaffold. Based on the phylogenetic analysis, the nine PbrNPR1-like proteins were divided into three clades (Clades I–III) had similar gene structure, domain and conserved motifs. We sorted the cis-acting elements into three clades, including plant growth and development, stress responses, and hormone responses in the promoter regions of PbrNPR1-like genes. The result of qPCR analysis showed that expression diversity of PbrNPR1-like genes in various tissues. All the genes were up-regulated after SA treatment in leaves except for Pbrgene8896. PbrNPR1-like genes showed circadian rhythm and significantly different expression levels after inoculation with Alternaria alternata. These findings provide a solid insight for understanding the functions and evolution of PbrNPR1-like genes in Chinese pear.


INTRODUCTION
Pears are threatened with various diseases which were caused by fungus, bacterial, viruses, nematodes and insect bites in environment (Robert-Seilaniantz, Grant & Jones, 2011). These diseases limit pear quality and generate severe economic losses. Changes in phytohormones concentration or sensitivity can be triggered under biotic and abiotic stress conditions (Pieterse et al., 2009). Salicylic acid (SA) signaling triggers the resistance against biotrophic and hemibiotrophic pathogens (Nie et al., 2017), whereas a combination of Jasmonates (JA) and Ethylene (ET) signaling activates the resistance against necrotrophic pathogens (Attaran & He, 2012). These two pathways are mostly antagonistic: elevated biotroph resistance is often correlated with increased necrotroph Focusing on the function of redox rhythm and NPR1 in the plant immunity and circadian clock  proved the complex relationship between plant immune response and circadian clock. Plant circadian clocks play an important role in regulating the growth-defense balance. The nocturnal stomatal closure and active defense in the morning are good examples of how the circadian clock stop the unsustainable energy shift to immunity. Through the interaction with the circadian redox rhythm of metabolism, the circadian clock shuts off immune induction to prevent conflict with growth at night (Dodd et al., 2005;Somers et al., 1998). The circadian clock combines environmental cues with temporal information to regulate the plant physiology and development.
NPR1-like genes play an important role in resistance disease, growth, and development of plant tissues and organs in many plants, for examples, tomato, strawberry, avocado, banana (Backer et al., 2015;Endah et al., 2008;Li et al., 2019;Wang et al., 2018). However, NPR1-like gene family was still unidentified and uncharacterized in Chinese pear. In this study, it was the first report about discovering and identifying the NPR1-like gene family in Chinese pear genome. Also, PbrNPR1-like genes were analyzed the phylogeny, gene structure, conserved motif, cis-elements, chromosomal location, and tissue-specific expression of P. bretschneideri Reld. Using qPCR technology, we indicated how NPR1-like genes response after SA treatment and A. alternata infection and expression changes within 24 h. This study provided a valuable foundation for further functional analysis of the PbrNPR1-like genes and pear genetic improvement.

Plant materials and treatment
The three-year-old P. bretschneideri Rehd (cv.Yali) trees were cultivated in the greenhouse (16/8 h light/dark, 24 C/20 C day/night and relative humidity of 25%) in the experimental orchard of Hebei Agricultural University, Baoding, Hebei, China. Yali leaves were collected at 0, 1, 3, 6, 12, 24, 48, and 72 h after 0.2 mM SA treatment. Samples from different tissue of 7-year-old Yali were gained at different stages of development. Flower buds were sampled 1 day before flowering, and then flower, young leaves, young stems, young fruit, immature leaves, old stems, mature fruits and seed were collected at 1, 15,15,29,115,115,159 and 159 days after flowering (DAF) in 2018, respectively. The samples were taken at 8 time points throughout the day, 9:00 am, 12:00 am, 15:00 pm, 18:00 pm, 21:00 pm, 24:00 am, 3:00 am and 6:00 am. The A. alternata was isolated from the diseased leaves, and grew at 28 C incubator.

Identification of PbrNPR1-like gene family member
The genome of P. bretschneideri Rehd (cv. Dangshan Suli) was downloaded from the Genome Database for Rosaceae (GDR) (http://www.rosaceae.org/) to identify the genome-wide NPR1-like gene family in Chinese pear. Using hmmsearch tool in the HMMER (v3.0) software package, we detected the PbrNPR1-like genes in the Protein database. HMMER file according to the ankyrin domain (PF00023) and BTB/POZ domain (PF00651) was downloaded from the Pfam protein database (http://pfam.xfam.org/). The CDS length, protein size, isoelectric point (IP), molecular weight (MW) of the deduced PbrNPR1-like genes were analyzed using the ExPASy website (https://web.expasy. org/compute_pi/).

Sequence alignment and phylogenetic analysis
The NPR1-like proteins sequences of rice, arabidopsis, grape, apple, pear, and other plants were obtained from NCBI and Ensamble plants (http://plants.ensembl.org/index. html) to perform alignments by using Clustal W 2.0 (Larkin et al., 2007). The phylogenetic tree was constructed through the neighbor-joining (NJ) method with the bootstrap value 1,000 replicates in the MEGA 6.0 software (Tamura et al., 2013).

Sequence analysis of PbrNPR1-like gene family
The individual PbrNPR1-like genes structure was displayed using Gene Structure Display Server 2.0 software (GSDS, http://gsds.cbi.pku.edu.cn/) the exon/intron structure based on the alignments of NPR1-like genes CDS and their genomic sequences in Chinese pear and Arabidopsis (Hu et al., 2015). The MEME software program (http://meme-suite.org/) was used to analyze the NPR1-like protein motif with the following parameters: (1) the width of optimum motif 6-200; (2) the maximum number of motif 10 (Bailey et al., 2009). The domain was analyzed by using a web CD-search tool (https://www.ncbi.nlm.nih.gov/ Structure/bwrpsb/bwrpsb.cgi) (Lu et al., 2020).

Cis-elements in the promoters of PbrNPR1-like genes
We analyzed the promoter region up to 2,000 bp except the Pbrgene12425 gene including 1,300 bp. The promoter sequences of the PbrNPR1-like genes were analyzed using PlantCARE databases (http://bioinformatics.psb.ugent.be/webtools/plantcare/html/).

Real-time quantitative PCR (qPCR)
The total RNA of samples was extracted using the plant RNA extraction kit (Takara, Dalian, China). First-stand cDNA was synthesized through the FastQuant RT Kit with gDNase (TianGen, Beijing, China). The qPCR was performed using the Top Green qPCR SuperMix Kit (TranGen, Beijing, China). Each 20.0 ml reaction mixture included 10.0 ml SYBR supermix, 2.0 ml cDNA template, 0.4 ml forward and reverse primer (10.0 mM) and 7.2 ml ddH 2 O. LightCyclerÒ 96 system (Roche, Germany) was used for qPCR using the following PCR parameters: pre-denaturation 94 C for 30 s, followed by 42 cycles of 94 C denaturation for 5 s, 55 C annealing for 15 s, 72 C extension for 10 s, melting-curve analysis was done by 95 C for 15 s, and then a constant increase from 58 C to 95 C at a 2% ramp rate. All gene-specific primers for PbrNPR1-like gene family were designed by using DNAMAN 8.0 software except Pbrgene33340, and listed in Table S2 except Pbrgene33340. The relative expression of the PbrNPR1-like genes was calculated using the 2 −ΔΔCT method (Livak & Schmittgen, 2001), and three technical replicates and three biological replicates were applied.

Statistical analysis
The results of analysis the variance (ANOVA) and duncan multiple comparison were analyzed using SPSS software (SPSS version 17.0, Chicago, IL, USA). Data was presented in graphs through OriginPro 9.1.0 software (Microcal Software, Inc., Northampton, MA, USA).

Identification of PbrNPR1-like genes
After removing the redundant entries, nine PbrNPR1-like genes with typical BTB/POZ and ankryin domain were identified to unevenly distribute on the chromosomes. Among the nine PbrNPR1-like genes, three PbrNPR1-like genes were detected on chromosome 5, the other 6 genes were identified on chromosome 2, 9, 10, 14, 17, and a scaffold respectively. The characteristics of PbrNPR1-like genes were listed in Table 1, including the chromosomal location, intron number, coding sequence (CDS) length, amino acid length, protein molecular weight (MW) and isoelectric point (PI). The lengths of the open reading frame (ORF) of PbrNPR1-like genes range from 1,494 bp (Pbrgene40077/34018) to 1,770 bp (Pbrgene8341). The length of the PbrNPR1-like proteins range from 498 amino acids (Pbrgene40077/34018) to 591 amino acids (Pbrgene8341). The molecular weight of these proteins range from 54.52 kDa (Pbrgene40077/34018) to 66.25 kDa (Pbrgene8341) and the PIs range from 5.7 (Pbrgene8896/33340) to 6.4 (Pbrgene8895). All genes contain one or three introns.

Phylogenetic and protein domain analysis of PbrNPR1-like gene family
The phylogenetic analysis of family determines the origin of the evolutionary process of the family, and the internal branch of the evolutionary tree reflects the distance of evolutionary relationships among different genes. To explore the evolutionary relationship of PbrNPR1-like gene family, NPR1-like protein sequences from 20 different monocotyledon and dicotyledon plant species were used to construct the phylogenetic tree (Table S1). According to the phylogenetic tree, these NPR1-like proteins were divided into three clades: clade I, II, III. As shown in Fig. 1, Pbrgene12425 was classified into the cluster of AtNPR1 and AtNPR2 as the positive regulator of SAR. However, the Pbrgene6286, Pbrgene8895, Pbrgene8896, Pbrgene2529, Pbrgene8341 and Pbrgene33340 involved in the clade II of AtNPR3 and AtNPR4 negative SAR regulation (Zhang et al., 2006). The other two protein Pbrgene34018 and Pbrgene40077 were classified into the clade III along with AtBOP1 and AtBOP2 to participated in the organ determinacy and symmetry (Hepworth et al., 2005;McKim et al., 2008).
Gene structure, motif composition, and domain of NPR1-like gene family The classification of PbrNPR1-like genes was presented in Fig. 2B, in which clade I and II both contained four exons and three introns, and the Clade III (Pbrgene34018 and Pbrgene40077) included two exons and one intron. The protein domain composition revealed that all NPR1-like proteins included a N-terminal BTB-POZ domain and a ankyrin repeat domain located central region (Fig. 2C). The genes in Clade I and the Clade II both have NPR1-like-C terminal region which has been proved to be an important part in NPR1-like gene family (Boyle et al., 2009;Rochon et al., 2006). The NPR1-like-C terminal protein usually existed in the eukaryotes with the size between 251 and 588 amino acids. The C-terminal region included a nuclear localization signal (NLS), a NIMINTERACTING (NIMIN) 1/2 protein binding site, and a conservative penta-amino acid motif (LENRV) . The presumed domain DUF3420 was found in eukaryotes with the length about 50 amino acids which was functionally uncharacterized. The Pbrgene2529 did not have the DUF3420 domain. Ten conserved motifs were identified in the PbrNPR1-like gene family (Fig. S1, Table S3). As shown in Fig. 2D, Clade I and Clade II harbored motif 1-10 and Clade III featured motifs 1, 2, 3, 5, 6, 7, 9, and 10. The same conserved motif ingredients in each clade supported the polygenetic classification in the PbrNPR1-like gene family.
Multiple sequence alignments were executed to examine the conservation of residues, domains, and motifs in PbrNPR1-like and known-function AtNPR1-like (AtNPR1 to AtNPR6) (Fig. 2E). We found that npr1-2 (Cys150Tyr), nim1-2 (His300Tyr), and npr1-1 (His334Tyr) mutations in AtNPR1 were completely conserved in all PbrNPR1-like proteins. Nim1-4 (Arg432Lys) in AtNPR1 and npr4-4D (Arg419Gln) in AtNPR4 mutant sites were conserved in Clade I and Clade II. C82, C216, and C156 cysteine residues in AtNPR1 participated in its oligomer-monomer transition were highly conservation in PbrNPR1-like proteins. The Arg432 residue in AtNPR1, Arg428 in AtNPR3, and Arg419 in AtNPR4 required for their expression of SA were conserved in Clade I and Clade II in PbrNPR1-like protein.

Cis-elements in the promoter regions of PbrNPR1-like genes
In this study, about 2,000 bp promoter region sequence was identified in the PbrNPR1-like gene family except for the Pbrgene12425 with 1,300 bp. The cis-elements were divided into three aspects through the PlantCARE database, including plant growth and development, stress responses, and hormone responses (Table S4). In the first aspect, MRE, G-box, GTI motif, GATA motif, I-box, AE-box, ATCT-motif, TCT motif, TCCC motif and Box 4 for light responsiveness, circadian for the circadian control, O 2 -site for zein metabolism regulation (Fig. 3) (Abdullah et al., 2018). In the second aspect, included a range of stress-related elements, such as, STRE and TC-rich repeats involved in stress . WUN-motif related to a wound-responsive element (Ni, Cui & Gelvin, 1996), LTR for low-temperature responsiveness (Xiang, Huang & Xiong, 2007). In the third aspect, we detected the ABRE related to ABA was the most ordinary motif (Kim et al., 2011), CGTAC motif and TGACG motif for MeJA responsiveness, GARE-motif, P-box, and TATC box for gibberellin responsiveness, TCA-element for SA responsiveness, the TGA-element for auxin-response (Abdullah et al., 2018). These results indicated that PbrNPR1-like gene family have the advantages for enhancing abiotic stress responses and hormones responsiveness and may respond to abiotic stress and hormones.
Expression patterns for the PbrNPR1-like genes in different plant tissues To obtain a first glance at the functions of PbrNPR1-like genes during different developmental stages of pear, the transcript accumulation levels were observed in 10 tissues, including the leaves (April, July), stems (April, July), flowers (4 April), flowers buds (8 April), fruits (May, September), seed (September) and bud (8 March) during the budding and reproductive stages. The result revealed that PbrNPR1-like genes were constitutively expressed in different tissue, and the expression level varied (Fig. 4). The expressions level of PbrNPR1-like genes in young fruits, mature leaves and mature stems were higher than that in mature fruits, young leaves and young stems except for Pbrgene6268/8896/12425. Among the eight genes, the expression level of Pbrgene6286 gene was relatively stable in 10 tissues/ organs, Pbrgene12425 gene was moderately expressed in leaves, Pbrgene8895 was highly expressed in leaf, stem, seed and bud tissue, Pbrgene8341/8896 were relatively highly  Differential expression pattern of the PbrNPR1-like genes in response to SA treatment The transcript levels of Pbrgene2529/6286/8341/8895/12425/34018/40077 were up-regulated among that the Pbrgene2529/8341/8895/12425 reached their highest expression levels at 12 h with 2.29-fold, 2.72-fold, 6.86-fold, and 3.96-fold comparing to the 0 h samples, respectively, and then slowly returned to the baseline levels after 12 h (Fig. 5). Comparing to the 0 h samples, the expression levels of Pbrgene6286, Pbrgene34018, Pbrgene40077 were the highest at 6 h with 4.45-fold, 3 h with 2.57-fold, 3 h with 1.97-fold, respectively. Whereas Pbrgene8896 was first down-regulated and then up-regulated at 24 h. We observed the expression of the Pbrgene6286/8341/12425 gene changed more than other PbrNPR1-like genes. In brief, the expression patterns of PbrNPR1-like genes were significantly affected by exogenous SA hormones.
Expression analysis of PbrNPR1-like genes in response to A. alternata these genes were reached a small peak at 12 hpi (Fig. 6). The expression of Pbrgene12425 was significantly higher than the other PbrNPR1-like genes at 96 hpi. The Pbrgene2529/ 8341/8896 were up-regulated and reached a peak at 12 hpi, and dramatically decreased except Pbrgene8896 was then up-regulated at 72 hpi.     Figure 6 Expression of PbrNPR1-like genes in Yali leaves after inoculation with A. alternata. The leaves were harvest at 0, 6, 12, 24, 48, 72, 96, 120 h after A. alternata treatment. Different letters associated with each time point indicate statistically significant differences at the 5% level. The same letters indicate that the statistics did not differ significantly at P < 0.05 according to Duncan's multiple range tests.
Full-size  DOI: 10.7717/peerj.12617/ fig-6 Differential expression profile of the PbrNPR1-like genes in the circadian rhythm The circadian rhythm expression patterns of PbrNPR1-like genes under long daytime were compared. The results showed that the expression of all the PbrNPR1-like genes had the phenomenon of circadian rhythm. Most of the eight genes had similar expressed patterns except for Pbrgene40077, which were expressed higher in the daytime than in the night (Fig. 7). The expression of Pbrgene40077 was the highest at 21:00 pm, and low in all other time. This circadian clock phenomenon was related to the function of the family members, the clade I, II involved in plant immune, and the clade III took part in the organ symmetry and determinacy during the leaf morphogenesis.

DISCUSSION
Nine PbrNPR1-like genes, Pbrgene12425/8895/8896/8341/6286/34018/2529/33340/40077 were identified in Chinese pear. PbrNPR1-like proteins have similar gene structures, domains, and conserved motifs and amino acid residues with NPR1-like sequences in Arabidopsis, suggesting that their orthologs probably display similar biological functions in Chinese pear, while the difference between the groups specified their function was diversified. Similar results had existed in Oryza sativa, Malus domestica, Populus trichocarpa (Shao et al., 2013;Yuan et al., 2007;Zhang et al., 2016). These findings indicated that the NPR1-like genes were highly conserved in numerous plant species. The result of phylogenetic tree suggested that potential functional distinctions between the PbrNPR1-like gene family were existed. A large number of stress response elements were existed in the promoter regions, such as ARE, LTR, MBS, MYB, MYC, STRE, TC-rich repeats, which suggested that PbrNPR1-like genes played an important role in the stress response process. In addition, G-box, I-box, and other optical response elements were found in the promoter regions, indicating that these genes may also be regulated by optical signals. Meanwhile, we detected the ABRE, CGTAC, TGACG, GARE, TATC box, P-box, TCA-element, and TGA-element in the nine genes, so we speculate these genes take part in the hormone-responsive and are induced by GA, SA, ABA, auxin, and MeJA hormone (Zhao et al., 2020). All results implied that different PbrNPR1-like genes played a role in special environment at different times, and thus the patterns of PbrNPR1-like genes responses to hormones were complex.
To further state the possible functions of the NPR1-like gene family in the development and growth of P. bretschneideri Rehd, the transcription profiles of PbrNPR1-like genes were studied through qPCR in 10 different tissues. The study displayed Pbrgene12425 was expressed in all the tissues with the highest expression in mature leaves. Pbrgene6268 and Pbrgene8896 were expressed more higher in flower bud than that in flower. The tissue-and stage-specific expression profile implied a specific function of the two genes in the early flower development (Shi et al., 2013). Pbrgene34018 and Pbrgene40077 in clade III displayed specific expression in the leaves, stem, and flower (Fig. 4), implying that they may participated in the organ symmetry and determinacy during the leaf morphogenesis which was similar to the study on Persea Americana. The expression of PaNPR1, PaNPR2, and PaNPR4 were seen in all sorts of tissues (Backer et al., 2015). The expression of PaNPR3 was expressed much higher in fruit and aerial tissues. The PaNPR5 was displayed in the roots. Therefore, database retrieval and functional prediction of PbrNPR1-like genes in different tissues and various stages of development indicated that PbrNPR1-like genes may play an important role in plant growth, and some PbrNPR1-like genes may have unique functions at specific stages of development.
The expression profile of the PbrNPR1-like genes showed obviously different after SA induction. The same results were obtained in Qinguan apple resistant disease cultivars, the MdNPR1/2/3/4/5/6/7 were up-regulated at 6 h, and MdNPR8 were expressed to the highest level at 12 h (Zhang et al., 2016). The time of reaching peak expression is different in the various plants, the possible reason is different defense system and the time of the infection and reaction varies. The expression of defense-related genes increased after 3 h in the morning with SA inducing in the morning or evening under constant light for 3 h while the expression of growth-related genes increased in the evening (Zhou et al., 2015). The NPR1-like genes would be induced after SA treatment and participated in the SA signal transduction which plays an important role in local defense and distal tissue of wide-spectrum SAR (Baldwin & Meldau, 2013).
Plant disease activated by the pathogens promote defense immunity to plant pathogens and the defense immunity is a difficult mechanism including the triggering of multiple immune mechanisms reaction . Several resistant genes and proteins associated with the pathogenesis have been isolated and they can be used to improve the plant defense against different diseases. The PR proteins play an important role in plant defense systems. NPR1 can interact with some members of the TGA family of bZIP trascription factors which combine with the as-1-like (TGACG) element in the PR gene promoter and take charge of PR gene expression (Zhang et al., 2003;Zhou et al., 2000). The results (Fig. 6) not only showed that there were significantly different in the expression of the family genes after infecting with A. alternata and some PbrNPR1-like genes were assciated with pear resistance to A. alternata in pear. PbrNPR1-like genes were pathogen-inducible and participated in the immune system of pear, enhancing the disease resistance.
Studies have found that the expression of NPR1 has circadian rhythm in A. thaliana (Zhou et al., 2015), which is balanced between the regulation of self-growth and immunity. The expression of PbrNPR1-like gene family showed circadian rhythm in natural conditions during a day, indicating that the PbrNPR1-like genes may be conserved in plants circadian rhythm. The stomata on the surface of plant leaves opens in the morning, making it vulnerable to pathogen invasion, the plant activates a defense mechanism which makes it more resistant to disease in the morning than at night, the plant grows normally at night (Korneli, Danisman & Staiger, 2014;Zhang et al., 2013). This result verified the conclusion that SA mediated immune responses in the morning helped to avoid the conflicts between SA immune responses and growth-related activities that need to transport water at night (Zhou et al., 2015). This change in plant defense response may reflect an adaptation to change physical conditions during the day and the temperatures and humidity are generally more conducive to pathogen challenges in the morning (Karapetyan & Dong, 2018).
The functions of the NPR1-like genes maybe much more intricate in Chinese pear than in A. thaliana. Overexpression of AtNPR1 can improve the disease resistance of A. thaliana, while AtNPR3 and AtNPR4 have redundant functions and opposite functions to the AtNPR1 (Ding et al., 2018). When the biotic and abiotic infection, SA will combine with its receptors (NPR1, NPR3, and NPR4) and induce the expression of PR protein, which leads to trigger the immune responses. AtNPR2 was more similar to the AtNPR1 than the other NPR1 orthologs and would play an important role in the SA perception and acted as an evolutive reservoir of the NPR1 (Castello et al., 2018). In pear, it was speculated there were more homologs would interact with each other to regulate the balance PbrNPR1. The PbrNPR1 homologs can be regulated by different transcripts or transcription factors.
This study provided the preliminary PbrNPR1-like gene family information and functional annotation of the nine discovered PbrNPR1-like genes from pear. Sequence structure, homology, and phylogenetic analysis suggested that seven PbrNPR1-like proteins might participate in the defense responses, the rest two genes were likely involved in tissue development. Hormone and expression in various tissues provided a support for this and allow future research to learn much more about the possible role of the PbrNPR1-like genes in SAR in pear. The future efforts will be focused on the localization and the intracellular interactions of defense-related PbrNPR1-like proteins as well as