RNA-seq Gene Profiling Reveals Transcriptional Changes in the Late Phase during Compatible Interaction between a Korean Soybean Cultivar (Glycine max cv. Kwangan) and Pseudomonas syringae pv. syringae B728a

Article information

Plant Pathol J. 2022;38(6):603-615
Publication date (electronic) : 2022 December 1
doi : https://doi.org/10.5423/PPJ.OA.08.2022.0118
1Department of Applied Bioscience, Dong-A University, Busan 49315, Korea
2Department of Molecular Genetics, Dong-A University, Busan 49315, Korea
3Department of Biological Sciences, Chonnam National University, Gwangju 61186, Korea
4Institute of Agricultural Life Science, Dong-A University, Busan 49315, Korea
*Co-corresponding authors: H. J. Park, Phone) +82-62-530-3410, FAX) +82-62-530-3409, E-mail) hjpark_bio@jnu.ac.kr. H. W. Jung, Phone) +82-51-200-7546, FAX) +82-51-200-7505, E-mail) hwjung@dau.ac.kr
†These authors contributed equally to this work.Handling Editor: Sang-Wook Han
Received 2022 August 25; Revised 2022 September 23; Accepted 2022 September 26.

Abstract

Soybean (Glycine max (L) Merr.) provides plant-derived proteins, soy vegetable oils, and various beneficial metabolites to humans and livestock. The importance of soybean is highly underlined, especially when carbon-negative sustainable agriculture is noticeable. However, many diseases by pests and pathogens threaten sustainable soybean production. Therefore, understanding molecular interaction between diverse cultivated varieties and pathogens is essential to developing disease-resistant soybean plants. Here, we established a pathosystem of the Korean domestic cultivar Kwangan against Pseudomonas syringae pv. syringae B728a. This bacterial strain caused apparent disease symptoms and grew well in trifoliate leaves of soybean plants. To examine the disease susceptibility of the cultivar, we analyzed transcriptional changes in soybean leaves on day 5 after P. syringae pv. syringae B728a infection. About 8,900 and 7,780 differentially expressed genes (DEGs) were identified in this study, and significant proportions of DEGs were engaged in various primary and secondary metabolisms. On the other hand, soybean orthologs to well-known plant immune-related genes, especially in plant hormone signal transduction, mitogen-activated protein kinase signaling, and plant-pathogen interaction, were mainly reduced in transcript levels at 5 days post inoculation. These findings present the feature of the compatible interaction between cultivar Kwangan and P. syringae pv. syringae B728a, as a hemibiotroph, at the late infection phase. Collectively, we propose that P. syringae pv. syringae B728a successfully inhibits plant immune response in susceptible plants and deregulates host metabolic processes for their colonization and proliferation, whereas host plants employ diverse metabolites to protect themselves against infection with the hemibiotrophic pathogen at the late infection phase.

Soybean (Glycine max (L) Merr.), an edible legume, is one of the most important crops as an excellent source of plant-based oils, vegetable proteins, and secondary metabolites. Wild soybean (G. soja) was domesticated to G. max in Northeast Asia, including Korea, China, and Japan. So far, over 60,000 varieties have been developed and cultivated in different countries (Carter et al., 2004; Li et al., 2020). Since a reference genome of cultivated accession Williams 82 was published, functional studies of soybean are accelerating, and the pan-genome of soybean was also analyzed (Liu et al., 2015b, 2020; Schmutz et al., 2010; Wang and Tian, 2015). A Korean soybean cultivar (cv.), Kwangan, also known as ‘Kwangankong’ (Suwon No. 159), was bred by crossing a high-protein Japanese cv. Dongsan (Touzan) 69 with a virus disease-resistant Korean cv. Hwangkeum (Suwon No. 97), and largely cultivated in the Central and Southern regions of Korea (Kim et al., 2011; Lee et al., 2015; Soybean Breeding Team, Upland Crop Div., Crop Experiment and Experiment Station, 1994). The soybean cv. Kwangan exhibits good crop traits such as high protein content and seed yield and relatively higher resistance against pests and pathogens than other Korean sprout soybean cultivars (Kim et al., 2011; Soybean Breeding Team, Upland Crop Div., Crop Experiment and Experiment Station, 1994).

Like other crops and vegetables, the productivity and quality of soybean rely on disease management. Indeed, the average soybean yield loss due to pests and pathogens was 21.4% (11.0–32.4%) from 2010 to 2014 globally (Savary et al., 2019). Worldwide, cyst nematode (Heterodera glycines), white mold (Sclerotinia sclerotiorum), and soybean rust (Phakopsora pachyrhizi) are the most destructive pathogens in soybean, and many fungi, oomycetes, and nematodes also cause severe yield loss of soybean (Savary et al., 2019). Additionally, bacterial blight and bacterial pustule, the major bacterial diseases in fields, are caused by Pseudomonas syringae pv. glycinea and Xanthomonas axonopodis pv. glycines, respectively (Faske et al., 2014; Hartman et al., 2015). Even though agricultural practices such as aggressive usage of agricultural chemicals, crop rotation, and tillage, have been performed to avoid primary and secondary inoculum from previously infected crop residues, disease management of soybean is primarily dependent on cultivating disease-resistant plants (Faske et al., 2014). Thus, considering eco-friendly and sustainable food supply, previously unidentified genetic resources have to be unveiled and applied to molecular breeding to increase disease resistance in soybean plants.

P. syringae pv. syringae among more than 60 pathovars of P. syringae exhibits polyphagous bacterium with a broad-spectrum host range (Kennelly et al., 2007). P. syringae pv. syringae can infect soybean leaves and pods, and consequently cause necrotic lesions surrounded by chlorotic zones (Gnanamanickam and Ward, 1982). The genome sequence of P. syringae pv. syringae B728a (PssB728a), a representative strain of P. syringae pv. syringae that was firstly isolated from a snap bean (Phaseolus vulgaris) leaflet in Wisconsin, USA, was complete (Feil et al., 2005; Loper and Lindow, 1987). PssB728a retains an effector repertoire secreted into its host plants through its type III secretion system (Vinatzer et al., 2006). In the field, this bacterium can multiply and survive for many generations on the (non)host leaf surfaces as an epiphyte (Ercolani et al., 1974; Hirano and Upper, 1990; Lindemann et al., 1984). The size of the epiphytic population increases upon heavy rains, promoting the bacteria’s rapid multiplication (Hirano et al., 1994, 1996). Thus, it is considered an intriguing source to study the molecular interaction between host plants and parasitic microbes (Vinatzer et al., 2006).

Recently, (molecular) components of soybeans associated with biotic stress response have been identified (Cooper et al., 2013; Delgado-Cerrone et al., 2018; Helm et al., 2019; Liu et al., 2015b; Pottinger et al., 2020; Russell et al., 2015; Wei et al., 2020; Yuan et al., 2020). Thus, we felt it is necessary to fill the knowledge gap in a Korean soybean accession. In this study, we evaluated the bacterial disease in the soybean cv. Kwangan inoculated with PssB728a and analyzed transcriptional change in soybean plants to identify genes involved in pathogenesis and immune response against a virulent Pseudomonas infection. In the infected soybean leaves, 8,902 and 7,786 genes were up- and down-regulated, respectively, at the late phase of infection. They appeared to be associated and involved in mainly biological pathways, including metabolic pathways, secondary metabolites, hormonal signal transduction, and plant-pathogen interaction upon PssB728a bacterial infection. Taken together, this pathosystem shows that the expression of soybean orthologs to well-known plant immune-related genes, especially in plant hormone signal transduction, MITOGEN-ACTIVATED PROTEIN KINASE (MAPK) signaling, and plant-pathogen interaction, is reduced during a compatible interaction. In contrast, soybean genes involved in the biosynthesis of primary and secondary metabolites are strongly expressed by PssB728a infection. This finding can help understand how PssB728a suppresses the plant immunity in the susceptible cultivar. We propose a few signaling pathways primarily affected by a virulent bacterial infection during a compatible interaction.

Materials and Methods

Plant material and growth conditions

G. max cv. Kwangan was used for this study. Soybean seeds were sterilized with chlorine gas (100 ml of 5.25% NaOCl supplemented with 3 ml of 12 N HCl) in a fume-hood overnight and sowed on soils (Nongwoo Bio, Suwon, Korea). Soybean plants grew under long-day conditions (25°C, 16 h light/8 h dark, 120 μmol/m2/s) in a walk-in growth chamber.

Bacterial strain and inoculation

The second trifoliate leaves of 3-week-old soybean plants were inoculated with PssB728a. We sprayed PssB728a suspension (OD600 = 0.01) in 10 mM MgSO4 containing 0.01% Silwet L-77 to the abaxial and adaxial side of leaves and kept the infected plants in the dark conditions with high humidity during the indicated periods. At least eight leaf discs were taken from infected leaves, and the number of bacteria in the leaf disc was counted using the typical serial dilution method. To analyze the epiphytic growth of PssB728a in plants, leaf discs of infected leaves with PssB728a were washed thoroughly with 10 mM MgSO4, in which bacterial growth was counted. The washed leaves were then soaked in 70% ethanol during several seconds for surface sterilization and were ground in 10 mM MgSO4 to measure the endophytic growth of PssB728a (Lee et al., 2012).

Transcriptome analysis

Trifoliate leaf tissues of soybean plants (~1 g) were harvested on the 5 days after spray-inoculation with 10 mM MgSO4 (2 of the control group) and PssB728a (OD600 = 0.01) (2 of treatment group) for the extraction of total RNA. RNA purification using RNeasy spin column (Qiagen, Hilden, Germany) was then conducted following the manufacturer’s instruction. After confirming the RNA quality (Bioanalyzer, Agilent, Palo Alto, CA, USA), the total RNA was subject to construct mRNA sequencing libraries using a TruSeq Stranded Total RNA LT Sample Prep Kit (Plant) according to the manufacturer’s instruction (Illumina, San Diego, CA, USA) (Martin and Wang, 2011). Each library from two biological replicates per treatment (total 4 libraries) was sequenced for 101 bp paired-end reads with a HiSeq 2500 platform (Illumina) serviced by Macrogen in Seoul, Korea.

The RNA-seq data were processed by TopHat2 and Bowtie2 (Langmead and Salzberg, 2012; Trapnell et al., 2009). After quality control of raw sequence reads, provided by Phred quality score and FastQC, adaptor sequences were trimmed using the Trimmomatic program to generate trimmed data (Bolger et al., 2014). The trimmed reads were aligned and assembled with the G. max v2.1 gene model (https://www.ncbi.nlm.nih.gov/assembly/GCA_000004515.4/) using the HISAT2 and StringTie (Kim et al., 2019; Pertea et al., 2015, 2016). The raw read counts, FPKM (fragments per kilobase of transcript per million mapping reads), TPM (transcripts per million), and simple description of the mapped transcripts to known soybean genes were presented in Supplementary Table 1.

Followed by the normalization by trimmed mean of M-values in edgeR, differentially expressed genes (DEGs) in infected plants were determined by two statistical outputs (|fold change (FC)| ≥ 2 and exactTest raw. P < 0.05) compared with the mock-inoculated plants (Chen et al., 2016). FC (infected samples/non-infected samples) and related statistical data of each DEG were listed in Supplementary Table 2.

Functional classification and enrichment analyses of DEGs

The DEGs were used for gene ontology (GO) enrichment analysis using Protein Analysis Through Evolutionary Relationships (PANTHER) (Mi et al., 2019; Thomas et al., 2003) and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway analysis (Kanehisa et al., 2007). The GO terms and KEGG pathways with P < 0.01 and Q < 0.01 were considered to be significantly enriched in DEGs. The list of enriched GO terms and KEGG pathways related to each DEG were shown in Supplementary Tables 3 and 4. We can share detailed analysis procedures under a request.

Results and Discussion

The soybean cv. Kwangan owns good agronomic features such as high resistance to soybean mosaic virus and high protein quality (Kim et al., 2011; Lee et al., 2015; Soybean Breeding Team, Upland Crop Div., Crop Experiment and Experiment Station, 1994). Therefore, the cv. Kwangan was broadly cultivated in Korea and used as a well-established parental plant to develop the Biotech soybean (Kim et al., 2012; Yeom et al., 2020). However, the molecular response of cv. Kwangan to bacterial disease, especially bacteria brown spot caused by P. syringae pv. syringae, has not been done. This study focused on the transcriptional changes in susceptible response in soybean plants upon PssB728a infection.

Growth of P. syringae pv. syringe B728a in trifoliate leaves of soybean plants

PssB728a exploits two different habitats (leaf surface and apoplast) as a phytopathogenic bacterium. The epiphytic PssB728a population acts as an inoculum to infect host plants successfully (Hirano and Upper, 2000). The epiphytic bacterial growth was counted every 2 days after a spray-inoculation of PssB728a. The number of bacteria on the leaf surface increased 2 days post inoculation (dpi), and the bacterial growth was saturated on the 4th day and 6th day (4 and 6 dpi) with populations of 106 colony forming unit (cfu) per leaf disc (Fig. 1A, left). Consistently, the endophytic bacterial population reached > 108 cfu per leaf disc at 4 dpi (Fig. 1A, right), and the lesion increased in size as time progressed (Fig. 1B). The bacterial growth assay reveals that PssB728a causes a compatible interaction with soybean cv. Kwangan.

Fig. 1

Disease symptom development and bacterial growth in leaves of soybean cv. Kwangan infected with Pseudomonas syringae pv. syrinage (Pss) B728a. (A) The epiphytic and endophytic bacterial growths were tested on infected soybean leaves. At least eight leaf discs were used to count the endophytic and epiphytic bacterial growth of PssB728a on the indicated day post inoculation (dpi). The bars present the mean ± standard error. (B) The second trifoliate leaves from 3-week-old soybean plants grown in long-day conditions were spray-inoculated with 10 mM MgSO4 (mock) and PssB728a (OD600 = 0.01) and incubated under dark conditions with high humidity and 25°C during the experiment. Disease symptoms were taken in photographs 0, 3, and 5 dpi.

DEGs during P. syringae pv. syringe B728a infection

Transcriptome sequencing libraries were made with total RNA from soybean leaves collected on day 5 after inoculation with 10 mM MgSO4 (non-infected sample, NI) and PssB728a (infected sample, IF). Among 60,506 known soybean transcripts, 31,608 genes, except 28,898 not detected in at least one sample, were used in the transcriptome analysis. The correlation matrix of Pearson’s coefficient among four samples (NI_1, NI_2, IF_1, and IF_2) indicated that two independent biological replicates were highly reproducible (Fig. 2A). In addition, multidimensional scaling analysis showed that NIs and IFs were separately clustered, and components I and II reflected a 97.4% and 2.1% difference among samples, respectively (data not shown). These statistical analyses show that the transcriptome-seq data can represent the transcriptional reprogramming occurring in infected soybean leaves with PssB728a at 5 dpi.

Fig. 2

Quality of RNA-seq analysis and number of differentially expressed genes (DEGs) in infected soybean leaves with Pseudomonas syringae pv. syringae B728a. (A) Correlation matrix of Pearson’s coefficient (–1 ≤ r ≤ 1) to determine the similarity of each sample. NI and IF mean samples from non-infected and PssB728a-infected leaves of soybean plants, respectively. A closer value to r = 1 presents a higher similarity between samples. (B) A volcano plot to verify DEGs in PssB728a-infected plants versus non-infected soybean. Blue and brown dots indicate genes displaying fold change (FC) ≤ −2 and FC ≥ 2 with P < 0.05, respectively. Vertical dotted lines showed the |FC| = 2 and a horizontal dotted line indicates P = 0.05. (C) Numbers of genes whose FC is up- and down-regulated on day 5 after PssB728a infection. (D) Heat map showing the expression pattern of DEGs in two biological replications of non-infected (NI) leaves and PssB728a-infected (IF) leaves. Genes are grouped according to hierarchical clustering analysis (Euclidean distance, complete linkage). Rows are clustered using distance and average linkage. The color gradient represents the Z-score for log2 FC-based normalized values.

In a volcano plot, the up-regulated and down-regulated genes with |FC| ≥ 2 and P < 0.05 were displayed on the right (brown dots) and the left (blue dots) sides, respectively (Fig. 2B). About 8,900 up-regulated genes (14.7% of known soybean transcripts) and 7,780 down-regulated genes (12.9%) were assigned as DEGs in infected leaves with PssB728a (Fig. 2C). Hierarchical clustering of expression levels of DEGs (up to 3,000 of 16,680 genes) indicated that two biological replicates for non-infected plants and infected plants had also been done reproducibly (Fig. 2D). Therefore, we suggest that DEGs’ information obtained from the transcriptome analysis helps us to understand the transcriptional changes at 5 dpi in susceptible soybean plants to PssB728a infection. Whole lists of DEGs were presented in Supplementary Table 2.

The enriched GO terms of DEGs

The cellular and physiological processes of plants are altered upon pathogen infection. Molecular networks composed of complex interactions between/among genes, proteins, and metabolic compounds are involved especially in plant immune signaling, transcription of immune-related genes, and biosynthesis of phytohormones and secondary metabolites, as the outcome of a plant-microbe interaction (Delplace et al., 2022). A compatible interaction between host and pathogen leads to successful infection and disease (Glazebrook, 2005). Different from an incompatible interaction accompanying a hypersensitive response due to intensified immune responses, the compatible interaction entails interruption of plant immunity by effectors and toxins derived from pathogens, resulting in suppression of immune responses during pathogenesis (Chisholm et al., 2006; Dangl and Jones, 2001; Dodds and Rathjen, 2010; Jones and Dangl, 2006). To analyze which biological processes and molecular functions of DEGs were implicated in the compatible interaction between cv. Kwangan and PssB728a, the GO enrichment analysis was accomplished, which explains the associations between gene products and GO terms using the PANTHER classification system.

First, all up-regulated and down-regulated genes were subjected to GO analysis (GO biological process) to point out biological processes in which DEGs are involved. The top 10 categories exhibiting the lowest P-values were listed with the fold enrichments determined by the ratio of the number of observed genes versus that of expected genes (Fig. 3A and B). The up-regulated DEGs were enriched in ‘small molecule metabolic process (GO:0044281)’, ‘child terms of organonitrogen compound metabolic process (GO:1901564) and cellular amide metabolic process (GO:0043603)’, and ‘generation of precursor metabolites and energy (GO:0006091)’ (Fig. 3A). On the other hand, certain kinds of genes involved in ‘regulation of metabolic process (GO:0019222)’, ‘macromolecule metabolic process (GO:0043170)’, and ‘protein modification process (GO:0006464)’ were significantly down-regulated by PssB728a at 5 dpi (Fig. 3B). In addition, large sets of genes associated with RNA biological process appeared to be influenced by pathogen infection, because the expression of those genes involved in RNA biosynthesis, RNA metabolic process, mRNA metabolic process, and regulation of RNA biosynthetic process decreased upon bacterial infection (Supplementary Table 3). In the same line with it, the number of observed genes related to transcription was significantly less in the infected plants compared with that of expected genes of the soybean genome (Fig. 3C). Supporting this, the expression of 99 genes which are implicated in the metabolism of purines (heterocyclic aromatic organic compound consisting of two rings fused, adenine and guanine here), essential building block components for DNA and RNA, was altered by PssB728a infection (Supplementary Table 4). Second, we also analyzed the molecular function of DEGs’ products (GO molecular function and GO-slim molecular function). Both up- and down-regulated DEGs were significantly enriched in 62 GO terms, including transferase activity, catalytic activity, kinase activity, and small molecule binding (Supplementary Table 3). These results suggest that primary and secondary metabolisms are mainly affected by PssB728a infection, whereas transcriptional activity decreased at 5 dpi.

Fig. 3

Enriched gene ontology (GO) terms of differentially expressed genes (DEGs) in infected soybean plants. (A, B) Top 10 GO biological process terms based on the number of observed genes in this analysis among up-regulated DEGs (A) and down-regulated DEGs (B). The numeric on the right side of the bars shows the exact numbers of observed genes. Fold enrichment was calculated by the ratio of observed gene versus expected gene. (C) Terms of GO-slim molecular function showing fold enrichment < 1, P < 0.005, and Q < 0.05 in up-regulated DEGs. The whole GO biological process, cellular component, and molecular function terms for each DEG were presented in Supplementary Table 3. FDR, false discovery rate.

Numerous genes encoding enzymes such as hydrolases, lyases, transferase, and metabolite interconversion enzymes were expressed highly upon bacterial infection (Fig. 4A). In contrast, a large number of genes encoding kinases, especially non-receptor serine/threonine protein kinases (about more than 150 genes) and protein modifying enzymes such as ubiquitin-protein ligases, expressed lower in infected soybean plants than mock-inoculated plants (Fig. 4B). In fact, KEGG pathway analysis revealed that the high population of DEGs was mostly associated with metabolic pathways (carbonate, secondary metabolite, starch and sucrose, amino sugar metabolism) and protein processing (Supplemental Table 4).

Fig. 4

Enriched gene ontology (GO) protein class terms of differentially expressed genes (DEGs) in infected soybean plants. (A, B) Top 10 GO protein class terms based on the number of observed genes in this analysis among up-regulated DEGs (A) and down-regulated DEGs (B). Up-regulated (8,902 genes, A) and down-regulated (7,786 genes, B) genes were applied to GO enrichment protein class analysis. For each protein, class GO term, the expected (empty) and observed gene (filled) numbers whose P-values are lower than 0.01 (statistical significance) are presented. No. of observed genes, number of DEGs associated with a GO term for protein class; no. of expected genes, number of genes expected for each GO term in the genome after normalization of the whole number of genes analyzed in this analysis (8,902 and 7,786).

The KEGG pathways related to DEGs

Transcriptional profiles of several crop plants after pathogen infections clarified particular metabolic and signaling pathways tightly related to compatible and incompatible interaction (Cregeen et al., 2015; Gyetvai et al., 2012; Meng et al., 2021; Sonah et al., 2016; van Esse et al., 2009; Wang et al., 2010).

Using 3,775 genes mapped in KEGG pathways (http://www.kegg.jp/kegg/pathway.html) among 16,688 DEGs, gene-set enrichments were analyzed to examine which pathways were activated or suppressed at 5 dpi in infected leaves (Kanehisa et al., 2022). Fig. 5 presents the top 15 pathways with P-values < 7.47E-45 and Q-values < 6.63E-44 as follows; metabolic pathways, biosynthesis of secondary metabolite, carbon metabolism, plant hormone signal transduction, MAPK signaling pathway, plant-pathogen interaction, photosynthesis, protein processing in the endoplasmic reticulum (ER), purine metabolism, pyruvate metabolism, etc. Many genes were mapped simultaneously into different pathways. For example, among 137 genes in the category of glycolysis [7], whose expressions were affected upon bacterial infection, 49 genes (about 36%) participate in the pyruvate metabolism [15]; 50% of genes (49 genes) in pyruvate metabolism work at glycolysis-related pathway as well.

Fig. 5

The top 15 Kyoto Encyclopedia of Genes and Genomes (KEGG) pathways in which differentially expressed genes (DEGs) are participating. Each protein ID of DEGs was converted to KEGG gene ID using a biological DataBase network (bioDBnet, https://biodbnet-abcc.ncifcrf.gov/) and the converted KEGG gene IDs were used for the KEGG pathway analysis (KEGG Mapper-Color). The top 15 pathways were presented here with the exact number of genes. The whole KEGG pathways related to each DEG were listed in Supplementary Table 4.

Cellular detoxification of reactive oxygen species (ROS) formed from different metabolic pathways and upon various unfavorable stress conditions is an essential biological process in living organisms to maintain homeostasis (Huang et al., 2019). Upon bacterial infection, expression of the biosynthetic genes of antioxidant compounds appeared to be increased as expected. For instance, the expression of genes encoding enzymes involved in ascorbate biosynthesis, such as mannose-6-phosphate isomerase 1 (EC: 5.3.1.8, Entrez_Gene_ID: 100815644), phosphomannomutase (EC: 5.4.2.8, Entrez_Gene_ID: 732605), and mannose-1-phosphate guanylyltransferase 1-like (EC: 2.7.7.13, Entrez_Gene_ID: 100777120), increased during infection (Supplementary Table 4). In addition, numerous glutathione S-transferases, which function as one of the major detoxifying enzymes by conjugating tripeptide (γ-Glu-Cys-Gly) glutathione to various hydrophobic and electrophilic substrates (Liu et al., 2015a), were actively transcribed in the leaves infected with PssB728a (Supplementary Table 4).

Furthermore, a variety of genes related to an ER stress response, which occurs when the capacity of the ER to fold proteins becomes saturated and leads to activation of stress-adaptive response called unfolded protein response (UPR) (Bao and Howell, 2017; Ruberti et al., 2015), showed the altered transcription patterns under PssB728a infection. Three soybean orthologs of an ER stress sensor, the ER-resident transmembrane protein inositol-requiring enzyme 1 (serine/threonine-protein kinase/endoribonuclease IRE1, EC: 2.7.11.1 or EC: 3.1.26, Entrez_Gene_IDs: 100785904, 100788389, 102661527), reduced their mRNA expression upon bacterial infection (Supplemental Table 4). An ortholog of another ER-resident transmembrane protein ER stress sensor, PERK (eukaryotic translation initiation factor 2-alpha kinase 4/eIF-2-alpha kinase GCN2 isoform X1 (EC:2.7.11.1, Entrez_Gene_ID: 100797996), which has been characterized as an ER stress sensor in animal cells but not identified experimentally in plant cells (Park and Park, 2019), also showed the reduced expression in infected leaves. In contrast, several ER-stress responsive genes and genes whose products participate in the protein secretion process exhibited higher expression upon bacterial infection (Supplemental Table 4). An ER stress-inducing compound, tunicamycin enhanced the susceptibility to P. syringae pv. tomato DC3000 in Arabidopsis while gene expression of PR1 and FLS2-mediated PTI such as callose deposition increased (Chakraborty et al., 2017). ire1 and bzip60 mutants, which are deficient in IRE1, an ER stress sensor activating UPR and its splicing target bZIP-transcription factor, bZIP60, respectively, were more susceptible to P. syringae pv. maculicola (Moreno et al., 2012). Consistently, knock-down tobacco plants of IRE1 and bZIP60 via virus-induced gene silencing were more sensitive to the necrotic fungal pathogen, Alternaria alternata (Xu et al., 2019). The transcriptome profile of soybean plants infected with PssB728a suggests that ER-mediated biological process is significantly affected in a compatible interaction (Saijo et al., 2009).

Among the top 15 pathways (Fig. 5), we analyzed further details about plant hormone signal transduction [5], MAPK signaling pathway-plant [10], and plant-pathogen interaction [11]. Genes participating in auxin- and cytokinin-mediated signaling were up-regulated in soybean leaves under PssB728a infection, whereas ethylene (ET)-, jasmonic acid (JA)-, and salicylic acid (SA)-mediated signaling seemed to be down-regulated in PssB728a-infected leaves (Supplementary Fig. 1). In the hormone-biosynthetic aspects, transcript levels of genes participating in the last half of α-linolenic acid metabolic pathway for JA biosynthesis decreased by more than 2-folds in the infected plant (Supplementary Fig. 2). On the other side, soybean orthologs encoding three critical enzymes for the ET biosynthesis, S-ADENOSYLMETHIONINE SYNTHETASE (SAM), 1-AMINOCYCLOPROPANE-1-CARBOXYLATE SYNTHASE (ACS), and AMINOCYCLOPROPANECARBOXYLATE OXIDASE (ACO) were actively transcribed under PssB728a infection (Supplementary Table 4). Additionally, a putative ortholog of Arabidopsis ISOCHORISMATE SYNTHASE1 (ICS1), ICS was strongly transcribed in the infected leaves (EC:5.4.4.2, Entrez_Gene_ID: 100799025) albeit genes involved in shikimate pathway did not show any trends on transcriptional change under infection (Supplementary Table 4). Note that we do not have any experimental clues on the dynamics of these hormones in the infected plants. To measure hormone levels during infection will be necessary to understand the compatible interaction. Taken together, we assume that PssB728a successfully suppresses plant immune signaling activated by ET, JA, and SA and employs growth-related hormone signaling for their colonization.

MAPK signal cascades consisting of three sequential serine/threonine kinases (MAPKKK → MAPKK → MAPK) transmit environmental and developmental signals upon a range of stresses, including pathogen infection, heat, freezing, ROS, drought, and wounding (Meng and Zhang, 2013). Various hormones such as auxin, abscisic acid (ABA), JA, SA, and ET affect the MAPK signaling pathway (Jagodzik et al., 2018). Upon PssB728a infection, ABA receptors PYRABACTIN RESISTANCE 1 (PYR1)/PYRABACTIN RESISTANCE 1-LIKE (PYL) and its downstream signaling components type 2C protein phosphatases (PP2C), negative regulators of ABA signaling, displayed higher mRNA expression. In contrast, the expressions of MAPK signaling components were reduced, probably leading to a non-efficient and unsuccessful stress-adaptive response with disease increase (Supplementary Fig. 3). MAPK pathway relays signals from membrane-residing immune receptors to downstream components for defense activation in plants (Thulasi Devendrakumar et al., 2018). Most of the MAPKKKs, MAPKKs, and MAPKs detected in this study were down-regulated in the transcript levels on day 5 after PssB728a infection (Supplementary Fig. 3). Besides MAPK signaling cascades, genes encoding extracellular (co)receptors, CYCLIC NUCLEOTIDE GATED CHANNELs, and Ca2+-signaling proteins were weakly expressed in the infected leaves of soybean (Supplementary Fig. 4). However, a soybean ortholog of FLAGELLIN-SENSITIVE 2 (FLS2) showed increased mRNA expression in the soybean leaves infected with PssB728a. However, its downstream pathways, MEKK1-MKK4/5-MPK3/6 and MPK4, reduced their transcription. In addition, soybean homologous genes of well-identified defense regulators in Arabidopsis and tomato, such as EDS1, RIN4, SGT1, and Pti1/5/6, displayed reduced mRNA levels in the infected leaves (Supplementary Fig. 4). This KEGG pathway analysis suggests that the plant immune responses of soybean cv. Kwangan are compromised by PssB728a.

In this study, we developed a pathosystem between soybean cv. Kwangan and PssB278a, a causal pathogen of bacterial brown spot disease, and presented the transcriptional changes in soybean leaves on day 5 after PssB728a infection. Unexpectedly, mRNA expression of homologous genes to well-identified immune-related genes in model and other crop plants decreased or was not dramatically altered in soybean leaves on day 5 after infection. In contrast, primary and secondary metabolisms underwent a change in PssB728a-infected conditions. Thus, these analyses represent the transcriptional changes at the late phase of a compatible interaction and reflect the suppression of plant immunity in susceptible plants by virulent pathogen infection. Furthermore, this study identifies large number of candidates for future functional analyses and provides a foundation to understand molecular mechanisms implicated in the compatible interaction of a susceptible Korean soybean cultivar and bacterial brown spot disease.

Acknowledgments

This work was supported by the New Breeding Technologies Development Program of the Rural Development Administration (PJ01477702 and PJ01653304) (Ho Won Jung), the Basic Science Research Program through the National Research Foundation of Korea (NRF) funded by the Ministry of Education, Republic of Korea (2020R1A6A1A03047729) (Ho Won Jung and Hee Jin Park), and the Green Fusion Technology Graduate School Program of the Ministry of Environment (Myoungsub Kim, Dohui Lee).

Notes

Conflicts of Interest

No potential conflict of interest relevant to this article was reported.

Electronic Supplementary Material

References

Bao Y., Howell S.H.. 2017;The unfolded protein response supports plant development and defense as well as responses to abiotic stress. Front. Plant Sci 8:344.
Bolger A.M., Lohse M., Usadel B.. 2014;Trimmomatic: a flexible trimmer for Illumina sequence data. Bioinformatics 30:2114–2120.
Carter T.E., Nelson R.L., Sneller C.H., Cui Z.. 2004. Genetic diversity in soybean. Soybeans: improvement, production, and uses 3rd edth ed. In : Shibles R.M., Harper J.E., Wilson R.F., Shoemaker R.C., eds. p. 303–416. American Society of Agronomy. Madison, WI, USA:
Chakraborty R., Macoy D.M., Lee S.Y., Kim W.-Y., Kim M.G.. 2017;Tunicamycin-induced endoplasmic reticulum stress suppresses plant immunity. Appl. Biol. Chem 60:623–630.
Chen Y., Lun A.T.L., Smyth G.K.. 2016;From reads to genes to pathways: differential expression analysis of RNA-Seq experiments using Rsubread and the edgeR quasi-likelihood pipeline. F1000Res 5:1438.
Chisholm S.T., Coaker G., Day B., Staskawicz B.J.. 2006;Host-microbe interactions: shaping the evolution of the plant immune response. Cell 124:803–814.
Cooper B., Campbell K.B., McMahon M.B., Luster D.G.. 2013;Disruption of Rpp1-mediated soybean rust immunity by virus-induced gene silencing. Plant Signal. Behav 8:e27543.
Cregeen S., Radisek S., Mandelc S., Turk B., Stajner N., Jakse J., Javornik B.. 2015;Different gene expressions of resistant and susceptible Hop cultivars in response to infection with a highly aggressive strain of Verticillium albo-atrum . Plant Mol. Biol. Rep 33:689–704.
Dangl J.L., Jones J.D.G.. 2001;Plant pathogens and integrated defence responses to infection. Nature 411:826–833.
Delgado-Cerrone L., Alvarez A., Mena E., Ponce de León I., Montesano M.. 2018;Genome-wide analysis of the soybean CRK-family and transcriptional regulation by biotic stress signals triggering plant immunity. PLoS ONE 13:e0207438.
Delplace F., Huard-Chauveau C., Berthom R., Roby D.. 2022;Network organization of the plant immune system: from pathogen perception to robust defense induction. Plant J 109:447–470.
Dodds P.N., Rathjen J.P.. 2010;Plant immunity: towards an integrated view of plant–pathogen interactions. Nat. Rev. Genet 11:539–548.
Ercolani G.L., Hagedorn D.J., Kelman A., Rand R.E.. 1974. Epiphytic survival of Pseudomonas syringae on hairy vetch in relation to epidemiology of bacterial brown spot of bean in Wisconsin. Phytopathology 64p. 1330–1339.
Faske T., Kirkpatrick T., Zhou J., Tzanetakis I.. 2014. Soybean diseases. Arkansas soybean production handbook - MP197 p. 1–18. The Soybean Commodity Committee of the Cooperative Extension Service, University of Arkansas. Fayetteville, AR, USA:
Feil H., Feil W.S., Chain P., Larimer F., DiBartolo G., Copeland A., Lykidis A., Trong S., Nolan M., Goltsman E., Thiel J., Malfatti S., Loper J.E., Lapidus A., Detter J.C., Land M., Richardson P.M., Kyrpides N.C., Ivanova N., Lindow S.E.. 2005;Comparison of the complete genome sequences of Pseudomonas syringae pv. syringae B728a and pv. tomato DC3000. Proc. Natl. Acad. Sci. U. S. A 102:11064–11069.
Glazebrook J.. 2005;Contrasting mechanisms of defense against biotrophic and necrotrophic pathogens. Annu. Rev. Phytopathol 43:205–227.
Gnanamanickam S.S., Ward E.W.B.. 1982;Characterization of Pseudomonas syringae strains causing disease symptoms on soybean. Can. J. Plant Pathol 4:233–236.
Gyetvai G., Sønderkær M., Göbel U., Basekow R., Ballvora A., Imhoff M., Kersten B., Nielsen K.L., Gebhardt C.. 2012;The transcriptome of compatible and incompatible interactions of potato (Solanum tuberosum) with Phytophthora infestans revealed by DeepSAGE analysis. PLoS ONE 7:e31526.
Hartman G.L., Rupe J.C., Sikora E.J., Domier L.L., Davis J.A., Steffey K.L.. 2015. Compendium of soybean diseases and pests 5th edth ed. American Phytopathological Society. St. Paul, MN, USA: p. 201.
Helm M., Qi M., Sarkar S., Yu H., Whitham S.A., Innes R.W.. 2019;Engineering a decoy substrate in soybean to enable recognition of the soybean mosaic virus NIa protease. Mol. Plant-Microbe Interact 32:760–769.
Hirano S.S., Baker L.S., Upper C.D.. 1996;Raindrop momentum triggers growth of leaf-associated populations of Pseudomonas syringae on field-grown snap bean plants. Appl. Environ. Microbiol 62:2560–2566.
Hirano S.S., Clayton M.K., Upper C.D.. 1994;Estimation of and temporal changes in means and variances of populations of Pseudomonas syringae on snap bean leaflets. Phytopathology 84:934–940.
Hirano S.S., Upper C.D.. 1990;Population biology and epidemiology of Pseudomonas syringae . Annu. Rev. Phytopathol 28:155–177.
Hirano S.S., Upper C.D.. 2000;Bacteria in the leaf ecosystem with emphasis on Pseudomonas syringae-a pathogen, ice nucleus, and epiphyte. Microbiol. Mol. Biol. Rev 64:624–653.
Huang H., Ullah F., Zhou D.-X., Yi M., Zhao Y.. 2019;Mechanisms of ROS regulation of plant development and stress responses. Front. Plant Sci 10:800.
Jagodzik P., Tajdel-Zielinska M., Ciesla A., Marczak M., Ludwikow A.. 2018;Mitogen-activated protein kinase cascades in plant hormone signaling. Front. Plant Sci 9:1387.
Jones J.D.G., Dangl J.L.. 2006;The plant immune system. Nature 444:323–329.
Kanehisa M., Araki M., Goto S., Hattori M., Hirakawa M., Itoh M., Katayama T., Kawashima S., Okuda S., Tokimatsu T., Yamanishi Y.. 2007;KEGG for linking genomes to life and the environment. Nucleic Acids Res 36:D480–D484.
Kanehisa M., Sato Y., Kawashima M.. 2022;KEGG mapping tools for uncovering hidden features in biological data. Protein Sci 31:47–53.
Kennelly M.M., Cazorla F.M., de Vicente A., Ramos C., Sundin G.W.. 2007; Pseudomonas syringae diseases of fruit trees: progress toward understanding and control. Plant Dis 91:4–17.
Kim D., Paggi J.M., Park C., Bennett C., Salzberg S.L.. 2019;Graph-based genome alignment and genotyping with HISAT2 and HISAT-genotype. Nat. Biotechnol 37:907–915.
Kim M.-J., Kim J.K., Kim H.J., Pak J.H., Lee J.-H., Kim D.-H., Choi H.K., Jung H.W., Lee J.-D., Chung Y.-S., Ha S.-H.. 2012;Genetic modification of the soybean to enhance the β-carotene content through seed-specific expression. PLoS One 7:e48287.
Kim Y.-J., Lee K.-W., Cho S.-K., Oh Y.-J., Shin S.-O., Paik C.-H., Kim K.-H., Kim T.-S., Kim K.-J.. 2011;Selection and quality evaluation of sprout soybean [Glycine max (L.) Merrill] variety for environment-friendly cultivation in southern paddy field. Korean J. Org. Agric 19:357–372. (in Korean).
Langmead B., Salzberg S.L.. 2012;Fast gapped-read alignment with Bowtie 2. Nat. Methods 9:357–359.
Lee C., Choi M.-S., Kim H.-T., Yun H.-T., Lee B., Chung Y.-S., Kim R.W., Choi H.-K.. 2015;Soybean [Glycine max (L.) Merrill]: importance as a crop and pedigree reconstruction of Korean varieties. Plant Breed. Biotechnol 3:179–196.
Lee J., Teitzel G.M., Munkvold K., del Pozo O., Martin G.B., Michelmore R.W., Greenberg J.T.. 2012;Type III secretion and effectors shape the survival and growth pattern of Pseudomonas syringae on leaf surfaces. Plant Physiol 158:1803–1818.
Li M.-W., Wang Z., Jiang B., Kaga A., Wong F.-L., Zhang G., Han T., Chung G., Nguyen H., Lam H.-M.. 2020;Impacts of genomic research on soybean improvement in East Asia. Theor. Appl. Genet 133:1655–1678.
Lindemann J., Arny D.C., Upper C.D.. 1984;Use of an apparent infection threshold population of Pseudomonas syringae to predict incidence and severity of brown spot of bean. Phytopathology 74:1334–1339.
Liu H.-J., Tang Z.-X., Han X.-M., Yang Z.-L., Zhang F.-M., Yang H.-L., Liu Y.-J., Zeng Q.-Y.. 2015a;Divergence in enzymatic activities in the soybean GST supergene family provides new insight into the evolutionary dynamics of whole-genome duplicates. Mol. Biol. Evol 32:2844–2859.
Liu J.-Z., Graham M.A., Pedley K.F., Whitham S.A.. 2015b;Gaining insight into soybean defense responses using functional genomics approaches. Brief. Funct. Genomics 14:283–290.
Liu Y., Du H., Li P., Shen Y., Peng H., Liu S., Zhou G.-A., Zhang H., Liu Z., Shi M., Huang X., Li Y., Zhang M., Wang Z., Zhu B., Han B., Liang C., Tian Z.. 2020;Pan-genome of wild and cultivated soybeans. Cell 182:162–176.e13.
Loper J.E., Lindow S.E.. 1987;Lack of evidence for the in situ fluorescent pigment production by Pseudomonas syringae pv. syringae on bean leaf surfaces. Phytopathology 77:1449–1454.
Martin J.A., Wang Z.. 2011;Next-generation transcriptome assembly. Nat. Rev. Genet 12:671–682.
Meng H., Sun M., Jiang Z., Liu Y., Sun Y., Liu D., Jiang C., Ren M., Yuan G., Yu W., Feng Q., Yang A., Cheng L., Wang Y.. 2021;Comparative transcriptome analysis reveals resistant and susceptible genes in tobacco cultivars in response to infection by Phytophthora nicotianae . Sci Rep 11:809.
Meng X., Zhang S.. 2013;MAPK cascades in plant disease resistance signaling. Annu. Rev. Phytopathol 51:245–266.
Mi H., Muruganujan A., Huang X., Ebert D., Mills C., Guo X., Thomas P.D.. 2019;Protocol update for large-scale genome and gene function analysis with the PANTHER classification system (v.14.0). Nat. Protoc 14:703–721.
Moreno A.A., Mukhtar M.S., Blanco F., Boatwright J.L., Moreno I., Jordan M.R., Chen Y., Brandizzi F., Dong X., Orellana A., Pajerowska-Mukhtar K.M., Polymenis M.. 2012;IRE1/bZIP60-mediated unfolded protein response plays distinct roles in plant immunity and abiotic stress responses. PLoS ONE 7:e31944.
Park C.-J., Park J.M.. 2019;Endoplasmic reticulum plays a critical role in integrating signals generated by both biotic and abiotic stress in plants. Front. Plant Sci 10:399.
Pertea M., Kim D., Pertea G.M., Leek J.T., Salzberg S.L.. 2016;Transcript-level expression analysis of RNA-seq experiments with HISAT, StringTie and Ballgown. Nat. Protoc 11:1650–1667.
Pertea M., Pertea G.M., Antonescu C.M., Chang T.-C., Mendell J.T., Salzberg S.L.. 2015;StringTie enables improved reconstruction of a transcriptome from RNA-seq reads. Nat. Biotechnol 33:290–295.
Pottinger S.E., Bak A., Margets A., Helm M., Tang L., Casteel C., Innes R.W.. 2020;Optimizing the PBS1 decoy system to confer resistance to Potyvirus infection in Arabidopsis and soybean. Mol. Plant Microbe-Interact 33:932–944.
Ruberti C., Kim S.-J., Stefano G., Brandizzi F.. 2015;Unfolded protein response in plants: one master, many questions. Curr. Opin. Plant Biol 27:59–66.
Russell A.R., Ashfield T., Innes R.W.. 2015; Pseudomonas syringae Effector AvrPphB suppresses AvrB-induced activation of RPM1 but not AvrRpm1-induced activation. Mol. Plant Microbe-Interact 28:727–735.
Saijo Y., Tintor N., Lu X., Rauf P., Pajerowska-Mukhtar K., Häweker H., Dong X., Robatzek S., Schulze-Lefert P.. 2009;Receptor quality control in the endoplasmic reticulum for plant innate immunity. EMBO J 28:3439–3449.
Savary S., Willocquet L., Pethybridge S.J., Esker P., McRoberts N., Nelson A.. 2019;The global burden of pathogens and pests on major food crops. Nat. Ecol. Evol 3:430–439.
Schmutz J., Cannon S.B., Schlueter J., Ma J., Mitros T., Nelson W., Hyten D.L., Song Q., Thelen J.J., Cheng J., Xu D., Hellsten U., May G.D., Yu Y., Sakurai T., Umezawa T., Bhattacharyya M.K., Sandhu D., Valliyodan B., Lindquist E., Peto M., Grant D., Shu S., Goodstein D., Barry K., Futrell-Griggs M., Abernathy B., Du J., Tian Z., Zhu L., Gill N., Joshi T., Libault M., Sethuraman A., Zhang X.-C., Shinozaki K., Nguyen H.T., Wing R.A., Cregan P., Specht J., Grimwood J., Rokhsar D., Stacey G., Shoemaker R.C., Jackson S.A.. 2010;Genome sequence of the palaeopolyploid soybean. Nature 463:178–183.
Sonah H., Zhang X., Deshmukh R.K., Borhan M.H., Fernando W.G.D., Bélanger R.R.. 2016;Comparative transcriptomic analysis of virulence factors in Leptosphaeria maculans during compatible and incompatible interactions with canola. Front. Plant Sci 7:1784.
Soybean Breeding Team, Upland Crop Div., Crop Experiment and Experiment Station. 1994;A new high seed-protein, small grain and high-yielding soybean variety “Kwangankong. Korean J. Breed. Sci 16:462.
Thomas P.D., Kejariwal A., Campbell M.J., Mi H., Diemer K., Guo N., Ladunga I., Ulitsky-Lazareva B., Muruganujan A., Rabkin S., Vandergriff J.A., Doremieux O.. 2003;PANTHER: a browsable database of gene products organized by biological function, using curated protein family and subfamily classification. Nucleic Acids Res 31:334–341.
Thulasi Devendrakumar K., Li X., Zhang Y.. 2018;MAP kinase signalling: interplays between plant PAMP- and effector-triggered immunity. Cell. Mol. Life Sci 75:2981–2989.
Trapnell C., Pachter L., Salzberg S.L.. 2009;TopHat: discovering splice junctions with RNA-Seq. Bioinformatics 25:1105–1111.
van Esse H.P., Fradin E.F., de Groot P.J., de Wit P.J.G.M., Thomma B.P.H.J.. 2009;Tomato transcriptional responses to a foliar and a vascular fungal pathogen are distinct. Mol. Plant Microbe-Interact 22:245–258.
Vinatzer B.A., Teitzel G.M., Lee M.-W., Jelenska J., Hotton S., Fairfax K., Jenrette J., Greenberg J.T.. 2006;The type III effector repertoire of Pseudomonas syringae pv. syringae B728a and its role in survival and disease on host and non-host plants. Mol. Microbiol 62:26–44.
Wang X., Liu W., Chen X., Tang C., Dong Y., Ma J., Huang X., Wei G., Han Q., Huang L., Kang Z.. 2010;Differential gene expression in incompatible interaction between wheat and stripe rust fungus revealed by cDNA-AFLP and comparison to compatible interaction. BMC Plant Biol 10:9.
Wang Z., Tian Z.. 2015;Genomics progress will facilitate molecular breeding in soybean. Sci. China Life Sci 58:813–815.
Wei Y., Balaceanu A., Rufian J.S., Segonzac C., Zhao A., Morcillo R.J.L., Macho A.P.. 2020;An immune receptor complex evolved in soybean to perceive a polymorphic bacterial flagellin. Nat. Commun 11:3763.
Xu Z., Song N., Ma L., Wu J.. 2019;IRE1-bZIP60 pathway is required for Nicotiana attenuata resistance to fungal pathogen Alternaria alternata . Front. Plant Sci 10:263.
Yeom W.W., Kim H.J., Lee K.-R., Cho H.S., Kim J.-Y., Jung H.W., Oh S.-W., Jun S.E., Kim H.U., Chung Y.-S.. 2020;Increased production of α-Linolenic acid in soybean seeds by overexpression of LesquerellaFAD3-1 . Front. Plant Sci 10:1812.
Yuan Y., Yang Y., Yin J., Shen Y., Li B., Wang L., Zhi H.. 2020;Transcriptome-based discovery of genes and networks related to R SC3Q -mediated resistance to Soybean mosaic virus in soybean. Crop Pasture Sci 71:987.

Article information Continued

Fig. 1

Disease symptom development and bacterial growth in leaves of soybean cv. Kwangan infected with Pseudomonas syringae pv. syrinage (Pss) B728a. (A) The epiphytic and endophytic bacterial growths were tested on infected soybean leaves. At least eight leaf discs were used to count the endophytic and epiphytic bacterial growth of PssB728a on the indicated day post inoculation (dpi). The bars present the mean ± standard error. (B) The second trifoliate leaves from 3-week-old soybean plants grown in long-day conditions were spray-inoculated with 10 mM MgSO4 (mock) and PssB728a (OD600 = 0.01) and incubated under dark conditions with high humidity and 25°C during the experiment. Disease symptoms were taken in photographs 0, 3, and 5 dpi.

Fig. 2

Quality of RNA-seq analysis and number of differentially expressed genes (DEGs) in infected soybean leaves with Pseudomonas syringae pv. syringae B728a. (A) Correlation matrix of Pearson’s coefficient (–1 ≤ r ≤ 1) to determine the similarity of each sample. NI and IF mean samples from non-infected and PssB728a-infected leaves of soybean plants, respectively. A closer value to r = 1 presents a higher similarity between samples. (B) A volcano plot to verify DEGs in PssB728a-infected plants versus non-infected soybean. Blue and brown dots indicate genes displaying fold change (FC) ≤ −2 and FC ≥ 2 with P < 0.05, respectively. Vertical dotted lines showed the |FC| = 2 and a horizontal dotted line indicates P = 0.05. (C) Numbers of genes whose FC is up- and down-regulated on day 5 after PssB728a infection. (D) Heat map showing the expression pattern of DEGs in two biological replications of non-infected (NI) leaves and PssB728a-infected (IF) leaves. Genes are grouped according to hierarchical clustering analysis (Euclidean distance, complete linkage). Rows are clustered using distance and average linkage. The color gradient represents the Z-score for log2 FC-based normalized values.

Fig. 3

Enriched gene ontology (GO) terms of differentially expressed genes (DEGs) in infected soybean plants. (A, B) Top 10 GO biological process terms based on the number of observed genes in this analysis among up-regulated DEGs (A) and down-regulated DEGs (B). The numeric on the right side of the bars shows the exact numbers of observed genes. Fold enrichment was calculated by the ratio of observed gene versus expected gene. (C) Terms of GO-slim molecular function showing fold enrichment < 1, P < 0.005, and Q < 0.05 in up-regulated DEGs. The whole GO biological process, cellular component, and molecular function terms for each DEG were presented in Supplementary Table 3. FDR, false discovery rate.

Fig. 4

Enriched gene ontology (GO) protein class terms of differentially expressed genes (DEGs) in infected soybean plants. (A, B) Top 10 GO protein class terms based on the number of observed genes in this analysis among up-regulated DEGs (A) and down-regulated DEGs (B). Up-regulated (8,902 genes, A) and down-regulated (7,786 genes, B) genes were applied to GO enrichment protein class analysis. For each protein, class GO term, the expected (empty) and observed gene (filled) numbers whose P-values are lower than 0.01 (statistical significance) are presented. No. of observed genes, number of DEGs associated with a GO term for protein class; no. of expected genes, number of genes expected for each GO term in the genome after normalization of the whole number of genes analyzed in this analysis (8,902 and 7,786).

Fig. 5

The top 15 Kyoto Encyclopedia of Genes and Genomes (KEGG) pathways in which differentially expressed genes (DEGs) are participating. Each protein ID of DEGs was converted to KEGG gene ID using a biological DataBase network (bioDBnet, https://biodbnet-abcc.ncifcrf.gov/) and the converted KEGG gene IDs were used for the KEGG pathway analysis (KEGG Mapper-Color). The top 15 pathways were presented here with the exact number of genes. The whole KEGG pathways related to each DEG were listed in Supplementary Table 4.