TRANSCRIPTOMIC ANALYSIS OF SYM28 AND SYM29 SUPERNODULATING MUTANTS OF PEA (PISUM SATIVUM L.) UNDER COMPLEX INOCULATION WITH BENEFICIAL MICROORGANISMS Текст научной статьи по специальности «Биологические науки»

FULL COMMUNICATIONS

SYMBIOGENETICS

Citation: Zhukov, V., Zorin, E., Zhernakov, A., Afonin, A., Akhtemova, G., Bovin, A., Dolgikh, A., Gorshkov, A., Gribchenko, E., Ivanova, K., Kirienko, A., Kitaeva, A., Kliukova, M., Kulaeva, O., Kusakin, P., Leppyanen, I., Pavlova, O., Romanyuk, D., Rudaya, E., Serova, T., Shtark, O., Sulima, A., Tsyganova, A., Vasileva, E., Dolgikh, E., Tsyganov, V., and Tikhonovich, I. 2021. Transcriptomic analysis of sym28 and sym29 supernodulating mutants of pea (Pisum sativum L.) under complex inoculation with beneficial microorganisms. Bio. Comm. 66(3): 181-197. https://doi.org/10.21638/spbu03.2021.301 Authors' information: Vladimir Zhukov, PhD, Head of Laboratory, orcid.org/0000-0002-2411-9191; Evgeny Zorin, PhD Student, Research Engineer, orcid.org/0000-0001-5666-3020; Aleksandr Zhernakov, Junior Researcher, orcid. org/0000-0001-8961-9317; Alexey Afonin, PhD Student, Research Engineer, orcid.org/0000-

0002-8530-0226; Gulnar Akhtemova, PhD, Senior Researcher, orcid.org/0000-0001-7957-3693; Andrej Bovin, PhD Student, Researcher, orcid. org/0000-0003-4061-435X; Aleksandra Dolgikh, Research Engineer, orcid.org/0000-0003-1845-9701; Artemii Gorshkov, PhD Student, Research Engineer, orcid.org/0000-0002-7296-0188; Emma Gribchenko, Technician, orcid.org/0000-0002-1538-5527; Kira Ivanova, Junior Researcher, orcid. org/0000-0002-9119-065X; Anna Kirienko, PhD, Researcher, orcid.org/0000-0002-2519-306X; Anna Kitaeva, PhD, Researcher, orcid.org/0000-0001-7873-6737; Marina Kliukova, Research Engineer, orcid.org/0000-0003-1119-5512; Olga Kulaeva, PhD, Senior Researcher, orcid.org/0000-

0003-2687-9693; Pyotr Kusakin, PhD Student, Research Engineer, orcid.org/0000-0003-2795-515X; Irina Leppyanen, PhD, Senior Researcher, orcid.org/0000-0002-2158-0855; Olga Pavlova, PhD, Researcher, orcid.org/0000-0003-0528-5618; Daria Romanyuk, PhD, Researcher, orcid. org/0000-0001-9576-1256; Elizaveta Rudaya, PhD Student, Researcher, orcid.org/0000-0002-3081-9880; Tatiana Serova, PhD, Researcher, orcid. org/0000-0002-7236-1003; Oksana Shtark, PhD, Leading Researcher, orcid.org/0000-0002-3656-4559; Anton Sulima, PhD, Researcher, orcid. org/0000-0002-2300-857X; Anna Tsyganova, PhD, Leading Researcher, orcid.org/0000-0003-3505-4298; Ekaterina Vasileva, Research Engineer, orcid.org/0000-0001-5599-0361; Elena Dolgikh, Dr. of Sci. in Biology, Head of Laboratory, orcid. org/0000-0002-5375-0943; Viktor Tsyganov,

Dr. of Sci. in Biology, Head of Laboratory, orcid. org/0000-0003-3105-8689; Igor Tikhonovich, RAS Academician, Dr. of Sci. in Biology, Professor, Director for Science of ARRIAM, Dean of Faculty, orcid.org/0000-0001-8968-854X Manuscript Editor: Pavel Skutschas, Department of Vertebrate Zoology, Faculty of Biology, Saint Petersburg State University, Saint Petersburg, Russia Received: November 26, 2020; Revised: April 18, 2021; Accepted: May 11, 2021. Copyright: © 2021 Zhukov et al. This is an open-access article distributed under the terms of the License Agreement with Saint Petersburg State University, which permits to the authors unrestricted distribution, and self-archiving free of charge. Funding: This study was supported by the Russian Science Foundation grant No 17-76-30016. Ethics statement: This paper does not contain any studies involving human participants or animals performed by any of the authors. Competing interests: The authors have declared that no competing interests exist.

Transcriptomic analysis of sym28 and sym29 supernodulating mutants of pea (Pisum sativum L.) under complex inoculation with beneficial microorganisms

Vladimir Zhukov1, Evgeny Zorin1, Aleksandr Zhernakov1, Alexey Afonin1, Gulnar Akhtemova1, Andrej Bovin1, Aleksandra Dolgikh1, Artemii Gorshkov1, Emma Gribchenko1, Kira Ivanova1, Anna Kirienko1, Anna Kitaeva1, Marina Kliukova1, Olga Kulaeva1, Pyotr Kusakin1, Irina Leppyanen1, Olga Pavlova1, Daria Romanyuk1, Elizaveta Rudaya1, Tatiana Serova1, Oksana Shtark1, Anton Sulima1, Anna Tsyganova1, Ekaterina Vasileva1, Elena Dolgikh1, Viktor Tsyganov1, and Igor Tikhonovich12

''Department of Biotechnology, All-Russia Research Institute for Agricultural Microbiology, Shosse Podbel'skogo, 3, Saint Petersburg, 196608, Russian Federation 2Faculty of Biology, Saint Petersburg State University, Universitetskaya nab., 7-9, Saint Petersburg, 199034, Russian Federation

Address correspondence and requests for materials to Vladimir Zhukov, vzhukov@arriam.ru

Abstract

The garden pea (Pisum sativum L.), like most members of Fabaceae family, is capable of forming symbioses with beneficial soil microorganisms such as nodule bacteria (rhizobia), arbuscular mycorrhizal (AM) fungi and plant growth promoting bacteria (PGPB). The autoregulation of nodulation (AON) system is known to play an important role in controlling both the number of nodules and the level of root colonization by AM via root-to-shoot signaling mediated by CLAVATA/ ESR-related (CLE) peptides and their receptors. In the pea, mutations in genes Sym28 (CLV2-like) and Sym29 (CLV1-like), which encode receptors for CLE peptides, lead to the supernodulation phenotype, i.e., excessive nodule formation. The aim of the present study was to analyze the response of pea cv. 'Frisson' (wild type) and mutants P64 (sym28) and P88 (sym29) to complex inoculation with rhizobia, AM fungi and PGPB, with regard to biomass accumulation, yield and transcriptomic alterations. The plants were grown in quartz sand for 2 and 4 weeks after inoculation with either rhizobia (Rh) or complex inoculation with Rh + AM, Rh + PGPB, and Rh+AM+PGPB, and the biomass and yield were assessed. Transcrip-tome sequencing of whole shoots and roots was performed using a modified RNAseq protocol named MACE (Massive Analysis of cDNA Ends). In the experimental conditions, P88 (sym29) plants demonstrated the best biomass accumulation and yield, as compared to the wild type and P64 (sym28) plants, whereas P64 (sym28) had the lowest rate of biomass and seed yield. The transcriptome analysis showed that both supernodulating mutants more actively responded to biotic and abiotic factors than the wild-type plants and demonstrated increased expression of genes characteristic to late stages of nodule development. The roots of P64 (sym28) plants responded to AM+Rh treatment with upregulation of genes encoding plastid proteins, which can be connected with the activation of carotenoid biosynthesis (namely, the non-mevalonate pathway that takes place in root plastids). The more active response to symbionts in p88 (sym29) plants, as compared to cv. 'Frisson', was associated with counterregulation of transcripts involved in chloroplast functioning and development in leaves, which accompanies successful plant development in symbiotic conditions. Finally, the effect of retardation of plant aging upon mycorrhization on a transcriptomic level was recorded for cv. 'Frisson' but not for P64 (sym28) and P88 (sym29) mutants, which points towards its possible connection with the AON system. The results of this

work link the plant's autoregulation with the responsiveness to inoculation with beneficial soil microorganisms. Keywords: RNAseq, transcriptomics, arbuscular my-corrhiza, nodule bacteria, complex inoculation, autoregulation of nodulation, garden pea

Introduction

Legume plants (family Fabaceae) are able to establish three types of beneficial symbiosis: with rhizobia (Rh), arbuscular mycorrhizal (AM) fungi and plant growth-promoting bacteria (PGPB) (Tikhonovich et al., 2015). In symbiosis with rhizobia, legume plants develop new organs — root nodules — where the rhizobia are converted to symbiotic forms called symbiosomes that fix atmospheric nitrogen (Tsyganova, Kitaeva and Tsyganov, 2018). During symbiosis with AM fungi, the plant allows the fungal hyphae to penetrate the roots and to form symbiotic structures called arbuscules, through which the plant acquires water and sparingly soluble phosphates in exchange for nutrients (Parniske, 2008). The interaction of legumes with PGPB implies the colonization of root surfaces and stimulation of plant growth by compounds excreted by PGPB; PGPB also perform biocontrol over pathogenic microorganisms by producing antibiotics and other compounds (Lugtenberg, Rozen and Kamilova, 2017). Since the formation of beneficial symbioses with microorganisms can improve nitrogen and phosphorus supply and, consequently, increase the yield and quality of grains (Smith, Jakobsen, Gronlund and Smith, 2011; Courty et al., 2015), studying the molecular mechanisms of formation and development of symbioses in legumes is a relevant task for modern genetics.

At the onset of both nitrogen-fixing symbiosis (NFS) and AM, mutual recognition of partners due to signal exchange occurs: plant secondary metabolites (flavonoids in the case of NFS and strigolactones — the carotenoid derivatives — in the case of AM) excreted from the roots attract microorganisms and stimulate them to produce chitin-derived signal molecules, Nod factors and Myc factors, respectively. These Nod and Myc factors are then perceived by specialized plant receptors (Zipfel and Oldroyd, 2017; Leppyanen et al., 2018; Müller et al., 2019). After successful recognition of a microsymbiont, the legume plant activates the signal cascade (also known as CSP, the common symbiosis pathway, because it is shared between NFS and AM), that results in transcriptional changes leading to subsequent development of the symbiotic structures (Old-royd, 2013; Larrainzar et al., 2015).

Another commonality between NFS and AM relates to the regulation of the nodule number and the AM colonization level in the root; these are controlled in legumes by a system known as AON (autoregulation of nodula-tion) (Reid et al., 2011; Wang, Reid and Foo, 2018). The

central role in the AON is played by the Leucine-rich repeat (LRR) receptor kinase homologous to Arabidopsis thaliana CLAVATA1 (CLV1), known as SYM29 in pea (reviewed in Wang et al., 2018). This protein acts as a shoot receptor of root-derived signals, CLAVATA/ESR-related (CLE) peptides. Some members of the CLE peptide family, represented by CLE12 and CLE13 in pea, are synthesized in emerging nodules and are transported through the xylem to shoots, where they serve as a signal of nodulation intensity (Mortier et al., 2010). In shoots, CLE peptides are perceived by receptor complexes that include CLV1 receptor kinase, CORYNE (CRN) receptor pseudokinase and CLAVATA2 (CLV2) receptor-like protein (Roy et al., 2020); these complexes generate a signal of a still unknown nature that suppresses further nodulation. Mutations in genes encoding the proteins that constitute the receptor complexes perceiving CLE peptides lead to supernodulation. An important role in AON is played by the TOO MUCH LOVE (TML1 and TML2) Kelch-repeat F-box proteins, which act as negative regulators of nodulation in roots in response to the shoot-derived signal (Magori et al., 2009).

The garden pea (Pisum sativum L.) is an important legume crop grown worldwide (FAOSTAT 2018) and the oldest object of genetics as a scientific discipline since Gregor Mendel's experiments. The pea has been used as a model object for studying genetic bases of NFS and AM for many years (Borisov et al., 2007; Tsyganov and Tsyganova, 2020). Several mutant lines of pea with different defects in nodulation and mycorrhization, for which the nucleotide sequence and point mutations leading to mutant phenotype have been identified, remain a useful and valuable material for genetic studies (Zhukov et al., 2016). Among them, lines carrying mutations in CLV1 (sym29) and CLV2 (sym28) are available in pea. They were characterized by a supernodulation phenotype, which is accompanied in the case of sym28 by shoot fasciation (which points towards the role of pea CLV2 gene in control of shoot meristem, in addition to nodulation) (Krusell et al., 2011). For CLV1 (sym29) the supermycorrhization phenotype was also described, which was not the case for the sym28 mutant (Morandi, Sagan, Prado-Vivant and Duc, 2000). So far, no studies of complex, multi-component plant-microbial systems in pea supernodulating mutants have been reported, although for some pea varieties and cultivars with normal nodule and mycorrhizal phenotype the effect of dual Rhizobium-AM fungal on proteome and metabolome has been studied (Desalegn, Turetschek, Kaul and Wienkoop, 2016). Here, we investigated the impact of pea AON mutations sym28 and sym29 on biomass, seed yield and transcriptome profiles of plants inoculated with NB, AM and PGPB, and found that the AON system influences the responsiveness of plants to complex interactions with beneficial microorganisms.

Materials and methods

Experimental setup

The seeds of cultivar (cv.) 'Frisson' (wild type) and mutant lines P64 (sym28) and P88 (sym29) (Sagan and Duc, 1996) were surface sterilized with concentrated sulfuric acid, rinsed with sterile water, germinated on wet ver-miculite for three days in darkness at 25 °C, planted in 5 L pots filled with quartz sand (five plants per pot), and inoculated with 150 ml of water suspension (106 CFU * l-1) of Rhizobium leguminosarum bv. viciae RCAM1026 (Afonin, Sulima, Zhernakov and Zhukov, 2017). In variants with AM inoculation, the AM fungus Rhizophagus irregularis strain BEG144, provided by the International Bank for the Glomeromycota (Dijon, France), was used to inoculate pea seedlings. Allium schoenoprasum L. was used as a host plant for R. irregularis cultivation. Fresh roots of A. schoenoprasum colonized by R. irregularis were surface-disinfected as described by Cranenbrouck et al. (2005), cut into 0.5-1 cm segments and used as AM fungal inoculum in this study (0.2 g of inoculum per plant into wells made in the sand before planting). Finally, in variants with PGPB inoculation, 2 ml of water suspension of Arthrobacter mysorens 7 strain (accession number in The Russian Collection of Agricultural Microorganisms (RCAM) (ARRIAM, St. Petersburg) RCAM 01094 (Afo-nin et al., 2021)) in concentration 5 * 106 CFU/ml were poured under each planted seedling. This strain is an active component of 'Mysorin' biopreparation that was previously described as possessing plant growth promoting qualities and was shown to enhance the effect of other biopreparations containing R. leguminosarum (Kozhe-myakov et al., 2015). Thus, the following combinations of inoculation were obtained: 1) rhizobia (Rh), 2) rhizobia + AM (Rh + Myc), 3) rhizobia + PGPB (Rh + PGPB), and 4) rhizobia + AM + PGPB (Rh + Myc + PGPB). Note that the inoculation with rhizobia was applied as an obligatory component in all treatments, since Rhizobium leguminosarum bv. viciae strains capable of forming nodules on pea plants are present in most soils of Russia, and the present experiment was designed to model field conditions.

Before planting, the weight of pots was adjusted to the same value. The plants were grown under non-controlled light, humidity, and temperature conditions in a greenhouse of the All-Russia Research Institute for Agricultural Microbiology, St. Petersburg (June-August 2018). The plants were harvested as follows: three pots per variant were harvested in 2 and 4 weeks post-inoculation (= planting; wpi), and six pots were harvested at the stage of mature seeds (approx. 3 months post inoculation). The dry weight of plants' aerial parts was determined at all time points, the nodule number and AM colonization were analyzed at 4 wpi only, and the weight of seeds and the total number of seeds per plant were recorded at full maturation of the

plants. Several lateral roots (30 cm length) from each pea root system were randomly selected and then subjected to analysis of AM development as described by Shtark et al. (2016). No intraradical mycelium was detected at this time point; however, extraradical mycelium was observed on the surface of the analyzed roots.

Statistical analysis

Statistical analysis of plant growth parameters for three genotypes under different inoculations was performed using Generalized Linear Models (GLMs). For each of the analyzed parameters, we defined the distribution type and the link function (see Table S1)1 which fitted the observed data in the best way based on the model diagnostic, which included visual analysis of residual plots, Cook's distance plots and Q-Q plots. Initial (full) models included all three analyzed factors: (i) Genotype, (ii) PGPB inoculation, (iii) AM inoculation, and their interactions.

P ~ Genotype*PGPB*AM

Using the Likelihood Ratio Test, we tested the significance of model terms and afterwards reduced full models to ones that included three main factors and interactions of genotype with PGPB and genotype with arbuscular mycorrhiza.

P ~ Genotype+PGPB+AM+ Genotype:PGPB+Genotype:AM

For comparisons between genotypes under different inoculation treatments and between treatments for each genotype, we formulated several hypotheses coded in Contrast Matrix that were tested with Wald test using the 'multicomp' R-package. The false discovery rate (FDR) approach was used for multiple comparisons correction in multiple hypothesis testing.

Transcriptomics

Three plants from each variant were removed from the pots, the roots were rinsed with water, and the whole root systems and shoots were immediately frozen in liquid nitrogen. The samples were then ground in liquid nitrogen, and total RNA was extracted with Qiagen Plant Mini kit (Qiagen, Hilden, Germany) and used for sequencing libraries preparation with a MACE kit (GenX-Pro GmbH, Frankfurt-am-Main, Germany) according to the manufacturer's protocol, described in detail by Zhernakov et al. (2019). The libraries were sequenced on an Illumina HiSeq 2500 in Macrogen (Seoul, South Korea). The sequencing results are deposited in the NCBI database (bioproject PRJNA658774, accession numbers SAMN15887753-SAMN15887848).

1 Supplemental material to the article is available at https://biocomm.spbu.ru/article/view/9796

p=0.002

p=0.018

p=0.002

p=0.267

p=0.002

p=0_054

p=0.002

r

p= 1.000

1

p=1.000p=1.000 I I I p=1.000p=1.000 I I I

p=0.489

r

p=1.000

1

p=1.000p=1.000 I I I p=1.000p=1.000 r ' '

T

1

*

r

p= 1.000

p=1.000p=1.000 I I I p=1.000p=1.000 I I I

w

P64 (sym28)

Shoot dry weight, g - final

i

Frisson Genotype

p=0.964799550

P88 (sym29)

p-0000000535

p=0.377024549

p=0.964799550

p=0.000000002

r

p=0.377024549

p=0.964799550

p=0.607337255

p=1.000

' p=1.000p=0.70^ p-0.724 p^l'ooo '

I I I

p=1.000

p=Q.685p=1.00b

p=1.000p=0.632

I I I

1

p=1.000

' p=0-616p=1-00l)

p'=i.ooop=o!ei6 ' I I • "I

4 wf

P64 (sym28)

Frisson Genotype

P88 (sym29)

0.6-

0 4-

p=0.0043

p=0.9999

p=0.0003

p=0.9999

p=0.02CM

p=0.9999

p=0.0003

p= 1.000

p=1.000p=1.000 I I I p=1.000p=1.000 I I I

E*

p=0.9999

p=1.000

p=1.000p=1.000 I I I p=1.000p=1000

p=1.000

I

p=1.000p=1.000 I '

p=1.000p=1.000

I I I

m

P64 (sym28)

Nodule number - 4 wpi

Frisson Genotype

750-

p=0 175

P88 (sym29)

p=0.047

r

p=0.314

p=0.121

p=0.121

p=0.047

p=0.175

p= 1,000

1 p=1.000p=1.00b

p=i.ooop=i.'ooo ' I—I-1

p-0-121

I-

p=1 000 1 p=1.000p=1.00l0

p=1.000p=1.'000 ' I—I-1

p= 1-000 1 p=1-000p=1.Q|j0

p=i.ooop=i!ooo ' '—I—•—I

P64 (sym28)

Frisson Genotype

P88 (sym29)

p= 1.00000000

p= 100000000

p=0.00003224

p= 1-00000000

P=0-00000006

r

p=1.00000000

p=1.000

I p=1-000p=a97|0 p=1.000p=l'000 ' r—I-1

p= 1-00000000

p=1.00000000

p= 1000 I p-0-523 p=ÏÔÔ|o p'=1.000p=0.'519 '

I I I

p= 1.000 j p=Q-523p=1-Qo[) p=1.000p=0.'785 '

I—r-1

p=1e+00

r

p=6e-10

r

p-8e-01

p=8e-13

p=8e-01

p=8e-01

r

p-3e-01

p= 1.000

1 p=1 000p=1.00b

p=1.000p=l'000 ' I I I

T

p-2e-01

p= 1.000

1 p=Q-330p=10oi) p'=1.000p=0'217 '

I I I

p= 1000

1 p=0 217 p=1 oob p'=1.000p=0.l217 '

I I i 1

I

Í

P84 (sym28)

Frisson Genotype

P88 (sym29) P64 (sym28)

Inoculation Epl Rh Epl Rh+PGPB E±3 Rh+AM E-3 Rh+PGPB+AM

Frisson Genotype

P88 (sym29)

Fig. 1. Growth parameters measured for plants of wild-type cv. 'Frisson' and mutants P64 (sym28) and P88 (sym29) Inoculated with rhizobia (Rh), plant growth promoting bacteria (PGPB) and arbuscular mycorrhlzal (AM) fungi In different combinations.

Brackets denote the compared groups; FDR- adjusted p-values correspond to Wald tests based on the generalized regression results.

The quality of the raw reads was evaluated in FastQC (Andrews et al., 2014); the trimming of polyA tails, adapters, and removal of low quality sequences was performed in BBDuk (Bushnell, 2018). Mapping the reads to the reference genome of pea cv. 'Frisson' (Afonin et al., unpublished data) was performed in STAR (ver. 2.7.3 a) (Dobin et al., 2013), with parallel quantification; at least 1.2 mln reads per sample were mapped, with >89.5 % of mapped reads. PCA (Principal Component Analysis) plots were constructed with the use of DESeq2 and ggplot2 (Wick-ham, 2016) packages. After analysis of PCA plots (Sup-pl. Fig. 1 A, B), five outlying samples belonging to the P64 (sym28) mutant were excluded. Differential gene expression analysis was performed using the DESeq2 package (ver 1.26.0) (Love, Huber and Anders, 2014); the genes were considered differentially expressed if the Wald test was passed with a q-value below 0.05 and log2 fold change below 0.5. GO terms for genes of interest were obtained using the Trinotate suite (Bryant et al., 2017). GO enrichment analysis (using the weight01 algorithm and Fisher's exact test) and further visualization were performed in the packages topGO (Alexa and Rahnenfuhrer, 2009) and ggplot2, respectively. To estimate the expression level of individual genes, the raw counts of aligned reads were normalized by the counts per million method (CPM), logarithmed, and transformed into a z-scale.

Results

Growth and yield parameters of P64 (sym28), P88 (sym29) and cv. 'Frisson' (wild type)

The shoot biomass accumulation was different in all tested lines. Both sym28 and sym29 mutants had less shoot weight at 2 wpi compared to the wild-type cv. 'Frisson'. Nevertheless, at 4 wpi the P88 (sym29) shoot biomass did not statistically differ from the shoot biomass of 'Frisson' plants under any type of inoculation. The final shoot biomass (10 wpi) of the P88 (sym29) mutant was higher than that of the wild-type 'Frisson'. That difference between wild type and the P88 mutant was statistically significant for both AM-inoculated groups, with and without PGPB inoculation. At the same time P64 (sym28) still gained less biomass in shoots in comparison with the wild-type cv. 'Frisson'.

The separate influence of AM and PGPB (with rhi-zobia as the background factor) on growth and yield parameters of the genotypes was assessed (Fig. 1). The comparisons of inoculated and non-inoculated plants showed no significant effect of PGPB under experimental conditions for all three genotypes. At the same time, the influence of AM fungi was considerable: AM slightly reduced the seed number and seed weight in cv. 'Frisson' but increased these parameters in the P88 (sym29) mutant (Fig. 1).

Linear approximation of shoot biomass accumulation at 2 and 4 wpi time points suggests that the P64 (sym28) mutant accumulated biomass more slowly, thus gaining less biomass at the end of the experiment than cv. 'Frisson' and P88 (sym29). P88 (sym29) demonstrated a slight delay as compared to cv. 'Frisson' at initial stages but had the highest shoot biomass at the end of the experiment. The difference between cv. 'Frisson' and P88 (sym29) was statistically significant for the mycorrhi-zal plants. Thus, P88 (sym29) plants in the experimental conditions had higher shoot growth rate at some stages, as compared to the wild type and the P64 (sym28) plants.

As expected, both P64 (sym28) and P88 (sym29) mutants formed an increased number of nodules in comparison with the wild-type cv. 'Frisson' at 4 wpi (106 % and 142 % more, respectively) (Fig. 1). Inoculation with PGPB and AM showed no significant effect on the number of nodules formed in all three lines. Furthermore, we did not observe any differences in the root biomass between lines or treatment groups.

PCA (Principal Component Analysis) of 3'MACE-seq data

The 3'MACE-seq reads were mapped onto the reference genome after the quality control and trimming step, and the resulting gene expression data was analyzed. The PCA plots were built for root and shoot samples separately (Fig. S1A, S1B). The first component (80 % of variance for roots and 60 % for shoots) clearly distinguished time points (2 and 4 wpi) and may therefore be considered as the time scale. The second component (12 % and 20 % for roots and shoots, respectively) distinguished mutants from wild type and, at the same time, P64 (sym28) mycorrhizal samples from non-mycorrhizal ones at 2 wpi. For the shoot samples of cv. 'Frisson' at 2 wpi, a shift of non-mycorrhizal samples towards the group of older samples (4 wpi), as compared to the mycorrhizal samples, was noticed.

Transcriptomic response to mono-inoculation with rhizobia (Rh)

In the variant of mono-inoculation with rhizobia, the expression profiles of two supernodulating mutants P64 (sym28) and P88 (sym29), as expected, differed significantly from that of the wild-type cv. 'Frisson' and shared a substantial amount of common differentially expressed genes (DEGs) (Fig. 2 and Fig. S2 A, B, C). The P88 (sym29) mutant demonstrated more unique DEGs, as compared to P64 (sym28), in shoots at 2 and 4 wpi time points, and in roots at the 2 wpi time point (Fig. S2 A, B, C). Functional annotation of DEGs showed that, in terms of biological process, the roots of mutants at 2 wpi share the increased response to several factors, as compared to the roots of the wild-type plants. Namely, the groups of

tyrosine catabolic process translation response to wounding response to water deprivation response to salt stress response to salicylic acid response to reactive oxygen species response to oxidative stress response to osmotic stress response to oomycetes response to molecule of bacterial origin response to light stimulus response to jasmonic acid response to hydrogen peroxide response to freezing response to cold response to chitin response to auxin response to abscisic acid regulation of transcription, DNA-templated regulation of systemic acquired resistance regulation of stomatal movement regulation of seedling development regulation of jasmonic acid mediated signaling pathway regulation of growth regulation of defense response reductive pentose-phosphate cycle protein ubiquitination protein polyubiquitination protein folding protein complex oligomerization protein autoubiquitination protein autophosphorylation protein-chromophore linkage pollen tube guidance plant-type secondary cell wall biogenesis plant-type hypersensitive response photosynthesis, light harvesting in photosystem I photosynthesis, light harvesting photosynthesis photoprotection

negative regulation of transcription, DNA-templated multicellular organism growth mucilage biosynthetic process involved in seed coat development

mRNA processing maturation of LSU-rRNA maturation of 5.8S rRNA L-phenylalanine catabolic process glucosinolate metabolic process glucosinolate biosynthetic process ethylene-activated signaling pathway defense response to fungus defense response to bacterium, incompatible interaction defense response to bacterium defense response by callóse deposition in cell wall defense response cytoplasmic translational termination cellular response to nickel ion cellular response to hypoxia cellular aromatic compound metabolic process cell morphogenesis amino acid transmembrane transport abscisic acid-activated signaling pathway

Common DEG

1

j

• m

m

*

*

*

*

m

1 W

w

P64 un que DEG

P88 un que DEG

m

m m

•

m

m

•

A

• m

V

W

A

V

9

A

•

A

W

w A

w

•

-•

p-va lue

Number of genes

• 25

• 50

• 75

A loo

Gene ratio

Fig. 2. Gene Ontology enrichment analysis of DEGs in roots of cv. 'Frisson', P64 (sym28) and P88 (sym29) in 2 wpi with rhizobia.

B

en <D o O

CD O O) O O

o

O

transcription, DNA-templated signaling signal transduction RNA biosynthetic process response to stress response to stimulus response to oxygen-containing compound response to other organism response to organic substance response to lipid response to inorganic substance response to hormone response to external biotic stimulus response to endogenous stimulus response to drug response to chemical response to biotic stimulus response to bacterium response to alcohol response to acid chemical response to abscislc acid response to abiotic stimulus regulation of transcription, DNA-templated regulation of RNÂ biosynthetic process regulation of response to stimulus regulation of primary metabolic process regulation of nucleobase-containing compound metabolic process regulation of nucleic acid-templated transcription regulation of nitrogen compound metabolic process regulation of macromolecule biosynthetic process regulation of gene expression regulation ofcellular process regulation of cellular metabolic process regulation of cellular macromolecule biosynthetic process regulation of biological process oxidation-reduction process organic cyclic compound biosynthetic process nucleobase-containing compound biosynthetic process nucleic acid-templated transcription multi-organism process hormone-mediated signaling pathway heterocycle biosynthetic process defense response to other organism defense response cellular response to oxygen-containing compound cell communication biosynthetic process biological regulation aromatic compound biosynthetic process

DEG iri Rh+AM inoculation compared to Rh inoculation in Frisson and P64 (S, 2 wpi)

2 O i—

O o

o >

Common DEG

Friss on uniqu e DEG

P6 4 unique DEG

jt

W m i

1 & m 1

a ' 9

W A

A w

w A

w

A ■

w ±

\ »

A •

w

A

• a

m 9

W »

V

A

W A

w

t

m ' P

w

A

A W

w

A

w

A

••

• t 1

r A

W

0 1 0 2 0 3 0 4 o.;

p-value

0.002

0.001

Number of genes

• 200

# 400

% 600

A 800

Gene ratio

Fig. 3. Gene Ontology enrichment analysis of DEGs In roots (A) and shoots (B) of cv. 'Frisson' and P64 (sym28) In 2 wpi with mycorrhlzal fungus and rhlzobla.

DEGs included responses to: 1) environmental factors, 2) biotic factors such as response to chitin (probably, to Nod factors), 3) phytohormones (auxin, salicylic acid, ab-scisic acid, jasmonic acid), and 4) the defense response to other organisms. Other biological processes over-represented in both mutants were related to regulation of transcription (both positive and negative), protein autoubiq-uitination, multicellular organism growth, and ethylene-and abscisic acid-activated signaling pathways, including genes related to induced systemic resistance.

A few DEGs unique to P64 (sym28) roots were related to protein folding, while in roots of P88 (sym29), more numerous unique DEGs corresponded to the same GO terms as those shared with P64 (sym28), pointing at a more active response to the factors (including the response to bacterium) and associated activation of the related hormone-mediated signaling pathway and cell morphogenesis. Remarkably, a large group of genes related to cellular aromatic compound metabolic process was found to be specific to P88 (sym29) roots.

At the first and second time points, 2 and 4 wpi, in roots, the genes related to nodulation and nodule morphogenesis were upregulated in both mutants, in line with the supernodulation phenotype which implies the increase of the expression of nodule-related genes (Sup-pl. Fig. 2A). In shoots of mutants in 2 wpi, similarly to the roots, the common DEGs were related to responses and regulation (Fig. S2B), and P64 (sym28) shoots, similarly to the roots, over-expressed genes related to protein folding. Shoots of P88 (sym29) showed more pronounced response to biotic and abiotic stimuli than in cv. 'Frisson' and P64 (sym28) in 2 and 4 wpi (Fig. S2C).

Transcriptomic response to

combined inoculation with Rh, PGPB and AM

The influence of factors Rh+AM and Rh+PGPB was assessed by a two-factor analysis. Similarly to the results of the plant growth parameter analysis, no significant changes in gene expression were detected after PGPB treatment in any of the studied variants. The effect of AM fungi was clearly seen at 2 wpi for cv. 'Frisson' and P64 (sym28) (Fig. 3), and was less pronounced in the case of P88 (sym29) (Fig. S3). At the second time point (4 wpi), unfortunately, the high variability and dispersion detected in the data (probably due to growth in non-controlled conditions) blurred the effect of AM inoculation so that fewer DEGs could be detected.

Transcriptomic response of P64 (sym28) to Rh+AM inoculation

The gene expression upon AM inoculation was specifically examined in cv. 'Frisson' and P64 (sym28) roots and shoots at 2 wpi. The functional annotation of DEGs

in roots highlighted the processes common for both genotypes, such as several types of responses and regulation of transcription and metabolism. The differences between cv. 'Frisson' and P64 (sym28) were mainly related to genes involved in response to another organism (specifically, to a fungus), response to lipids and nitrogen compounds, generation of precursor metabolites and energy, and oxidation-reduction process (Fig. 3A), which indicates a more active response to AM fungus in P64 (sym28). Also, some genes annotated as encoding the chloroplastic products were significantly upreg-ulated in P64 (sym28) roots upon mycorrhization; this may be connected with carotenoid biosynthesis, since some reactions of this process proceed in plastids in both shoots and roots (Strack and Fester, 2006). Interestingly, in roots of P88 (sym29), the same genes had a high expression level regardless of inoculation type, and at 4 wpi, both mutants also showed upregulation of these genes compared to the wild type.

In shoots of these genotypes at 2 wpi, the pattern of DEGs was similar to that obtained for the roots, with common DEGs related to responses, regulation and the hormone-mediated signaling pathway (Fig. 3B). The shoots of P64 (sym28) were more responsive to AM inoculation, as there were more DEGs, which were distributed among the same groups of biological processes. The unique DEGs were related to response to lipids, organic substances, alcohol, abscisic acid, and the oxidation-reduction process.

Transcriptomic response of P88 (sym29) and cv. 'Frisson' to triple inoculation

The plants' response to triple inoculation (Rh + AM + PGPB) was studied using the cv. 'Frisson' and P88 (sym29) data (because some P64 (sym28) samples had been discarded after assessment by PCA, and the remaining number of samples was insufficient for this analysis). Remarkably, the response to triple inoculation, as compared to mono-inoculation with rhizobia, was different for wild type and supernodulating mutant P88 (sym29) in roots in 2 wpi (Fig. S4 A). In roots, in 4 wpi the response was weaker and less specific (Fig. S4 B). In shoots, in 2 wpi the response was almost indistinguishable (Fig. S4 C), but in 4 wpi P88 (sym29) showed a significant response, which, intriguingly, was opposite to that of cv. 'Frisson' (Fig. S4 D). In total, 40 genes were found to be counterregulated: 30 were upregulated in P88 (sym29) and 10 were downregulated. The upregu-lated genes are mainly involved in chloroplast functioning and development (e.g., those encoding RUBISCO, CURVATURE THYLAKOID and chloroplastic thiore-doxin and ferredoxin), which is in agreement with active growth of plant shoots at this time point. Interestingly, these genes have low expression level in shoots of

Absolute expression values (in CPM) of putative markers (S, 4 wpi)

--162.27 "19? 7R 1W RR 107? IgSeS 1¿LCI 1R9 94 rinSI p.

max

157.15 | 85.45 172.4 177.32 130.26

241.99

139.88 153.25 102.38 564.46 585.05 463.58

123.78 135.88 104.63 107.2

66.38 | 87.45 58.19 68.81

88.46 101.27 105.35 97.38

136.38 130.5 150.33 157.9 19.7 | 34.35 | 13.28 | 21.5

369.08 342.89 365.24 373.31

115.19 125.07 85.07 119.31

94.96 | 70.02

70.97 84.32~ 86.19 89.56

14.46 15.8~_i_

253.15 247.51 194.49 295.51

192.3 192.02 152.5 248.18

69.44 139.14 83.01 112.03

13744.334168.14 3732.5 114.17 125.55 92.11

235.86 234.87 165.42 1606.781527.98 1342.9 1235.22 996.99 1041.84 904.47

341.25 448.61 383.58 454.84 2345.363240.542849.873229.51 71.31 99.54 97.34 86.74 153.11 1191.79 | 162.1 | 172.81 1274.171493.561319.481694.04 752.8 1114.07 938.44 990.42

L

1068.59 890.68 909.09 777.24

| 123.97 173.31 127.71

| 224.77 246.26 220.86 J 1442.51 1519.631304.24

| 324.05 316.36 216.85

| 90.79 95.13 101.64

279.16 387.13 359.07 211.1 232.57 274.57 299.35

144.61 168.94 153.71 235.72 363.8

688.56 925.94 851.85 715.83

74.02 87.13 85.46 97.26

148.24 194.81 144.93 189.53

756.19 977.9 827.29 986.87

256.31 | 4 316.36 238.96

101.64 | 9 102.31 95.27

I

385.68 418.14 627.95

445.73 380.79 507.44

445.74 379.48 545.9 206.81 170.86 251.43 29.84 31.08 60.75 14.54 13.67 33.55

£

O a

§

it

.w it

GDSL esterase/lipase APG

Alkane hydroxylase MAH1 {ECO:0000303|PubMed:17905869} Defensin-like protein Defensin-like protein 230

Very-long-chain aldehyde decarbonylase CER1 {ECO:0000305}

1-aminocyclopropane-1-carboxylate oxidase 3

B-box zinc finger protein 25 {ECO:0000303|PubMed:19920209}

Peroxidase 4 {ECO:0000250|UniProtKB:Q42578}

Probable non-specific lipid-transfer protein AKCS9

Cysteine protease XCP2 {ECO:0000305}

Triacylglycerol lipase 2

Chalcone synthase

Fasciclin-like arabinogalactan protein 11 MYB-like transcription factor ETC1 Thioredoxin M-type, chloroplastic

Ribulose bisphosphate carboxylase small chain, chloroplastic Defensin-like protein 39

Protein CURVATURE THYLAKOID 1B, chloroplastic

Germin-like protein subfamily 3 member 1

Photosystem I reaction center subunit IV A, chloroplastic

Oxygen-evolving enhancer protein 3-2, chloroplastic

pEARLI1-like lipid transfer protein 3

Proline-rich protein 4

Ferredoxin-1, chloroplastic

Chlorophyll a-b binding protein AB80, chloroplastic

Organ-specific protein S2

Probable aquaporin PIP-type 7a

Probable GTP diphosphokinase RSH3, chloroplastic

Uncharacterized protein At4g22758

S-adenosylmethionine decarboxylase proenzyme

Arginine decarboxylase

Protein DEK

Ribosomal protein S3, mitochondrial Eukaryotic translation initiation factor 5B DEAD-box ATP-dependent RNA helicase 42 DnaJ homolog subfamily C member 2

I

min

Fig. 4. The expression profile of the top 25 genes, which were counterregulated in P88 (sym29) shoots, as compared to cv. 'Frisson', in 4 wpi.

E

(/>

<D

99.!

26.98 13.92 5.67

i<>

Absolute expression values (in CPM) of Sym genes (R, 2 wpi)

89.23 84.9 3.5 0.96 21.95

31.63

26.4

262.05

110.64 114.83 120.26

3.17 3.06 3.44

0 0 0.11

28.74 35.49 33.45

1.87 2.04 2.81

0 0 3.46

Jo

81.64

112.19

0 0.16

2.85 3.29

0 0

34.02 33.96

38.84

38.02 21.35

8.94 5.43

89.65 104.65

134.73 176.2

5.51 6.42

1.9 1.64

36.87 33.75

153.18 170.99

51.7 56.61

112.16 93.01

14.64 12.74

128.43

max

108.26

11.59 8.77

2.66 1.32

12.22 13.11 7.24

3.75 3.42 3.07

34.17 33.25 27.48

26.59 24.33

14.06

3.87

Sucrose synthase Nodulation protein H

Protein NODULATION SIGNALING PATHWAY 2 {ECO:0000303|PubMed: 15961668}

Chitin elicitor receptor kinase 1 {ECO:0000303|PubMed:21070404}

LysM domain receptor-like kinase 3 {ECO:0000303|PubMed:12947035}

Phosphoethanolamine N-methyltransferase

Myb family transcription factor IPN2 {ECO:0000305}

Probable ion channel SYM8

Nodulin-26

Ethylene-responsive transcription factor ERN1 {ECO:0000305} Probable esterase D14L

Ethylene-responsive transcription factor ERN2 {ECO:0000305}

Putative E3 ubiquitin-protein ligase LIN—1 {ECO:0000250|UniProtKB:D1FP53}

Probable receptor-like protein kinase At1g30570

Early nodulin-12A

Nodulation receptor kinase

Early nodulin-93

Uricase-2

Serine/threonine receptor-like kinase NFP {ECO:0000303|PubMed:16723404} Probable 2-isopropylmalate synthase

Glycerol-3-phosphate acyltransferase RAM2 {ECO:0000305} Nodulin-13

Protein CYCLOPS {ECO:0000303|PubMed:19074278} NF_YA1

Nod factor hydrolase protein 1 {ECO:0000303|PubMed:24082029}

GRAS family protein RAMI {ECO:0000305}

Nodulin-25

Flotillin—like protein 1

Early nodulin-70

Early nodulin-16

Early nodulin-75

Early nodulin-5

Early nodulin-12B

Nodulin-1

NF_YA2

Putative E3 ubiquitin-protein ligase LIN-2 {ECO:0000303|PubMed:19776163} Protein CNGC15a {ECO:0000303|PubMed:27230377} Calcium and calcium/calmodulin-dependent serine/threonine-protein kinase U-box domain-containing protein 1 {ECO:0000303|PubMed:20971894} Flotillin—like protein 4

BTB/POZ domain and ankyrin repeat-containing protein NOOT2 {ECO:0000305}

I

I

min

B

Absolute expression values (in CPM) of Sym genes (R, 4 wpi)

max

i

r{ -{

\

i \

IE

C

c (/) f* <D

(/)

180.24 2.72 23.3 51.06

-X

2.6

15 372.68 22.9 184.72 1838.76 54.95 1762.04 1033.7 282.65 4880.9 141.12 17.37

15.72 339.04 19.94 186.76 1976 47.56 1787.46 1142.72 289.44 3561.39 134.98 15.37 28.14 2243.01 237.15 78.45 7092.35 210.32 52.25

26.59 206.76 2410.34 70.76 1844.17 1283.25 326.17 4773.74 166.03 19.74 25.81 2307.67 302.43 75.23 5549.72 188.59 49.84 105.14 33.58 0

14.86 429.54 17.14 169 1800.94 51.19 1890.54 1096.15 272.71 3264.44 135.36 16.02 25.02 1708.17 206.95 89.01 6191.85 185.31 33.96 104.63 22.27 0

76.45 0.09 Ai 259.74

0.25 0

0.18 0.16

0.7

| 161.28 156.7

1.16 1.38

12.14 15.58

25.96

5.29

177.85 173.62 176 |

14.47 10.7

19.09 16.44

44.64 43.12

8.97 8.88

<f* x^ rt* <r

Putative E3 ubiquitin-protein ligase LIN-1 {ECO:0000250|UniProtKB:D1FP53} Nodulin-13

U-box domain-containing protein 1 {ECO:0000303|PubMed:20971894}

Uricase-2

Early nodulin-12B

Early nodulin-70

Early nodulin-5

Early nodulin-16

Nodulin-1

Early nodulin-75

Probable 2-isopropylmalate synthase

BTB/POZ domain and ankyrin repeat-containing protein NOOT2 {ECO:0000305} Flotillin—like protein 1 Early nodulin-12A NF_YA1

Calcium and calcium/calmodulin-dependent serine/threonine-protein kinase Nodulin-25

Protein CYCLOPS {ECO:0000303|PubMed: 19074278}

Nod factor hydrolase protein 1 {ECO:0000303|PubMed:24082029}

Phosphoethanolamine N-methyltransferase

Nodulin-26

Flotillin—like protein 4

Early nodulin-93

Protein CNGC15a {ECO:0000303|PubMed:27230377}

Putative E3 ubiquitin-protein ligase LIN-2 {ECO:0000303|PubMed:19776163}

NF_YA2

| Nodulation protein H Probable ion channel SYM8

Serine/threonine receptor-like kinase NFP {ECO:0000303|PubMed:16723404} Glycerol-3-phosphate acyltransferase RAM2 {ECO:0000305} Ethylene-responsive transcription factor ERN1 {ECO:0000305} GRAS family protein RAMI {ECO:0000305}

Protein NODULATION SIGNALING PATHWAY 2 {ECO:0000303|PubMed:15961668} Probable receptor-like protein kinase At1g30570 Sucrose synthase

Ethylene-responsive transcription factor ERN2 {ECO:0000305} | Nodulation receptor kinase Myb family transcription factor IPN2 {ECO:0000305} Chitin elicitor receptor kinase 1 {ECO:0000303|PubMed:21070404} Probable esterase D14L

LysM domain receptor-like kinase 3 {ECO:0000303|PubMed:12947035}

I

I

min

Fig. 5. Expression profiles of known symbiotic genes In plant roots in Rh and Rh+AM variants In 2 wpl (A) and 4 wpi (B).

P64 (sym28) (Fig. 4), in line with poor development of these plants. Thus, this set of shoot-expressed genes can be suggested as markers of successful development of pea plants under triple inoculation.

Differential expression of known symbiotic genes

The expression of genes related to nodulation and to AON was examined in detail with the use of MACE-seq data. In 2 and 4 wpi, roots of sym28 and sym29 mutants, as compared to the wild type, demonstrated elevated expression of genes specific to nodules (i.e., encoding nod-ulins Enod12A, Enod12B, Enod5, Enod16 and Enod75, as well as transcription factors CYCLOPS/IPD3, NF-YA1 and NF-YA2, which corresponds to the supernodu-lation phenotype of these mutants (Fig. 5A, B). These genes are related to infection and nodule development, and their upregulation in the roots reflects the loss of the ability to limit the extent of symbiosis via a negative feedback loop in AON mutants. At the same time, at 4 wpi, Sym genes involved in signal transduction at early stages of the legume-rhizobial symbiosis were less expressed in root systems of the mutants than in that of cv. 'Frisson' (Fig. 5B). This decrease in early symbiotic gene expression may be connected with AON-independent local control of nodulation in the roots of sym28 and sym29 mutants that already formed an excessive number of nodules.

The expression of genes related to autoregulation of nodulation was also examined. At 4 wpi, TML1 and TML2 encoding negative regulators of nodulation were strongly expressed in roots of wild-type plants, but were decreased in the roots of P64 (sym28) and P88 (sym29) mutants defective in CLV1 and CLV2 components of AON (Fig. S5A). Similarly, NNC1 encoding a transcriptional repressor that disturbs CLE gene expression was decreased in the roots of both mutants at 4 wpi. At the same time, CEP1 suppressing TML1, TML2 expression via CRA2 receptor kinase was upregulated in the roots of P64 (sym28) and P88 (sym29) mutants at 2 and 4 wpi.

The effect of mycorrhization was also estimated in wild-type and mutants. At 2 wpi, mycorrhization caused downregulation of CRE1 and CEP1 in roots of all lines, and CRA2 and CORYNE in the roots of AON mutants only (Fig. S5A). At 4 wpi, CRE1 and CEP1 were strongly expressed in P88 (sym29) roots and were downregu-lated upon mycorrhization, and NNC1, CORYNE and CRA2 were upregulated in roots of cv. 'Frisson' but not in AON mutants, except for CRA2 in P88 (sym29) (Suppl. Fig. 5B). In shoots, at 2 and 4 wpi the expression of CLV1, regardless of AM inoculation, was decreased in P88 (sym29) and strongly increased in P64 (sym28), as compared to wild type (Fig. S5C, D). Expression of CRA2 and CORYNE was higher in wild-type shoots than in shoots of both AON mutants.

Discussion

In this study, plants of P64 (sym28) and P88 (sym29) were compared with plants of parental genotype cv. 'Frisson' (wild type) after growth in quartz sand under inoculation with beneficial microorganisms in the following combinations: 1) rhizobia (Rh), 2) rhizobia + AM (Rh + Myc), 3) rhizobia + PGPB (Rh + PGPB), and 4) rhizobia + AM + PGPB (Rh + Myc + PGPB). In previous experiments, the growth parameters of supernodu-lating mutants obtained from cv. 'Frisson' did not exceed that of the parental line (Sagan, Ney and Duc, 1993; Mo-randi, Sagan, Prado-Vivant and Duc, 2000; Salon et al., 2001; Bourion et al., 2007), whereas in the present study P88 (sym29) plants showed superiority over cv. 'Frisson' and P64 (sym28). This may have happened because the experimental conditions of the present study, which appeared to be favorable for P88 (sym29), did not replicate the conditions of previous experiments (where double inoculation with Rh+AM was studied in field conditions, and the effect of only mono-inoculation with AM or Rh was evaluated in sand or another substrate). The advantage of P88 (sym29) may be linked to the nature of the mutation, which causes amino acid change L290F in one of the LRR repeats (Krusell et al., 2002), which negatively affects the receptor functioning but, possibly, doesn't fully block its activity. In turn, P64 (sym28) carries a strong allele with a mutation causing a preliminary stop codon (W456Stop) in the LRR region (Krusell et al., 2011), which leads to a lack of transmembrane domain and, apparently, to a total dysfunction of the protein.

The supernodulation phenotype and corresponding transcriptomic changes in roots and shoots of P64 (sym28) and P88 (sym29) mutants impaired in the CLV receptors reflect the breakdown of the autoregulation mechanism. Indeed, the nodule number in P64 (sym28) and P88 (sym29) was 2 and 2.5 times higher, respectively, than in wild type, and this ratio was not significantly affected by any inoculation type. The presence of excessive nodules was clearly seen on a tran-scriptomic level: the nodule-specific genes had a higher expression level in roots of both mutants than in wildtype roots. At the same time, the roots of the mutants showed a decrease in expression level of several genes involved in AON and encoding negative regulators of nodulation such as TML1, TML2 and NNC1, which may indicate continuous nodule formation in sym28 and sym29 mutants due to loss of the system control of nod-ulation. Indeed, the loss of balanced control of nodula-tion may trigger upregulation of genes controlling plant response to environmental and biotic factors and the defense response to other organisms, as it was shown for AON mutants.

The elevated response to plant hormones such as salicylic acid, abscisic acid and jasmonic acid, and

downregulation of CRE1 and CEP1 in roots of all lines, including the AON mutants, which occur upon mycor-rhization, point at the activity of alternative regulatory systems that may compensate for the defects in CLV-based autoregulation of nodulation in both mutants. The fact that the genes CRE1 and CEP1 were downregulated upon mycorrhization in P88 (sym29) roots also suggests that the local control over microsymbionts is not affected in this AON mutant.

Transcriptome profiles show that both supernodu-lating mutants more actively responded to biotic and abiotic factors as compared to wild type, and they demonstrated a more pronounced defense response. This active response at early time points was associated with slower biomass accumulation in both mutants. However, at the final time point, the P88 (sym29) genotype overtook the others, indicating that the high responsiv-ity to inoculation (or, the strategy of relying on symbioses rather than refusing them) may, in some conditions, be advantageous.

The effect of inoculation with AM fungi and PGPB on top of rhizobia appeared to be minor in our experiment. With regard to PGPB, this can be explained by the use of a suboptimal concentration of PGPB strain in the inoculum and/or due to incompatibility of this strain with pea genotype cv. 'Frisson' and the corresponding mutant lines. Also, the set-up of the experiment in non-controlled temperature and humidity conditions obviously resulted in high variation between the samples, which might have masked the slight effects of the tested PGPB strain on gene expression, if it indeed showed any. The effect of AM inoculation, however, was more pronounced (although at early time points the internal mycelium had not yet formed) so that we detected alterations in gene expression connected with AM+Rh treatment. Despite the fact that there was no intraradical mycelium in the roots, the extraradical mycelium could influence the plant. It is known that plants are capable of perceiving molecular signals sent by the germinated spores and external mycelium of AM fungi at the pre-symbiotic stage of mycorrhiza development (Nadal and Paszkowski, 2013). In addition, fungal mycelium can release substances that stimulate plant growth (Felten et al., 2009; Splivallo et al., 2009).

Interestingly, the strongest effect of AM+Rh treatment was characteristic for the P64 (sym28) mutant, in the roots of which we detected the upregulation of genes encoding plastid proteins. Plastid proteins are known to play a significant role in mycorrhization, since some stages of carotenoid biosynthesis (namely, the non-mevalonate pathway) occur in plastids (Walter, Fester and Strack, 2000; Strack and Fester, 2006). In shoots of P64 (sym28) we detected a similarly strong effect of AM+Rh treatment, which may reflect a different developmental status of mycorrhizal and non-mycorrhizal

plants (for example, retardation of transition to flowering).

Intriguingly, on the PCA plots the non-mycorrhizal samples of shoots of the wild-type plants at 2 wpi were shifted towards the 4 wpi samples, as compared to my-corrhizal samples, indicating a delay in plant development for cv. 'Frisson' under mycorrhization. The same effect of retardation of plant aging upon mycorrhization was described in our previous work with the use of me-tabolomics (Shtark et al., 2019). It is tempting to further investigate the transcriptomics data from this experiment in order to detect genes whose expression can serve as a marker of plant's rejuvenation under influence of AM, but we refrain from doing so until more samples from repeated experiments are obtained. Interestingly, this effect was not seen in the P88 (sym29) mutant and was scarcely detected in the P64 (sym28) mutant as well, which points at its possible connection with the autoregulation system.

Conclusion

This study added new information on phenotyp-ic characterization of pea supernodulating mutants P64 (sym28) and P88 (sym29) in symbiosis with beneficial soil microorganisms. Namely, it was found that in some conditions AON defects leading to excessive nodule formation might be advantageous; however, this effect should be studied in detail in the controlled environment of a growth chamber. AON mutants appeared to be more responsive to inoculation than the wild-type plants, and the strongest reaction to mycorrhization was described for P64 (sym28), for which the activation of plastid metabolism in roots and modulation of plant immune reactions in shoots was detected. Particular groups of genes differentially expressed in AON mutants regardless of the inoculation type were also described, such as genes related to protein folding in P64 (sym28), and genes related to secondary metabolite production in P88 (sym29). In general, the results of this work link the activity of the plant autoregulation system with the plant's responsiveness to inoculation with beneficial soil microorganisms.

References

Afonin, A., Sulima, A., Zhernakov, A., and Zhukov, V. 2017. Draft genome of the strain RCAM1026 Rhizobium legumi-nosarum bv. viciae. Genomics Data 1 1:85-86. https://doi. org/10.1016/j.gdata.2016.12.003 Afonin, A.M., Gribchenko, E.S., Akhtemova, G.A., Lak-tionov, Y. V., Kozhemyakov, A. P., and Zhukov, V.A. 2021. Complete genome sequence of the bacterial component of mysorin biopreparation. Microbiology Resource Announcements 10(11): e01287-20. https://doi.org/10.1128/ MRA.01287-20

Alexa, A. and Rahnenfuhrer, J. 2009. Gene set enrichment analysis with topGO. Bioconductor Improv 27.

Andrews, S.F., Krueger, F., Seconds-Pichon, A., Biggins, F., and Wingett, S.F. 2014. A quality control tool for high throughput sequence data. Babraham Bioinformatics. Borisov, A.Y., Danilova, T.N., Koroleva, T.A., Kuznetsova, E.V., Madsen, L., Mofett, M., Naumkina, T.S., Nemankin, T.A., Ovchinnikova, E.S., Pavlova, Z.B., Petrova, N.E., Pinaev, A.G., Radutoiu, S., Rozov, S.M., Rychagova, T.S., Shtark, O.Y., Solovov, 1.1., Stougaard, J., Tikhonovich, I.A., Topunov, A. F., Tsyganov, V. E., Vasil'chikov, A. G., Voroshi-lova, V. A., Weeden, N. F., Zhernakov, A. I., and Zhukov, V. A. 2007. Regulatory genes of garden pea (Pisum sativum L.) controlling the development of nitrogen-fixing nodules and arbuscular mycorrhiza: A review of basic and applied aspects. Applied Biochemistry and Microbiology 43(3):237-243. https://doi.org/10.1134/S0003683807030027 Bourion, V., Laguerre, G., Depret, G., Voisin, A.-S., Salon, C., and Duc, G. 2007. Genetic variability in nodulation and root growth affects nitrogen fixation and accumulation in pea. Annals of Botany 100(3):589-598. https://doi. org/10.1093/aob/mcm147 Bryant, D.M., Johnson, K., DiTommaso, T., Tickle, T., Coug-er, M.B., Payzin-Dogru, D., Lee, T.J., Leigh, N.D., Kuo, T.-H., and Davis, F. G. 2017. A tissue-mapped axolotl de novo transcriptome enables identification of limb regeneration factors. Cell Reports 18(3):762-776. https:// doi.org/10.1016/j.celrep.2016.12.063 Bushnell, B. 2018. BBTools: a suite of fast, multithreaded bioinformatics tools designed for analysis of DNA and RNA sequence data. Jt Genome Inst. Courty, P. E., Smith, P., Koegel, S., Redecker, D., and Wipf, D.

2015. Inorganic nitrogen uptake and transport in beneficial plant root-microbe interactions. Critical Reviews in Plant Sciences 34(1-3):4-16. https://doi.org/10.1080/073 52689.2014.897897

Cranenbrouck, S., Voets, L., Bivort, C., Renard, L., Strullu, D.-G., and Declerck, S. 2005. Methodologies for in vitro cultivation of arbuscular mycorrhizal fungi with root organs. In Declerck, S., Fortin, J. A., and Strullu, DG. (eds) In Vitro Culture of Mycorrhizas. Soil Biology, vol.4. Springer, Berlin, Heidelberg. https://doi.org/10.1007/3-540-27331-X_18 Desalegn, G., Turetschek, R., Kaul, H.-P., and Wienkoop, S.

2016. Microbial symbionts affect Pisum sativum pro-teome and metabolome under Didymella pinodes infection. Journal of Proteomics 143:173-187. https://doi. org/10.1016/j.jprot.2016.03.018

Dobin, A., Davis, C.A., Schlesinger, F., Drenkow, J., Zaleski, C., Jha, S., Batut, P., Chaisson, M., and Gingeras, T. R. 2013. STAR: ultrafast universal RNA-seq aligner. Bioinformatics 29(1):15-21. https://doi.org/10.1093/bioinformatics/ bts635

Felten, J., Kohler, A., Morin, E., Bhalerao, R.P., Palme, K., Martin, F., Ditengou, F.A., and Legué, V. 2009. The ec-tomycorrhizal fungus Laccaria bicolor stimulates lateral root formation in poplar and Arabidopsis through auxin transport and signaling. Plant Physiology 151(4):1991-2005. https://doi.org/10.1104/pp.109.147231 Kozhemyakov, A. P., Laktionov, Y.V, Popova, T. A., Orlova, A.G., Kokorina, A.L., Vaishlya, O.B., Agafonov, E.V, Guzhvin, S.A., Churakov, A.A., and Yakovleva, M.T. 2015. The scientific basis for the creation of new forms of microbial biochemicals. Sel'skokhozyaistvennaya biologiya 50(3):369-376. https://doi.org/10.15389/agrobiology.2015.3369eng Krusell, L., Madsen, L.H., Sato, S., Aubert, G., Genua, A., Szczyglowski, K., Duc, G., Kaneko, T., Tabata, S., and de Bruijn, F. 2002. Shoot control of root development and nodulation is mediated by a receptor-like kinase. Nature 420(6914):422-426. https://doi.org/10.1038/na-ture01207

Krusell, L., Sato, N., Fukuhara, I., Koch, B. E. V., Grossmann, C., Okamoto, S., Oka-Kira, E., Otsubo, Y., Aubert, G., and Na-kagawa, T. 2011. The Clavata2 genes of pea and Lotus japonicus affect autoregulation of nodulation. The Plant Journal 65(6):861-871. https://doi.org/10.1111/j.1365-313X.2010.04474.X Larrainzar, E., Riely, B.K., Kim, S.C., Carrasquilla-Garcia, N., Yu, H.-J., Hwang, H.-J., Oh, M., Kim, G.B., Surend-rarao, A.K., and Chasman, D. 2015. Deep sequencing of the Medicago truncatula root transcriptome reveals a massive and early interaction between nodulation factor and ethylene signals. Plant Physiology 169(1):233-265. https://doi.org/10.1104/pp.15.00350 Leppyanen, I.V, Shakhnazarova, V.Y., Shtark, O.Y., Vishnevs-kaya, N.A., Tikhonovich, I.A., and Dolgikh, E.A. 2018. Receptor-like kinase LYK9 in Pisum sativum L. is the CERK1-like receptor that controls both plant immunity and AM symbiosis development. International Journal of Molecular Sciences 19(1):8. https://doi.org/10.3390/ ijms19010008

Love, M. I., Huber, W., and Anders, S. 2014. Moderated estimation of fold change and dispersion for RNA-seq data with DESeq2. Genome Biology 15(12):550. https://doi. org/10.1186/s13059-014-0550-8 Lugtenberg, B., Rozen, D.E., and Kamilova, F. 2017. Wars between microbes on roots and fruits. F1000Research 6:343. https://doi.org/10.12688/f1000research.10696.1 Magori, S., Oka-Kira, E., Shibata, S., Umehara, Y., Kouchi, H., Hase, Y., Tanaka, A., Sato, S., Tabata, S., and Kawagu-chi, M. 2009. TOO MUCH LOVE, a root regulator associated with the long-distance control of nodulation in Lotus japonicus. Molecular Plant-Microbe Interactions 22(3):259-268. https://doi.org/10.1094/MPMI-22-3-0259 Morandi, D., Sagan, M., Prado-Vivant, E., and Duc, G. 2000. Influence of genes determining supernodulation on root colonization by the mycorrhizal fungus Glomus mosseae in Pisum sativum and Medicago truncatula mutants. Mycorrhiza 10(1):37-42. https://doi.org/10.1007/ s005720050285 Mortier, V., Den Herder, G., Whitford, R., Van de Velde, W., Rombauts, S., D'haeseleer, K., Holsters, M., and Goor-machtig, S. 2010. CLE peptides control Medicago truncatula nodulation locally and systemically. Plant Physiology 153(1):222-237. https://doi.org/10.1104/pp.110.153718 Müller, L. M., Flokova, K., Schnabel, E., Sun, X., Fei, Z., Frugoli, J., Bouwmeester, H.J., and Harrison, M.J. 2019. A CLE-SUNN module regulates strigolactone content and fungal colonization in arbuscular mycorrhiza. Nature Plants 5(9):933-939. https://doi.org/10.1038/s41477-019-0501-1 Nadal, M. and Paszkowski, U. 2013. Polyphony in the rhi-zosphere: presymbiotic communication in arbuscular mycorrhizal symbiosis. Current Opinion in Plant Biology 16(4):473-479. https://doi.org/10.1016/j. pbi.2013.06.005 Oldroyd, G.E.D. 2013. Speak, friend, and enter: signalling systems that promote beneficial symbiotic associations in plants. Nature Reviews Microbiology 1 1(4):252-263. https://doi.org/10.1038/nrmicro2990 Parniske, M. 2008. Arbuscular mycorrhiza: the mother of plant root endosymbioses. Nature Reviews Microbiology 6(10):763-775. https://doi.org/10.1038/nrmicro1987 Reid, D.E., Ferguson, B.J., Hayashi, S., Lin, Y.-H., and Gress-hoff, P.M. 2011. Molecular mechanisms controlling legume autoregulation of nodulation. Annals of Botany 108(5):789-795. https://doi.org/10.1093/aob/mcr205 Roy, S., Liu, W., Nandety, R. S., Crook, A., Mysore, K. S., Pislar-iu, C. I., Frugoli, J., Dickstein, R., and Udvardi, M. K. 2020. Celebrating 20 years of genetic discoveries in legume

nodulation and symbiotic nitrogen fixation. The Plant Cell 32(1):15-41. https://doi.org/10.1105/tpc.19.00279 Sagan, M. and Duc, G. 1996. Sym28 and Sym29, two new genes involved in regulation of nodulation in pea (Pisum sativum L.). Symbiosis 20:229-245. Sagan, M., Ney, B., and Duc, G. 1993. Plant symbiotic mutants as a tool to analyse nitrogen nutrition and yield relationship in field-growth peas (Pisum sativum L.). Plant and Soil 153(1):33-45. https://doi.org/10.1007/BF00010542 Salon, C., Munier-Jolain, N., Duc, G., Voisin, A.-S., Grandgi-rard, D., Larmure, A., Emery, R., and Ney, B. 2001. Grain legume seed filling in relation to nitrogen acquisition: a review and prospects with particular reference to pea. Agronomie 21(6-7):539-552. https://doi.org/10.1051/ agro:2001143

Shtark, O.Y., Puzanskiy, R.K., Avdeeva, G.S., Yurkov, A.P., Smolikova, G. N., Yemelyanov, V.V., Kliukova, M. S., Sha-varda, A.L., Kirpichnikova, A.A., Zhernakov, A.I., Afo-nin, A.M., Tikhonovich, I.A., Zhukov, V.A., and Shisho-va, M.F. 2019. Metabolic alterations in pea leaves during arbuscular mycorrhiza development. PeerJ 7:e7495. https://doi.org/10.7717/peerj.7495 Shtark, O.Y., Sulima, A.S., Zhernakov, A.I., Kliukova, M.S., Fedorina, J.V., Pinaev, A .G., Kryukov, A.A., Akhtemo-va, G.A., Tikhonovich, I.A., and Zhukov, V.A. 2016. Arbuscular mycorrhiza development in pea (Pisum sativum L.) mutants impaired in five early nodulation genes including putative orthologs of NSP1 and NSP2. Symbiosis 68(1-3):129-144. https://doi.org/10.1007/s13199-016-0382-2

Smith, S. E., Jakobsen, I., Granlund, M., and Smith, F. A. 2011. Roles of arbuscular mycorrhizas in plant phosphorus nutrition: interactions between pathways of phosphorus uptake in arbuscular mycorrhizal roots have important implications for understanding and manipulating plant phosphorus acquisition. Plant Physiology 156(3):1050-1057. https://doi.org/10.1104/pp.111.174581 Splivallo, R., Fischer, U., Gobel, C., Feussner, I., and Karlo-vsky, P. 2009. Truffles regulate plant root morphogenesis via the production of auxin and ethylene. Plant Physiology 150(4):2018-2029. https://doi.org/10.1104/ pp.109.141325

Strack, D. and Fester, T. 2006. Isoprenoid metabolism and plastid reorganization in arbuscular mycorrhizal roots.

New Phytologist 172(1):22-34. https://doi.org/10.1111/ j.1469-8137.2006.01837.x Tikhonovich, I.A., Andronov, E. E., Borisov, A.Y., Dolgikh, E. A., Zhernakov, A.I., Zhukov, V.A., Provorov, N.A., Roumi-antseva, M. L., and Simarov, B.V. 2015. The principle of genome complementarity in the enhancement of plant adaptive capacities. Russian Journal of Genetics 51(9):831-846. https://doi.org/10.1134/S1022795415090124 Tsyganov, V.E. and Tsyganova, A.V. 2020. Symbiotic regulatory genes controlling nodule development in Pisum sativum L. Plants 9(12):1741. https://doi.org/10.3390/ plants9121741

Tsyganova, A.V., Kitaeva, A.B., and Tsyganov, V.E. 2017. Cell differentiation in nitrogen-fixing nodules hosting sym-biosomes. Functional Plant Biology 45(2):47-57. https:// doi.org/10.1071/FP16377 Walter, M. H., Fester, T., and Strack, D. 2000. Arbuscular mycorrhizal fungi induce the non-mevalonate methyleryth-ritol phosphate pathway of isoprenoid biosynthesis correlated with accumulation of the 'yellow pigment'and other apocarotenoids. The Plant Journal 21(6):571-578. https://doi.org/10.1046/j.1365-313x.2000.00708.x Wang, C., Reid, J.B., and Foo, E. 2018. The art of self-con-trol-autoregulation of plant-microbe symbioses. Frontiers in Plant Science 9:988. https://doi.org/10.3389/ fpls.2018.00988 Wickham, H. 2016. Ggplot2: elegant graphics for data analysis. Springer. https://doi.org/10.1007/978-3-319-24277-4 Zhernakov, A.I., Shtark, O.Y., Kulaeva, O.A., Fedorina, J.V., Afonin, A.M., Kitaeva, A.B., Tsyganov, V.E., Afonso-Grunz, F., Hoffmeier, K., Rotter, B., Winter, P., Tikhonovich, I.A., and Zhukov, V.A. 2019. Mapping-by-sequencing using NGS-based 3'-MACE-Seq reveals a new mutant allele of the essential nodulation gene Sym33 (IPD3) in pea (Pisum sativum L.). PeerJ 7:e6662. https://doi. org/10.7717/peerj.6662 Zhukov, V.A., Shtark, O.Y., Nemankin, T.A., Kryukov, A.A., Borisov, A. Y., and Tikhonovich, I.A. 2016. Genetic mapping of pea (Pisum sativum L.) genes involved in symbiosis. Sel'skokhozyaistvennaya biologiya 51(5):593-601. https://doi.org/10.15389/agrobiology.2016.5.593eng Zipfel, C. and Oldroyd, G.E.D. 2017. Plant signalling in symbiosis and immunity. Nature 543(7645):328-336. https:// doi.org/10.1038/nature22009