Identification of key genes and transcription factors in ageing Arabidopsis papilla cells by transcriptome analysis

Identification of key genes and transcription factors in ageing Arabidopsis papilla cells by transcriptome analysis

Plant Physiology and Biochemistry 147 (2020) 1–9 Contents lists available at ScienceDirect Plant Physiology and Biochemistry journal homepage: www.e...

2MB Sizes 0 Downloads 34 Views

Plant Physiology and Biochemistry 147 (2020) 1–9

Contents lists available at ScienceDirect

Plant Physiology and Biochemistry journal homepage: www.elsevier.com/locate/plaphy

Research article

Identification of key genes and transcription factors in ageing Arabidopsis papilla cells by transcriptome analysis

T

Hong Yee,1, Fei Rena,1, Haoyu Guoc, Liping Guod, Jianfang Baib,∗, Yukun Wange,∗∗ a

School of Agricultural Science and Engineering, Shaoguan University, 288 Daxue Road, Zhenjiang District, Shaoguan, 512000, PR China Beijing Engineering Research Center for Hybrid Wheat, Beijing Academy of Agriculture and Forestry Sciences, Beijing, 100097, PR China c College of Life Science, Capital Normal University, Beijing, 100048, PR China d School of Landscape Architecture, Beijing Forestry University, Beijing, 100083, PR China e Division of Biological Science, Nara Institute of Science and Technology, Ikoma, 6300192, Japan b

A R T I C LE I N FO

A B S T R A C T

Keywords: Programmed cell death Papilla cell Floral organ Transcriptome DEGs

Programmed cell death (PCD) play essential roles in plant growth and development. Stigmatic papilla cells form an indispensable organ for plant reproduction. The lifetime of papilla cells is tightly controlled, and the developmental PCD (dPCD) process is involved in papilla cell death. Hence, papilla cell death is a good model for studying on PCD process. In this study, the dPCD signal was visualized in dying papilla cells by detecting the GUS signal of the PCD-related reporter gene BIFUNCTIONAL NUCLEASE 1 (BFN1). We found that the GUS was not expressed at young stage, but strongly expressed in papilla cells at the ageing stage, indicating the PCD process was triggered to terminate the papilla cell fate. Given this, the RNA-Seq data set, which covered the information of the whole lifespan of papilla cells, was analyzed aiming to understand which genes and pathways were involved in papilla cell death. 37 differential expressed genes (DEGs) were isolated. Moreover, the pathways related to energy production and transportation, autophagy, and plant hormone signal transduction were considered as the key pathways involved in the papilla cell death. 9 types, total of 104 transcriptional factors (TFs) were identified as well. Finally, a putative working model of papilla cell death was integrated. The findings herein will enrich the knowledge of the dPCD-mediated pathway in regulating plant organ/tissue growth, development, senescence, and death. Our study will provide some referential gene resources for studying on the dPCD in other plant organs or tissues.

1. Introduction Stigma, an essential floral organ, plays an important role in pollination in plants. The morphological features of stigma are huge different, some stigmas show single-cell forms, while some display multicellular forms. Meanwhile, the papilla cells, which existed on the top of stigma or the surface of stigma, have unique structures and functions on pollination process. Floral longevity is a species-specific process, and it shows a closely relationship with the reproductive strategy in plants (Shibuya et al., 2016). During the process of sexual reproduction in flowering plants, pollen cells will germinate after conglutination and hydration. After the success pollination, the floral organs will be suffered post-pollination induced senescence, especially the papilla cells on stigma will die rapidly. However, the senesce of floral organs can also be triggered in the absence of pollination. Such kind of species-

specific lifespan of floral organs can range from few hours to several weeks (Williams, 1965). The senescence of most floral organs now is considered to be related to the programmed cell death (PCD) process. The lifespan of these floral organs is tightly controlled by plant itself (Van Doorn and Woltering, 2008). The PCD process can be considered as the actively controlled cessation of a cell's vital functions (Dickman et al., 2017; Gao et al., 2018), indicating the termination of the senescence process (Tomas, 2013; Gao et al., 2018). PCD in plants has been classified into two main types, one is known as environmental PCD (ePCD), which is a response of various abiotic and biotic stresses; the other one is defined as developmental PCD (dPCD) that shows the self-control ability of plants (Daneva et al., 2016). It is a very important topic to uncover the natural manners that plants control their lifespans. Therefore, the dPCD process, as an indispensable regulation pathway of plant growth and



Corresponding author. Corresponding author. E-mail addresses: [email protected] (J. Bai), [email protected] (Y. Wang). 1 Contributed equally to this work. ∗∗

https://doi.org/10.1016/j.plaphy.2019.12.008 Received 14 August 2019; Received in revised form 4 December 2019; Accepted 6 December 2019 Available online 09 December 2019 0981-9428/ © 2019 Elsevier Masson SAS. All rights reserved.

Plant Physiology and Biochemistry 147 (2020) 1–9

H. Ye, et al.

approach for controlling the false discovery rate (FDR). The expression level of each DEG should show an expression Fold Change≥2, and show a FDR < 0.01. Gene function was annotated based on the following databases: Nr (NCBI non-redundant protein sequences), Nt (NCBI non-redundant nucleotide sequences), Pfam (Protein family), KOG/COG (Clusters of Orthologous Groups of proteins), and SWISSPROT (Tombuloglu et al., 2013, 2015). K-means clustering method was carried out to do gene cluster analysis.

successful sexual reproduction (Huysmans et al., 2017), should be understood deeply. During the past years, a large set of dPCD-associated genes have been isolated. Among these genes, BIFUNCTIONAL NUCLEASE1 (BFN1), PUTATIVE ASPARTIC PROTEASE A3 (PASPA3), RIBONUCLEASE3 (RNS3), CYSTEIN ENDOPEPDITASE 1 (CEP1), DUF 679 MEMBRANE PROTEIN4 (DMP4), and EXITUS1 (EXI1) are key components. The expression levels of these genes are up-regulated during the dPCD process (Olvera-Carrillo et al., 2015). Besides, the plant hormone signaling pathways are also considered as important factors in the dPCD regulation (Ueda and Kusaba, 2015). Moreover, the transcriptional factors, such as NAC, MYB, bHLH and so on, are also involved in dPCD process in plants. Recently, Gao and colleagues used the model plant Arabidopsis to uncover the senescent mechanism of papilla cells without pollination (Gao et al., 2018). Their results revealed that two NAC transcriptional factors (TFs), KIRA1 (ANAC074) and ORE1 (ANAC092), showed the sufficient effects on regulating the senescence of papilla cells, and also demonstrated that the dPCD process happened during papilla cell death (Gao et al., 2018). The Arabidopsis stigma has a tightly regulated functional lifespan, and it will be terminated through the dPCD process (Gao et al., 2018). In order to gain more information during the dPCD process in aging papilla cells in Arabidopsis stigma, in this study, we performed a comparative transcriptome analysis. We mainly focus on the pathways, including the genes involved in, and TFs. The results may provide us some new insights into the regulatory mechanism of dPCD, and may give us an estimable reference about the importance of dPCD in plant growth regulation.

2.4. Gene ontology enrichment and KEGG pathway analysis Gene Ontology (GO) enrichment analysis of the DEGs was implemented by the GOseq R packages based on Wallenius non-central hyper-geometric distribution, which can adjust for gene length bias in DEGs. Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway database was used to enrich DEGs in different pathways (http://www. genome.jp/kegg/). We used KO-Based Annotation System software to test the statistical enrichment of differential expression genes in KEGG pathways. 2.5. Differentially expressed genes validation Total RNA was extracted from Arabidopsis papilla cells using an RNeasy Plant Mini Kit (Qiagen). Genomic DNA was digested using an RNase-free DNase Set (Qiagen) before cDNA synthesis. Reverse transcriptional reaction was performed using the PrimeScript RT Master Mix (Takara) according to manufacturer's protocol. The quantitative reverse transcriptional PCR (qRT-PCR) experiments were carried out on an Eco Real-Time PCR System, the program followed as the description previously (Wang et al., 2019). The 2−ΔΔCT method was used to calculate gene expression. Each experiment was replicated three times, and EIF4 (AT3G13290) was used as the internal reference.

2. Materials and methods 2.1. Transcriptome data source The RNA-sequencing data sets were downloaded from the public database ArrayExpress (https://www.ebi.ac.uk/arrayexpress), and the accession number was E-MTAB-6279 (Gao et al., 2018). In these data sets, wild type Arabidopsis (Col-0) flowers were emasculated at flower stage 12c (the petals rapidly reach the height of the medial stamens). Total RNA was extracted from papilla cells in three independent biological replicates at stages 1, 2 and 3, representing the young, mature and senescent stage, respectively (Gao et al., 2018).

2.6. Vector construction and transgene To construct proBFN1::GUS-GFP, the 2.0 Kb BFN1 promoter was amplified using genomic DNA from Col and sub-cloned into pENTR/DTOPO vector according to the manufacturer's protocol (Thermo Fisher Scientific). Then the resulting BFN1 promoter fragment were Gatewaycloned into pBGWFS7 vector by the LR reaction. The recombinant construct was transformed into Agrobacterium tumefaciens strain GV3101 by using the freeze–thaw method. The transgenic Arabidopsis plants were produced by using the Agrobacterium-mediated floral dip method. T1 seeds were collected and screened using antibiotic Basta. More than 30 T1 plants were obtained, and representative lines were chosen for further characterization.

2.2. Mapping reads to the reference genome FASTX-Toolkit and Trimmomatic programs were used to remove reads containing adapter, reads containing ploy-N, and low quality reads from raw data. In this step, clean data(clean reads) were obtained by removing reads containing adapter, reads containing ploy-N and low quality reads from raw data. At the same time, Q30, GC-content, and sequence duplication level of the clean data were calculated. All the downstream analyses were based on clean data. These clean reads were then mapped to the reference genome sequence. Only reads with a perfect match (match ratio is 100%) or one mismatch were further analyzed and annotated based on the reference genome. Hisat2 tools (http://ccb.jhu.edu/software/hisat2/index.shtml) were used to map with reference genome of Arabidopsis.

2.7. GUS staining and imaging GUS staining was performed as described previously (Wang et al., 2019). Representative phenotypes were imaged using a Leica MZ16F stereomicroscope (Leica Microsystems, Wetzlar, Germany). 3. Results 3.1. Phenotypic observation of aging stigma First, we observed the changes of papilla cells during their whole lifespan, and used the key dPCD-related marker gene BFN1 as the reporter to monitor the dPCD process patterns. At stage 1 (defined as young stage), papilla cells started to develop, and initiated to elongate (Fig. 1A and B). Meanwhile, the GUS signal could not be detected at this time point, indicating that the dPCD process had not been started yet (Fig. 1G). At stage 2 (defined as mature stage), it was clear that large amount of papilla cells approached the mature phase, and the cell length and size were at the peaks (Fig. 1C and D). Very few of papilla

2.3. Identification of differentially expressed genes and gene clustering analysis Gene expression levels were estimated by fragments per kilobase of transcript per million fragments mapped (FPKM). Differential expression analysis was performed using the DEseq. DEseq provide statistical routines for determining differential expression in digital gene expression data using a model based on the negative binomial distribution. The resulting P-values were adjusted using the Benjamini-Hochberg 2

Plant Physiology and Biochemistry 147 (2020) 1–9

H. Ye, et al.

Fig. 1. Morphological observation of stigmatic papilla cells during their lifespan. A, C, and E, lateral view of papilla cells at three developmental stages; B, D, and F, top view of papilla cells at three developmental stages; G, H, and I, GUS staining images of papilla cells at three developmental stages. proBFN1::GUS-GFP reporter line was used to reveal whether programmed cell death process was happened. The representative images were shown. Scale bars: 0.5 mm.

genes (DEGs) in different sample, we established two comparative groups, one group was stage 1 vs stage 2, another one was stage 2 vs stage 3. In stage 1 vs stage 2 group, total 1235 genes showed differentially expressed patterns, and 528 DEGs were up-regulated, while 707 DEGs were down-regulated (Fig. 2A). In the group of stage 2 vs stage 3, there were 2347 genes displayed differentially expressed profiles. Among these DEGs, 683 DEGs were up-regulated, and 1664 DEGs showed down-regulated expression patterns (Fig. 2A). Then, we did a Venn diagram of all DEGs in these two comparative groups. As shown in Fig. 2B, 541 DEGs could be detected in both two comparative groups, and 694 and 1806 DEGs were expressed specifically in mature stage (stage 2 vs stage 1) and senescent stage (stage 3 vs stage 2), respectively. Given this, we filtered out 541 common DEGs, and put our eyesight on 1806 DEGs in senescent stage to find out the key genes.

cells showed weak GUS signal, suggesting that most of papilla cells were at mature stage, and would entrance into the senescent stage (Fig. 1H). At stage 3 (defined as senescent stage), a mass of papilla cells displayed wilted situations, only few of papilla cells were still at mature stage (Fig. 1E and F). After GUS staining, very strong GUS signal could be detected in the senescent papilla cells, revealing that almost all papilla cells were at the senescent stage (Fig. 1I). 3.2. Overview of the RNA-Seq data In order to check the quality of RNA-Seq data, the Pearson's Correlation Coefficient (PCC) among three biological replications of each sample was calculated. As shown in Fig. S1, the PCC values in stage 1 and stage 2 samples were higher than 0.97. In stage 3 samples, the PCC value showed a little lower, but still over 0.8. These results indicated that the three biological replications of each samples were highly interrelated. After filtering out adaptors and unknown or low-quality bases, more than 4,055,186,400 clean bases were obtained from each sample (Table S3). For each sample, the GC content was higher than 44.33%, and the Q30 value was higher than 89.16% (Table S3). Among 9 samples, the match ratio higher than 92.14%. Besides, the unique mapped read ratio was higher than 89.91% in each sample (Table S4).

3.4. GO enrichment of DEGs In order to identify the candidate DEGs involved in regulating the senescent process in papilla cells, first, we monitored the expression profiles of all DEGs. Based on the expression heatmap, we found that some DEGs showed very significant expression patterns (Fig. 3A). Therefore, we isolated 4 DEG clusters from all DEG clusters, and each DEG cluster displayed a clearly differential expression pattern between two different stages (stage 2 and stage 3) (Fig. 3B). In order to understand the functions of these DEGs, GO enrichment analysis was carried out. The top five GO terms of each cluster were shown (Fig. 4). During the senescent stage, DEGs in these 4 clusters mainly gathered into PCD-

3.3. Differentially expressed genes analysis To understand the expression patterns of differentially expressed 3

Plant Physiology and Biochemistry 147 (2020) 1–9

H. Ye, et al.

Fig. 2. Overview of DEGs in different comparable groups. A, the numbers of up- and down-regulated DEGs in comparable groups stage 1 vs stage 2 and stage 2 vs stage 3; B, venn diagram analysis of unique and common DEGs in comparable groups stage 1 vs stage 2 and stage 2 vs stage 3.

6 members were known as autophagy-related protein encoded genes, two members were annotated as serine/threonine-protein kinase ATG genes, and the last one was known as ubiquitin-like modifier-activating enzyme ATG gene. All DEGs showed up-regulated expression patterns (Table 2). Last, a sugar metabolism-related KEGG pathway, “starch and sucrose metabolism” (ko00500), were labeled out (Fig. S2). In order to learn more about this pathway, 15 DEGs were identified (Table S5).

related, phytohormone biosynthetic and energy-related GO terms, such as “positive regulation of programmed cell death” (GO: 0043068), “ethylene biosynthetic process” (GO: 0009693), “cadmium-transporting ATPase activity” (GO: 0015434), “regulation of salicylic acid biosynthesis” (GO: 0080142), and “indoleacetic acid biosynthetic process” (GO: 0009684) (Fig. 4). Based on the results of GO analysis, 13 DEGs, which showed the highest expression change folds in each term, were isolated, and the expression fold change of each DEG was over four folds. Only one gene encoded “aminocyclopropanecarboxylate oxidase” (AT1G77330) displayed a down-regulated expression trend, the rest of 12 DEGs were all up-regulated (Table 1).

3.6. Analysis of transcription factors (TFs) TFs have already been confirmed that involved in regulating PCD process of papilla cells (Gao et al., 2018). To reveal the types and amount of TFs involved in regulating PCD process of papilla cells, the common and important TFs were isolated. Totally, 9 types of TFs, including AP2/ERF (gene number: 18), ARF (gene number: 3), Dof (gene number: 3), C2H2 (gene number: 11), bHLH (gene number: 12), bZIP (gene number: 8), MYB (gene number: 14), WRKY (gene number: 17) and NAC (gene number: 18), were isolated from 1806 DEGs (Fig. 6A). During the senescent stage of papilla cells, most of TFs were up-regulated. For example, 16 out of 18 AP2/ERF were un-regulated; for ARF and MYB, 50% of these TFs showed increase expression trends; for C2H2, bHLH, bZIP, NAC, over 60% of these TFs were expressed increasingly; all of Dof and WRKY showed up-regulated expression profiles (Fig. 6B). The information of all TFs were listed in Table S2.

3.5. KEGG pathway analysis of DEGs To gain more information into these 1806 DEGs, KEGG pathway analysis was performed. The top 20 KEGG pathways were shown in Fig. S2. In accord with the results of GO analysis, the first enriched KEGG pathway was “plant hormone signal transduction” (ko04075) (Fig. S2). In order to understand which genes were involved in this pathway, we monitored the expression profiles of DEGs participated in auxin, ethylene (ET), jasmonic acid (JA), and salicylic acid (SA) signaling pathways (Fig. 5A). In auxin signaling pathway, four auxin-responsive gene (IAAs) and one auxin response factor (ARF) gene were identified, and their expression trends were changed clearly during the senescent period (Fig. 5B). In ET signaling pathway, two genes, which encoded the ethylene receptor (ETR) and ethylene-responsive transcription factor (ERF), respectively, displayed differently expression patterns (Fig. 5B). For JA signaling pathway, a Jasmonate ZIM domain (JAZ) gene, JAZ5, was considered as a candidate gene (Fig. 5B). In SA signal transduction pathway, two NONEXPRESSER OF PR GENES 1-like (NPR1like) genes and one SA-related transcription factor (TGA) showed the significant expression trends (Fig. 5B). Then, we noticed the KEGG pathway “regulation of autophagy” (ko04140), which was highly closed with dPCD process (Fig. S2). Totally, 9 DEGs, which annotated as autophagy-related genes, were confirmed (Table 2). Among these DEGs,

3.7. DEG validation by qRT-PCR In order to verify the reliability of the RNA-Seq data, qRT-PCR was utilized on 6 randomly selected genes from 104 TFs (Table S2). The log2Fold Change (log2FC) values, which represented the RNA-Seq and qRT-PCR results, were shown in Fig. S3. The expression trend of each selected TF in qRT-PCR assay was similar as that in RNA-Seq analysis (Fig. S3). These results indicated that the RNA-Seq data in this study were reliable. 4

Plant Physiology and Biochemistry 147 (2020) 1–9

H. Ye, et al.

thought plant hormone should be involved in papilla cell death. So far, there are no studies show some evidences that auxin is involved in regulating the PCD of papilla cells. The studies of auxinmediated PCD are mainly focused on root development. In Arabidopsis, the growth of lateral root primordium (LRP) and the remodeling of LRPoverlying cells rely to auxin (Porco et al., 2016). The PCD of lateral root cap (LRC) cells was shown to result in release of auxin (Xuan et al., 2016). During the chilling stress, the protective cell death of columella stem cell daughters will be induced, and this process is closely related with auxin distribution in root stem cells (Hong et al., 2017; Tombuloglu, 2019). In this study, we isolated 4 IAA genes and one ARF gene using KEGG pathway analysis (Fig. 5B). These genes are essential components of auxin signaling pathway (Fig. 5A), and they showed clearly differential expression patterns during papilla cell senescence (Fig. 5B). Given this, we thought that auxin might play a role in regulating papilla cell death. Ethylene (ET) regulates a variety of processes in plants, such as seed germination, ripening of fruits, stress responses, tissue senescence, and so on (Bleecker and Kende, 2000). In Arabidopsis leaf senescence, a number of transcription factors that mediate ethylene responses are upregulated (Chen et al., 2002). Among the ethylene-responsive element binding factors (ERFs), AtERF1, 2, and 5 are upregulated in response to exogenous ethylene. Besides, the relationship between ET and the senescence of flower petals has been established many years ago (Reid and Wu, 1992). In present study, we identified one ERF gene and one ETR gene showing differential expression profiles in ET signaling pathway (Fig. 5). This result indicated that ET, as an important senescent plant hormone, might be involved in regulating papilla cell death. JA is known as a PCD-associated phytohormone recently. Some previous studies have already revealed that JA induces leaf senescence, which a kind of PCD in leaf lifespan, and altered the expression profiles of various senescence-associated genes in some plant species (Reinbothe et al., 2009). Based on the KEGG analysis, we isolated a jasmonate-ZIM-domain protein encoded gene, JAZ5, which displayed a clear differential expression pattern (Fig. 5). This result revealed that JA might take part in the stigmatic papilla cell death by JA signaling pathway. SA is implicated in the induction of PCD-associated with pathogen defense responses due to SA levels increase in response to PCD-inducing infections (Brodersen et al., 2005). SA has been shown to be required for spontaneous cell death in several lesion-mimic mutants (Brodersen et al., 2005). NONEXPRESSOR OF PATHOGENESIS-RELATED GENES1 (NPR1), acting as a SA receptor, can promote SA-induced defense gene expression and pathogen resistance (Fan and Dong, 2002). In our work, two NPR1-like genes were isolated, and they displayed altered expression levels during papilla cell death process (Fig. 5). Moreover, we found a downstream gene of NPR1, named as TGA, also showed a changed expression profile (Fig. 5). These results demonstrated that SA might be involved in the PCD process of papilla cells. In conclusion, auxin, JA, ET and SA were considered to be involved in stigmatic papilla cell death. It is better to use the specific reporters to visualize the spatial-temporal distribution of phytohormone or monitor the expression manners of key components in plant hormone signaling pathways. Anyway, our eyesight on roles of plant hormone in dPCD should be extended into the floral organ papilla cells.

Fig. 3. Expression patterns of 1806 DEGs in comparable group stage 2 vs stage 3. A, expression heatmap of 1806 DEGs in comparable group stage 2 vs stage 3; B, expression heatmap of four clusters of DEGs showed significant expression trends.

4. Discussion As a huge mystery in the journey of life, death is known as a complex process in species, especially in plants. Recent years, researchers are beginning to study some of the processes involved in developmental programmed cell death (dPCD) in plants, but we are still on the road now. Stigma is the carrier of stigmatic papilla cells, and is an essential organ for the pollination process. Acting as the receptor of pollen, papilla cells play crucial role in the determination of reproductive success. Therefore, uncovering of the regulatory mechanism of the stigmatic PCD process is of great significance for the seed production and yield. Although a previous study has already revealed two key players, KIRA1 and ORE1, in the PCD of stigma (Gao et al., 2018), it is still necessary to gain more details for the perfection of the stigmatic PCD. In this study, we further obtained some information of the stigmatic PCD using transcriptome analysis, and would like to discuss as follows:

4.2. The roles of autophagy in regulating papilla cell death

4.1. Phytohormone signals play important roles in regulating papilla cell lifespan

Autophagy is known as a kind of PCD, and it is a major cellular degradation pathway in eukaryotes. Recent studies have revealed the importance of autophagy in many aspects of plant life, such as seedling establishment, plant development, stress resistance, metabolism, reproduction, and so on (Michaeli et al., 2016). The autophagy-related ATG1/13 kinase complex has been characterized in Arabidopsis (Li et al., 2014). ATG1 was shown to interact with ATG13 and ATG8. Although ATG proteins are conserved in higher eukaryotes, there are

The functions of plant hormone on the regulation of plant growth and development are undoubted. Recently, the functions of plant hormone on modulating plant tissues senescence or death are uncovered. Previous studies demonstrated that the main plant hormones, such as ET, auxin, JA and SA, play crucial roles on the tissue-specific PCD process. In our KEGG pathway analysis results, the “plant hormone signal transduction” (ko04075) raised our concern. Therefore, we 5

Plant Physiology and Biochemistry 147 (2020) 1–9

H. Ye, et al.

Fig. 4. Gene ontology (GO) classification of four clusters of DEGs. Top five GO terms were shown for each cluster of DEGs. The numbers of DEGs clustered in each GO term were shown with black values.

Table 1 Information of DEGs isolated from GO enrichment analysis. Gene ID

Log2FCa

Expression pattern

FDR

Annotation

AT1G51800 AT1G52200 AT1G72520 AT1G77330 AT1G79320 AT1G79330 AT2G19110 AT2G32020 AT2G32030 AT3G01830 AT3G11840 AT5G26920 AT5G27420

3.88 3.58 4.79 −3.70 4.88 4.79 4.34 2.79 4.40 4.74 3.52 3.12 3.67

Up Up Up Down Up Up Up Up Up Up Up Up Up

4.18E-07 5.60E-09 3.53E-28 2.26E-10 2.45E-07 4.30E-06 3.80E-04 1.52E-04 5.68E-12 3.99E-03 1.68E-07 5.17E-07 4.57E-09

LRR receptor-like serine/threonine-protein kinase IOS1 PLANT CADMIUM RESISTANCE 8 Lipoxygenase 4 Aminocyclopropanecarboxylate oxidase Metacaspase 6 Metacaspase 5 Putative cadmium/zinc-transporting ATPase HMA4 Acyl-CoA N-acyltransferases (NAT) superfamily protein Acyl-CoA N-acyltransferases (NAT) superfamily protein Probable calcium-binding protein CML40 E3 ubiquitin-protein ligase PUB24 Calmodulin-binding protein 60 G E3 ubiquitin-protein ligase ATL31

a

Log2FC = Log2 Fold Change.

4.3. The network of TFs is essential for dPCD regulation of stigma longevity

limited studies working on the roles of autophagy in plants. Recently, a study showed that mutations in AtATG4 gene block formation of autophagic bodies and result in the arrest of root growth under nutrient starvation (Yoshimoto et al., 2004). However, the functions of ATG genes and autophagy in stigmatic papilla cell death has still not been reported. Here, based on the KEGG pathway analysis results, we noticed another pathway “regulation of autophagy” (ko04140). Meanwhile, we isolated 9 autophagy-related genes listed in Table 2. The expression levels of these genes showed clearly up-regulated trends during papilla cell death (Table 2), indicating that some autophagy-related genes took part in regulating papilla cell death directly or indirectly. Gao and colleagues have already showed the key factors, KIRA1 and ORE1, in death process of papilla cells (Gao et al., 2018), it was possible that these 9 genes isolated in this study were downstream targets of KIR1 and ORE1 proteins. It should be worthy to observe the autophagy process during papilla cell death in future work.

There is no doubt that TFs participate in the regulation of PCD in plants. Based on different tissues in plants, the participant TFs are different correspondingly. In the tapetal PCD in Arabidopsis thaliana, MYB80 is required for pollen and tapetum development (Phan et al., 2012). In rice, the bHLH142 mutant shows retarded meiosis and defects in tapetal PCD. bHLH164, bHLH5, and bHLH141 are known to function in rice pollen development. bHLH141 positively regulates the expression of AP37 and AP25, which induce tapetal PCD death (Ko et al., 2014). The active PCD process is integral to Arabidopsis root cap differentiation. During this process, ANAC033/SOMBRERO makes a decisive role (Fendrych et al., 2014). During the endoplasmic reticulum (ER)-stress-induced PCD in Arabidopsis, NAC089 is the key factor that can promote this process. bZIP28 and bZIP60 directly control the expression pattern of NAC089 (Yang et al., 2014). In low temperaturedependent cell death of Arabidopsis, a C2H2 TF, LESION SIMULATING 6

Plant Physiology and Biochemistry 147 (2020) 1–9

H. Ye, et al.

Fig. 5. DEGs related to plant hormone signal transduction pathway. A, diagram of auxin, ethylene, jasmonic acid, and singicylic acid signaling transduction pathways; B, information and expression patterns of DEGs involved in auxin, ethylene, jasmonic acid, and singicylic acid signaling transduction pathways.

DISEASE RESISTANCE 1 (LSD1), plays an important role (Huang et al., 2010). In this study, 9 types, total of 104 TFs were isolated (Fig. 6 and Table S2). Most of these TFs were up-regulated, indicating they were promoted during papilla cell death. Besides, we found that the numbers of NAC and WRKY were the highest among 9 types of TFs, revealing that these two kinds of TFs might play important roles in regulating papilla cell death. Consistent with this, Gao and colleagues found two NAC TFs, KIRA1 and ORE1, play core roles in PCD of papilla cells (Gao et al., 2018). Given these evidences, we speculated that the TFs isolated in this study might be the targets of KIRA1 and ORE1, and they might be form a network with KIRA1 and ORE1 to control papilla cell lifetime.

can be induced by various factors, including developmental signals and environmental cues (Xu and Hanson, 2000). Some intracellular cues, such as the production and transportation of sucrose, ATP generation, can also induce the PCD process (Ruan, 2012). In present study, we isolated 15 sugar metabolism-related DEGs, including UDP-glucuronate 4-epimerase, sucrose synthase, endoglucanase, and so on (Table S5). Meanwhile, two acyl-CoA N-acyltransferase-encoded genes and one cadmium/zinc-transporting ATPase-encoded gene were also isolated by GO enrichment analysis (Table 1), indicating that intracellular energyrelated pathway might be involved in the PCD of papilla cells.

4.4. Other probable factors in stigma senescent process

PCD in stigmatic papilla cells provides us a good model system to study the molecular aspects of tissue (or organ) senescence. Compared with other tissues, Arabidopsis papilla cells have a relatively shorter lifespan. Benefit from the number and shape, papilla cell death is easy to monitor. In the study reported herein, we try to build a panoramic network for the regulation of papilla cell death. By combining the information of DEGs, GO terms, KEGG pathways gained from RNA-Seq analysis, a supposed working model of papilla cell death was showed in Fig. 7. The possible way was the network linked by plant hormone signals and TFs could regulate the intracellular energy-related processes and PCD/autophagy pathways (Fig. 7). We speculated that the plant hormone levels were unbalanced, the energy support was lacking, and the TF-mediated PCD process, such as autophagy, was triggered in the

5. Conclusions

Excepting the factors mentioned above, we also found some other players in the DEG list (Table 1). First, we noticed two matacaspaseencoded genes. Despite the absence of canonical caspases in plants, dying tissue cells show an increased proteolytic caspase-like activity (Bozhkov et al., 2004). Studies have proved that matacaspase is essential in plant PCD (Aravind and Koonin, 2002; Bostancioglu et al., 2018). Metacaspase-dependent PCD is crucial for plant embryogenesis (Suarez et al., 2004). Given this, we thought that these two metacaspase-encoded genes, which showed increase and high log2FC values (Table 1), might be involved in the regulation of papilla cell death. PCD is also an integral part of the life cycle of plants. This process Table 2 Information of DEGs involved in autophagy pathway. Gene ID

Log2FCa

Expression pattern

FDR

Annotation

AT2G45170 AT3G06420 AT3G15580 AT3G49590 AT3G53930 AT3G61960 AT4G04620 AT4G21980 AT5G45900

2.59 2.13 1.09 1.82 1.96 2.46 1.19 2.29 2.68

Up Up Up Up Up Up Up Up Up

3.20E-04 4.43E-05 7.27E-04 3.00E-04 3.01E-06 6.35E-06 2.86E-04 4.24E-04 4.47E-03

Autophagy-related protein 8e Autophagy-related protein 8 h Autophagy-related protein 8i Autophagy-related protein 13a Serine/threonine-protein kinase ATG1b Serine/threonine-protein kinase ATG1a Autophagy-related protein 8b Autophagy-related protein 8a Ubiquitin-like modifier-activating enzyme ATG7

a

Log2FC = Log2 Fold Change. 7

Plant Physiology and Biochemistry 147 (2020) 1–9

H. Ye, et al.

Fig. 7. A putative network in regulation of papilla cell death. Based on the results in this study, phytohormone signals and TFs are interacted each other. The integration of phytohormone signals and TFs will regulate the energy-related pathways, such as sucrose biosynthesis or transduction and ATP biosynthesis or metabolism. The final phenomenon is the PCD process (including autophagy), which regulated by plant hormones, TFs, and energy-related pathways, was triggered in senescent papilla cells.

manuscript. HY, FR, YW and JB revised the manuscript. All authors read and approved the final manuscript. Declaration of competing interest The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper. Acknowledgements We apologize to the researchers for not citing their nice works, due to the limited layout of a printed page. This study was supported by the National Science Foundation for Young Scientists of China (31901490) and the Beijing Academy of Agriculture and Forestry Sciences Foundation for Youths (QNJJ201916). Appendix A. Supplementary data Supplementary data to this article can be found online at https:// doi.org/10.1016/j.plaphy.2019.12.008. References Aravind, L., Koonin, E.V., 2002. Classification of the caspasehemoglobinase fold: detection of new families and implications for the origin of the eukaryotic separins. Proteins 46, 355–367. https://doi.org/10.1002/prot.10060. Bleecker, A.B., Kende, H., 2000. Ethylene: a gaseous signal molecule in plants. Annu. Rev. Cell Dev. Biol. 16 (1), 1–18. https://doi.org/10.1146/annurev.cellbio.16.1.1. Bostancioglu, S.M., Tombuloglu, G., Tombuloglu, H., 2018. Genome-wide identification of barley MCs (metacaspases) and their possible roles in boron-induced programmed cell death. Mol. Biol. Rep. 45 (3), 211–225. https://doi.org/10.1007/s11033-0184154-3. Bozhkov, P.V., Filonova, L.H., Suarez, M.F., Helmersson, A., Smertenko, A.P., Zhivotovsky, B., Von Arnold, S., 2004. VEIDase is a principal caspase-like activity involved in plant programmed cell death and essential for embryonic pattern formation. Cell Death Differ. 11, 175–182. https://doi.org/10.1038/sj.cdd.4401330. Brodersen, P., Malinovsky, F.G., Hématy, K., Newman, M.A., Mundy, J., 2005. The role of salicylic acid in the induction of cell death in Arabidopsis acd11. Plant Physiol. 138, 1037–1045. https://doi.org/10.1104/pp.105.059303. Chen, W., Provart, N.J., Glazebrook, J., Katagiri, F., Chang, H.S., Eulgem, T., Mauch, F., Luan, S., Zou, G., Whitham, S.A., Budworth, P.R., Tao, Y., Xie, Z., Chen, X., Lam, S., Kreps, J.A., Harper, J.F., Si-Ammour, A., Mauch-Mani, B., Heinlein, M., Kobayashi, K., Hohn, T., Dangl, J.L., Wang, X., Zhu, T., 2002. Expression profile matrix of Arabidopsis transcription factor genes suggests their putative functions in response to

Fig. 6. TFs involved in papilla cell death. A, types and numbers of TFs related to papilla cell death; B, the ratios of up- and down-regulated TFs in each type.

ageing papilla cells. The expression levels of genes isolated in this study were changed significantly, indicating these genes, which were involved in the pathway or network mentioned above, were important factors in regulating papilla cell death. One word, our study find out some useful candidate genes and TFs, which can be utilized to perfect the regulatory mechanism of papilla cell death. Author's contributions YW and JB conceived the study and designed the experiments. HG and LG cultivated and prepared the plant materials; HY and FR carried out the experiments and analyzed the data; HY and YW wrote the 8

Plant Physiology and Biochemistry 147 (2020) 1–9

H. Ye, et al.

37–43. https://doi.org/10.1007/BF00024431. Reinbothe, C., Springer, A., Samol, I., Reinbothe, S., 2009. Plant oxylipins: role of jasmonic acid during programmed cell death, defence and leaf senescence. FEBS J. 276, 4666–4681. https://doi.org/10.1111/j.1742-4658.2009.07193.x. Ruan, Y.L., 2012. Signaling role of sucrose metabolism in development. Mol Plant. Jul 5 (4), 763–765. https://doi.org/10.1093/mp/sss046. Shibuya, K., Yamada, T., Ichimura, K., 2016. Morphological changes in senescing petal cells and the regulatory mechanism of petal senescence. J. Exp. Bot. 67, 5909–5918. https://doi.org/10.1093/jxb/erw337. Suarez, M.F., Filonova, L.H., Smertenko, A., Savenkov, E.I., Clapham, D.H., von Arnold, S., Zhivotovsky, B., Bozhkov, P.V., 2004. Metacaspase-dependent programmed cell death is essential for plant embryogenesis. Curr. Biol. 14 (9), R339–R940. https:// doi.org/10.1016/j.cub.2004.04.019. Tomas, H., 2013. Senescence, ageing and death of the whole plant. New Phytol. 197, 696–711. https://doi.org/10.1111/nph.12047. Tombuloglu, H., Kekec, G., Sakcali, M.S., Unver, T., 2013. Transcriptome-wide identification of R2R3-MYB transcription factors in barley with their boron responsive expression analysis. Mol. Genet. Genom. 288 (3–4), 141–155. https://doi.org/10.1007/ s00438-013-0740-1. Tombuloglu, G., Tombuloglu, H., Sakcali, M.S., Unver, T., 2015. High-throughput transcriptome analysis of barley (Hordeum vulgare) exposed to excessive boron. Gene 557 (1), 71–81. https://doi.org/10.1016/j.gene.2014.12.012. Tombuloglu, H., 2019. Genome-wide analysis of the auxin response factors (ARF) gene family in barley (Hordeum vulgare L.). J. Plant Biochem. Biotechnol. 28, 14. https:// doi.org/10.1007/s13562-018-0458-6. Ueda, H., Kusaba, M., 2015. Strigolactone regulates leaf senescence in concert with ethylene in Arabidopsis. Plant Physiol. 169 (1), 138–147. https://doi.org/10.1104/ pp.15.00325. Van Doorn, W.G., Woltering, E.J., 2008. Physiology and molecular biology of petal senescence. J. Exp. Bot. 59 (3), 453–480. https://doi.org/10.1093/jxb/erm356. Wang, Y.K., Duan, W.J., Bai, J.F., Wang, P., Yuan, S.H., Zhao, C.P., Zhang, L.P., 2019. Constitutive expression of a wheat microRNA, TaemiR167a, confers male sterility in transgenic Arabidopsis. Plant Growth Regul. 88, 227. https://doi.org/10.1007/ s10725-019-00503-4. Williams, R.R., 1965. The effect of summer nitrogen applications on the quality of apple blossom. J. Hortic. Sci. 40, 31–41. https://doi.org/10.1080/00221589.1965. 11514118. Xu, Y., Hanson, M.R., 2000. Programmed cell death during pollination-induced petal senescence in petunia. Plant Physiol. 122, 1323–1333. https://doi.org/10.2307/ 4279215. Xuan, W., Band, L.R., Kumpf, R.P., Van Damme, D., Parizot, B., De Rop, G., Opdenacker, D., Möller, B.K., Skorzinski, N., Njo, M.F., De Rybel, B., Audenaert, D., Nowack, M.K., Vanneste, S., Beeckman, T., 2016. Cyclic programmed cell death stimulates hormone signaling and root development in Arabidopsis. Science 351 (6271), 384–387. https:// doi.org/10.1126/science.aad2776. Yang, Z.T., Wang, M.J., Sun, L., Lu, S.J., Bi, D.L., Sun, L., Song, Z.T., Zhang, S.S., Zhou, S.F., Liu, J.X., 2014. The membrane-associated transcription factor NAC089 controls ER-stress-induced programmed cell death in plants. PLoS Genet. 10 (3), e1004243. https://doi.org/10.1371/journal.pgen.1004243. Yoshimoto, K., Hanaoka, H., Sato, S., Kato, T., Tabata, S., Noda, T., Ohsumi, Y., 2004. Processing of ATG8s, ubiquitin-like proteins, and their deconjugation by ATG4s are essential for plant autophagy. Plant Cell 16, 2967–2983. https://doi.org/10.1105/ tpc.104.025395.

environmental stresses. Plant Cell 14 (3), 559–574. https://doi.org/10.1105/tpc. 010410. Daneva, A., Gao, Z., Van Durme, M., Nowack, M.K., 2016. Functions and regulation of programmed cell death in plant development. Annu. Rev. Cell Dev. Biol. 32, 441–468. https://doi.org/10.1146/annurev-cellbio-111315-124915. Dickman, M., Williams, B., Li, Y., Figueiredo, P., Wolpert, T., 2017. Reassessing apoptosis in plants. Nat. Plants 3, 773–779. https://doi.org/10.1038/s41477-017-0020-x. Fan, W., Dong, X., 2002. In vivo interaction between NPR1 and transcription factor TGA2 leads to salicylic acid-mediated gene activation in Arabidopsis. Plant Cell 14, 1377–1389. https://doi.org/10.1105/tpc.001628. Fendrych, M., Van Hautegem, T., Van Durme, M., Olvera-Carrillo, Y., Huysmans, M., Karimi, M., Lippens, S., Guérin, C.J., Krebs, M., Schumacher, K., Nowack, M.K., 2014. Programmed cell death controlled by ANAC033/SOMBRERO determines root cap organ size in Arabidopsis. Curr. Biol. 24 (9), 931–940. https://doi.org/10.1016/j.cub. 2014.03.025. Gao, Z., Daneva, A., Salanenka, Y., Van Durme, M., Huysmans, M., Lin, Z., De Winter, F., Vanneste, S., Karimi, M., Van De Velde, J., Vandepoele, K., Van de Walle, D., Dewettinck, K., Lambrecht, B.N., Nowack, M.K., 2018. KIRA1 and ORESARA1 terminate flower receptivity by promoting cell death in the stigma of Arabidopsis. Nat Plants. Jun 4 (6), 365–375. https://doi.org/10.1038/s41477-018-0160-7. Hong, J.H., Savina, M., Du, J., Devendran, A., Kannivadi Ramakanth, K., Tian, X., Sim, W.S., Mironova, V.V., Xu, J., 2017. A sacrifice-for-survival mechanism protects root stem cell niche from chilling stress. Cell 170 (1), 102–113. https://doi.org/10.1016/j. cell.2017.06.002. Huang, X.Z., Li, Y.S., Zhang, X.Y., Zuo, J.R., Yang, S.H., 2010. The Arabidopsis LSD1 gene plays an important role in the regulation of low temperature-dependent cell death. New Phytol. 187 (2), 301–312. https://doi.org/10.2307/40792379. Huysmans, M., Lema, A.S., Coll, N.S., Nowack, M.K., 2017. Dying two deaths –programmed cell death regulation in development and disease. Curr. Opin. Plant Biol. 35, 37–44. https://doi.org/10.1016/j.pbi.2016.11.005. Ko, S.S., Li, M.J., Ku, M.S.B, Ho, Y.C., Lin, Y.J., 2014. The bHLH142 transcription factor coordinates with TDR1 to modulate the expression of EAT1 and regulate pollen development in rice. Plant Cell 26, 2486–2504. https://doi.org/10.1105/tpc.114. 126292. Li, F., Chung, T., Vierstra, R.D., 2014. AUTOPHAGY-RELATED11 plays a critical role in general autophagy- and senescence-induced mitophagy in Arabidopsis. Plant Cell. Feb 26 (2), 788–807. https://doi.org/10.4161/auto.29320. Michaeli, S., Galil, i G., Genschik, P., Fernie, A.R., Avin-Wittenberg, T., 2016. Autophagy in plants–What's new on the menu? Trends Plant Sci. 21 (2), 134–144. https://doi. org/10.1016/j.tplants.2015.10.008. Olvera-Carrillo, Y., Van Bel, M., Van Hautegem, T., Fendrych, M., Simaskova, M., Van Durme, M., Buscaill, P., Rivas, S., Coll, N.S., Coppens, F., Maere, S., Nowack, M.K., 2015. A conserved core of programmed cell death indicator genes discriminates developmentally and environmentally induced programmed cell death in plants. Plant Physiol. 169, 2684–2699. https://doi.org/10.1016/10.1104/pp.15.00769. Phan, H.A., Li, S.F., Parish, R.W., 2012. MYB80, a regulator of tapetal and pollen development, is functionally conserved in crops. Plant Mol. Biol. Report. 78, 171–183. https://doi.org/10.1007/s11103-011-9855-0. Porco, S., Larrieu, A., Du, Y., Gaudinier, A., Goh, T., Swarup, K., Kuemoers, B., 2016. Lateral root emergence in Arabidopsis is dependent on transcription factor LBD29 regulating auxin influx carrier LAX3. Development 143, 3340–3349. https://doi.org/ 10.1242/dev.136283. Reid, M.S., Wu, M.J., 1992. Ethylene and flower senescence. Plant Growth Regul. 11,

9