# Identification of transcriptional modules linked to the drought response of *Brassica napus* during seed development and their mitigation by early biotic stress

Grégoire Bianchetti<fnref n="1" /> ID

Vanessa Clouet<fnref n="1" /> ID

Fabrice Legeai<fnref n="1" /> ID

Cécile Baron<fnref n="1" /> ID

Kévin Gazengel<fnref n="1" /> ID

Benoit Ly Vu<fnref n="2" />

Sébastien Baud<fnref n="3" /> ID

Alexandra To<fnref n="3" />

Maria J. Manzanares-Dauleux<fnref n="1" /> ID

Julia Buitink<fnref n="2" />† ID

Nathalie Nesi<fnref n="1" />† ID

### Correspondence

Nathalie Nesi,

Email: nathalie.nesi@inrae.fr

### Funding information

Institut National de Recherche pour l'Agriculture, l'Alimentation et l'Environnement, Grant/Award Number: INRAE-BAP-IB2017-SQUAL; Région Bretagne; Agence Nationale de la Recherche, Grant/Award Number: ANR-14-JFAC-007-01

Edited by R.M. Rivero

### Abstract

In order to capture the drought impacts on seed quality acquisition in *Brassica napus* and its potential interaction with early biotic stress, seeds of the ‘Express’ genotype of oilseed rape were characterized from late embryogenesis to full maturity from plants submitted to reduced watering (WS) with or without pre-occurring inoculation by the telluric pathogen *Plasmodiophora brassicae* (Pb + WS or Pb, respectively), and compared to control conditions (C). Drought as a single constraint led to significantly lower accumulation of lipids, higher protein content and reduced longevity of the WS-treated seeds. In contrast, when water shortage was preceded by clubroot infection, these phenotypic differences were completely abolished despite the upregulation of the drought sensor RD20. A weighted gene co-expression network of seed development in oilseed rape was generated using 72 transcriptomes from developing seeds from the four treatments and identified 33 modules. Module 29 was highly enriched in heat shock proteins and chaperones that showed a stronger upregulation in Pb + WS compared to the WS condition, pointing to a possible priming effect by the early *P. brassicae* infection on seed quality acquisition. Module 13 was enriched with genes encoding 12S and 2S seed storage proteins, with the latter being strongly upregulated under WS conditions. Cis-element promotor enrichment identified *PEI1/TZF6*, *FUS3* and *bZIP68* as putative regulators significantly upregulated upon WS compared to Pb + WS. Our results provide a temporal co-expression atlas of seed development in oilseed rape and will serve as a resource to characterize the plant response towards combinations of biotic and abiotic stresses.

<figure><img src="image_1.png" /></figure>

<footnote n="1">¹IGEPP, INRAE, Institut Agro Rennes-Angers, Université de Rennes, Le Rheu, France</footnote>

<footnote n="2">²IRHS, INRAE, Institut Agro Rennes-Angers, Université d'Angers, France</footnote>

<footnote n="3">³JPB, INRAE, AgroParisTech, Versailles, France</footnote>

<footnote>† These authors contributed equally to this work</footnote>
# 1 | INTRODUCTION

With the world's human population constantly increasing and expected to reach 9.7 billion in 2050 (UN projections), seeds are even more central levers regarding food security and the conservation of genetic resources. However, the seed market can be easily and unpredictably compromised by challenges facing the production, including the need to reduce chemical inputs, preserve natural resources and cope with environmental fluctuations. Therefore, to face the challenges of adapting agriculture to climate change, seed production should have increased resilience to abiotic stress factors in the context of an expected proliferation of pathogen infections. This issue concerns both seed nutritional quality and vigour, two crucial factors for seedling establishment and yield.

Seed qualities are acquired during its development that is divided into three phases well described for a number of plant species (Goldberg et al., 1994; Baud et al., 2002; Le et al., 2010; Bewley et al., 2013): (1) embryogenesis that leads to the formation of the embryo, the endosperm and the integuments, (2) maturation, characterized by a massive biosynthesis and accumulation of storage reserves, when germination ability and desiccation tolerance are acquired, and (3) late maturation, when longevity is acquired and chlorophyll is lost. The duration of these phases greatly varies between plant species (Leprince et al., 2017) and is modulated by the growing conditions (Righetti et al., 2015). Assimilates provided by the maternal tissues fuel the seed development process that is governed by the interplay of tight hormonal control (reviewed in Locascio et al., 2014) and a complex molecular network (Santos-Mendoza et al., 2008). In particular, the so-called LAFL network consisting of master regulatory genes (LEC1, LEC2, ABI3, FUS3, LEC1-LIKE) orchestrates many aspects of the seed maturation program including the accumulation of storage compounds as well as the acquisition of desiccation tolerance together with the induction of dormancy (Kotak et al., 2007; Santos-Mendoza et al., 2008; Verdier and Thompson, 2008; Fatihi et al., 2016).

How the environment influences the genetic and physiological determinisms of seed qualities remains to be elucidated. Among the numerous environmental constraints, drought stress is and will be one of the most common limitations encountered by crops and other plants in current and future climate scenarios, as it results in significant productivity losses (Li et al., 2009). Furthermore, drought stress during seed development is generally known to lower crop seed qualities by reducing the time of the seed developmental phases (Barnabas et al., 2008; Stagnari et al., 2016). Nevertheless, drought effects vary with plant species, timing, intensity and duration of the applied constraint. In addition, in nature, plants are simultaneously or sequentially exposed to a variety of stresses that limit crop productivity.

Numerous studies now support that stress combinations (i.e. stresses occurring concurrently with some degree of overlapping) trigger complex plant responses at the transcriptome, cellular and physiological levels that significantly differ from the responses to individual stresses (Rasmussen et al., 2013), leading to the conclusion that they must be considered as entirely new constraints (Mittler, 2006). For instance, the combination of heat and drought often results in much

greater damage to crops than either stress individually (Mittler, 2006; Barnabas et al., 2008). When examining combinations of biotic and abiotic stresses, the interactions are even more complex to predict, with either positive or negative effects of one stress to the other, depending on the nature and modalities of application of each single stress (reviewed in Atkinson and Urwin, 2012; Suzuki et al., 2014; Pandey et al., 2017). Especially, 'drought-pathogen' stress interactions are probably among the most common and important stress combinations in agriculture with drought either reducing or stimulating plant susceptibility to pathogens (reviewed in Pandey et al., 2017). On the other hand, pathogen infection can actually confer advantages to plants exposed to drought, which can be explained by the pathogen-mediated accumulation of specialized metabolites such as osmoprotectants in plant tissues (Xu et al., 2008). Within the large body of literature regarding plant stress responses, information on the response of reproductive tissues to stress combinations, especially those that lead to positive effects on crop yield and quality, is rather scarce (Creus et al., 2004).

With a world production of 71 Mt of seeds per year, oilseed rape (OSR, Brassica napus L.) is the third largest oleaginous crop in the world. In addition to producing oil and proteins for food, feed and green chemistry, OSR also fulfils various ecosystem services (e.g., diversification of crop systems, improvement of soil quality, limitation of weeds and pathogens, promotion of insect biodiversity). Seed oil (45-50% w/w) and proteins (20-25% w/w) are primarily accumulated in the cotyledon cells during seed development (Barthet et al., 2011; Borisjuk et al., 2013). Average seed yield of OSR crops ranges from 1.5-2 t ha⁻¹ in Canada and Eastern Europe (extensive production of spring varieties) to 3-3.5 t/ha in Western Europe (intensive production of winter-type cultivars) (Nesi et al., 2008). However, this production is more and more threatened by a number of stresses that occur as single, multiple individuals or combinations throughout the growing season (Brunel-Muguet et al., 2015; Elferjani and Soolanayakanahally, 2018). In OSR, drought severely reduces seed yield but also impacts seed qualities by reducing oil and increasing proteins in seeds while modulating seedling vigour traits (Bouchereau et al., 1996; Champolivier and Merrien, 1996; Hatzig et al., 2018). In addition, clubroot disease, caused by the oblige telluric biotroph protist *Plasmodiophora brassicae*, leads to severe yield losses for Brassica crops worldwide (Dixon, 2009). Clubroot induces swellings or galls on the roots of infected plants that decrease the flow of water and nutrients within the plant throughout the growth cycle, thus mimicking a water shortage. Infected plants have been shown to develop pathogen tolerance by using specialized cellular programs (PCD) in addition to the biosynthesis of osmoprotectants or protective proteins such as the HEAT SHOCK PROTEINS (HSP; Summanwar et al., 2021). In addition, the pathogen partly hijacks the plant metabolism for its own benefit (Ludwig-Müller et al., 2009). Depending on the plant genotype, *P. brassicae* infection can lead to reduced straw yield, oil content of the grain, and grain mass (Pageau et al., 2006). Regarding stress combinations, previous studies highlighted that OSR response to clubroot can be modified by nutrient availability (Laperche et al., 2017;
Aigu et al., 2018, 2022). Likewise, in Arabidopsis, the expression of genetic and epigenetic factors controlling clubroot resistance was dependent on the level of irrigation (moderate drought and flooding) as well as the temperature (Gravot et al., 2016; Liégard et al., 2019). However, it is not known whether and how early *P. brassicae* infection interacts with OSR drought response during the reproductive phase.

This prompted us to investigate the molecular bases of single and combined stresses on seed qualities in OSR. Our goal was to address (1) whether a pre-occurring biotic stress such as the *P. brassicae* infection can pre-dispose the plants' acquisition of seed quality towards an upcoming abiotic stress and (2) what are the underlying molecular pathways modulated by the single and combined stresses. To this purpose, we chose to combine infection by *P. brassicae* at the seedling stage with drought stress applied at the onset of flowering, given that these two constraints are critical for OSR crops under field conditions. We designed an experiment to capture the dynamics of 72 transcriptomes of developing OSR seeds grown under either drought conditions, the presence of *P. brassicae* alone, or a combination of both stresses. This dataset was first used to generate a robust regulatory coexpression network of seed development in OSR using the weighted gene co-expression network analysis (WGCNA; Zhang and Horvath, 2005) based on the so-called "guilt-by-association" approach (Wolfe et al., 2005). The integration of seed physiological traits acquired under the different treatments, as well as information on functions known to regulate seed quality traits, helped to pinpoint gene modules and candidate genes as potential key regulators of stress resistance mechanisms in seeds. Our data will contribute to a better understanding of the crosstalk between the plant responses to biotic-abiotic combinations and how it modulates seed quality.

## 2 | MATERIALS AND METHODS

### 2.1 | Plant material, growth conditions and physiological analyses

Plants of winter oilseed rape (cv. 'Express') were cultivated over the growing season 2017-2018 in an open, covered but unheated tunnel, thereby subjecting them to local ambient temperatures and undergoing vernalization comparable to field conditions, as described in Bianchetti et al. (2021a). Each biological replicate consisted of a tank of ~1 m³ with 42 plants regularly spaced and individually irrigated. The water potential of the substrate was followed in each tank using a ceramic tensiometer (Tensionic®, SDEC), as reported by Bianchetti et al. (2021b). Plants were grown under four different treatments, as presented previously (Bianchetti et al., 2021a, 2021b). The "Pb" treatment involved inoculating 7-day-old seedlings with 1 mL of a 10⁶ resting spore suspension of *P. brassicae* (eH isolate; P1 pathotype), according to the method described by Manzanares-Dauleux et al. (2000a). The "WS" condition for which the water potential of the substrate was lowered and maintained around -600 mbar from the onset of flowering (i.e., BBCH60 developmental stage, according to

Lancashire et al., 1991). The "Pb + WS" condition combined clubroot inoculation and water shortage from BBCH60 onwards. Finally, the control "C" tanks where the water potential of the substrate was maintained above -200 mbar throughout the growth cycle (with nutritive solution replaced by clear water after application of the water shortage in the WS treatment). The full experimental design included eight independent biological replicates for the C and WS treatments and four for the Pb and Pb + WS treatments (Table S1). Clubroot symptoms were recorded according to Manzanares-Dauleux et al. (2000b) at 8 or 30 weeks after *P. brassicae* inoculation and at BBCH60 (2 plants per replicate at each timepoint), as well as at the end of the crop cycle (30 plants). Pods on the main raceme were tagged at BBCH60, and immature seeds were harvested from 3-5 of the 20 central plants of each tank at different thermal times intervals, as follows (in GDD, growing degree-days calculated from BBCH60 with a base temperature of 0°C): 280 GDD (T1), 400 (T2), 520 (T3), 675 (T4), 785 (T5), 860 (T6), 985 (T8), and 1130 (T9) (Bianchetti et al., 2021a). Mature, dry seeds were harvested when the vegetative parts were fully senescent, and the seeds were completely dark and hard (BBCH90). Water contents were assessed gravimetrically for triplicate samples of five seeds by determination of the fresh weight and subsequent dry weight (DW) after drying at 105°C overnight. Water contents are expressed on a dry weight basis.

### 2.2 | Fatty acid and protein contents of seeds

For determination of the seed fatty acid content and composition, 6 mg of crushed lyophilized material were incubated in 3 mL of methanol/toluene/sulphuric acid (100/2.5/30; v/v/v) at 94°C for 90 min after addition of C17:0 fatty acid as an internal standard. Fatty acyl methyl esters were then extracted into 1 mL of hexane following the addition of 3.39 mL of water. After vigorous shaking, the upper organic phase was analyzed by gas chromatography, as described by Li et al. (2006). For determination of the seed protein content, 2 mg of crushed lyophilized material was grounded in 200 µL of extraction buffer (Hepes 50 mM pH 7.5, SDS 3% [w/v], glycerol 10% [v/v], NaCl 150 mM, DTT 20 mM, protease inhibitor cocktail [Roche]) using a bead beating grinder (3 × 45 s at 6.5 m s⁻¹). The mixture was then centrifuged for 40 min at 15000 g at room temperature. Quantification of the total protein content was performed on the supernatant using the 2-D Quant kit (GE Healthcare Europe GmbH) according to the manufacturer's instructions.

### 2.3 | Seed physiology

Seeds were harvested at different developmental stages and rapidly dried at 42% relative humidity (RH) under an airflow. Percentage of final germination, measured as protrusion of the radicle through the surrounding seed layers, was determined on dried seeds by imbibition at 20°C in the dark for 7 days. For accelerated ageing assays, seeds were placed over a saturated solution of NaCl (75%
RH) at 35°C in hermetically sealed equilibration boxes and viability was evaluated by germination at regular time intervals on 60 seeds for each timepoint. Longevity was expressed as P50, defined as the time (days) at which the stored seed lost 50% viability during storage, derived from survival curves of the percentage of final germination over storage time.

## 2.4 | RNA sequencing, data processing and read mapping

RNA extraction and sequencing (RNAseq) were performed on seed samples collected at all developmental stages for the four experimental conditions. However, due to a limitation in the number of available samples at stages T4 and T6, RNAseq was conducted only for the control and water shortage (WS) conditions. RNA extraction, sequencing as well as read mapping are described in detail by Bianchetti et al. (2021a). Briefly, the cDNA libraries for RNA-seq were prepared using the QuantSeq 3' mRNA-Seq Library Prep Kit REV (Lexogen). In total, 72 libraries were sequenced using a paired-end read protocol on the Illumina Hiseq™4000 sequencing platform at the Beijing Genomics Institute (BGI, Hong Kong; https://www.bgi.com). The raw sequence data were deposited in the NCBI Sequence Read Archive (Accession No. PRJNA738318). Processing of the RNA-seq data was carried out using the adapted nf-core RNA-seq pipeline version v1.4.2 (Ewels et al., 2020; https://github.com/nf-core/rnaseq) and the clean reads were mapped onto the reference genome Darmor-bzh-v4 (Chalhoub et al., 2014) with Hisat2 v2.1.0 (Kim et al., 2019). Transcript quantification was determined using the feature-Counts software v1.6.4 (Liao et al., 2014) in relationships to the Darmor-bzh-v5 gene annotation (Chalhoub et al., 2014). Read counts were converted to counts per million (CPM) values (cpm in edgeR) and CPMs were calculated by normalizing the read counts by the total counts per sample. Based on previous characterization of the whole transcriptomic dataset (refer to Bianchetti et al., 2021a), uniquely mapped read pairs assigned to an annotated gene of the reference genome represented a total of 75177 transcripts (≈ 64% of the clean reads) of which the 42467 genes with ≥1 CPM in all of the biological replicates of at least one treatment per seed stage were retained for analyses (refer to Table S1 from Bianchetti et al., 2021a).

## 2.5 | Construction of a gene co-expression network and quantification of gene expression-trait associations

Weighted gene co-expression network analysis (Zhang and Horvath, 2005) was performed using the WGCNA package (v1.68) in R (Langfelder and Horvath, 2008; https://github.com/cran/WGCNA). From the initial set of 42,467 genes with a minimum expression level of 1 CPM, two filters were applied for data refinement. Initially, the top 25% of genes with the highest standard deviation in transcript levels across biological replicates were eliminated to reduce noise from inconsistently expressed genes. Subsequently, to focus on genes exhibiting marked expression variations across different

stages and/or treatments, the top 50% of genes, based on average transcript levels computed across all 72 samples, were retained. Besides, we ran the analysis of network topology to select the most optimal soft threshold that balanced scale independence and mean connectivity of the WGCNA (power 22 was chosen for network construction, as shown in Figure S4A). Gene expression values of filtered genes were imported into the package, and the modules were determined through automatic one-step network construction with the following parameters: power (soft threshold) = 22, minModuleSize = 10, maxBlockSize = 16,000, mergeCutHeight = 0.05, TOMType = 'unsigned'. The module eigengene (ME) of a given module, defined as the first principal component of the module standardized expression profiles (Langfelder and Horvath, 2008), was considered the best representative gene expression of a given module. Finally, the modules were clustered into meta-modules according to the ME distances using an agglomerative nesting hierarchical clustering analysis (agnes) with Ward's criterion. The Euclidean distances between modules were exported, and the eigengene network was visualized using the ggplot2 (v3.3.3) package in R (Wickham, 2009). Module-trait associations were searched using Pearson correlation coefficients (PCC) between ME and physiological trait values. Hub genes for a trait of interest were identified based on both high gene significance (GS, i.e. PCC between individual gene expression level and trait value) and high module membership (MM, i.e. PCC between individual gene expression level and ME). All correlations were considered statistically significant for FDR-adjusted p-value <0.05.

## 2.6 | Identification and analysis of differentially expressed genes, bioinformatics

Transcriptomic dataset initial global characterization and variance partitioning, in relation to the primary factors of our study, were respectively conducted using both Principal Component Analysis (PCA) with FactoMineR R package (v2.4; Lê et al., 2008) and Redundancy Analysis (RDA) with vegan R package (v2.5-7; Oksanen et al., 2015). For gene expression analyses, the AskoR R library was used (Alves Carvalho et al., unpub.; https://github.com/askomics/askoR). Analyses of differentially expressed genes (DEG) were conducted using a generalized linear model with likelihood ratio tests (glmLRT in edgeR), and the Trimmed Mean of M-values (TMM) normalization method was applied. Genes exhibiting an FDR-adjusted p-value less than 0.05, following the Benjamini-Hochberg correction, were identified as differentially expressed in the comparisons. These comparisons were made for each seed developmental stage and for each specific stress condition, contrasting seeds produced under controlled conditions with those produced under stress conditions. The gene ontology (GO) enrichment analyses of selected genes compared to the reference genome (GO terms of the Darmor-bzh-v5 gene annotation; Chalhoub et al., 2014) were run with topGO R package (v2.40.0; Alexa et al., 2006) using default parameters with a Fisher's exact test (p-value <0.01). For analysis of cis-element enrichment in promoter sequences of target genes, the SeqEnrich software was used (Becker et al., 2017; https://github.com/nagreme/SeqEnrich) that gathers information on B. napus
transcription factors, promoter motifs and GO to produce predictive regulatory networks. Analyses were run using default parameters with a hypergeometric test, and DNA motifs with a *p*-value <0.05 were considered as enriched within the query list.

## 2.7 | RT-qPCR analysis

For reverse transcription, 1 µg of DNA-free RNA extract was converted to first-strand cDNA using the M-MLV Reverse Transcriptase (Promega Corp.) and oligo(dT)₁₂–₁₈ primers according to the manufacturer's instructions. For quantification of transcript levels, samples of cDNA (10 ng) were subjected to quantitative real-time PCR (qRT-PCR) using the SYBR Green I Master (Roche) on a LightCycler® 480 instrument (Roche). The following thermal profile was used: 95°C for 5 min, 45 cycles of 95°C for 15 s, and 60°C for 40s. Amplicon dissociation curves were recorded after 45 cycles (95°C for 30s, 60°C for 30s and 97°C with a ramp rate of 0.11°C/s) to check for PCR product integrity. Genes *BnaC.RibS3* (*BnaC04g13040D*) and *BnaC.UBQ11* (*BnaC04g09440D*) were used as reference genes in our conditions to normalize target gene expression levels through the geNorm method (Vandesompele et al., 2002). The relative gene expression levels were determined using the Efficiency Method described by Tellman (2006). Primers were designed with Primer3web version 4.1.0 and are listed in Table S2. Values are based on three biological replicates.

## 2.8 | Accession numbers

Gene IDs are those from AGI or GenBank/EMBL databases for Arabidopsis and from the Darmor-bzh-v5 reference genome (Chalhoub et al., 2014; https://www.genoscope.cns.fr/brassicanapus/) for B. napus. The correspondence between the genes studied in the paper for Darmor-bzh-V5 and Darmor-bzh-V10 (Rousseau-Gueutin et al., 2020) is available in Tables S2 and S3.

## 2.9 | Statistical analyses of plant and seed traits

The experimental design is presented in Table S1. Data were expressed as the mean ± standard error (SE). All statistical analyses were conducted using R software (R Development Core Team; https://www.r-project.org/). Comparisons between two groups (e.g., stress vs. control at a specific sampling time) were performed using Student's t-test, with *p*-values below 0.05 considered statistically significant.

# 3 | RESULTS

## 3.1 | Efficient application of multiple and combined stresses during *B. napus* seed development

In order to capture and analyse the impact of drought on seed development in *B. napus* and its potential interaction with early biotic stress

applied from seedling growth onwards, seeds of the 'Express' genotype of winter OSR were characterized from late embryogenesis to full maturity (Figure 1A). Mother plants were submitted to reduced watering (WS) during the whole reproductive phase with or without pre-occurring inoculation by the pathogen *P. brassicae* (Pb + WS or Pb respectively) at the seedling stage, and compared to the control conditions (C, normal watering and no Pb inoculation).

Stress applications were monitored regularly throughout the plant cycle, and stress-associated indicators were scored at critical time points to evaluate the efficiency of the treatments. Soil water potential rapidly and significantly decreased after the start of water shortage (Figure 1B) and was maintained at a level of -600 mbar to imitate a moderate drought. Transcript accumulation of a marker gene for plant drought response, RD20, was evaluated after one month of water deficit. RD20 was shown to be upregulated in mature leaves of plants submitted to a reduced watering regime compared to their well-watered counterparts, irrespective of the previous inoculation with clubroot (Figure 1C). Recording clubroot symptoms on Pb and Pb + WS plants at the end of the growth cycle showed about 40% of plants with clear root galls, with no significant effect of the watering regime (Figure 1D). Altogether, these results demonstrated that biotic and abiotic stress were efficiently applied and sensed by the plants in our experimental design.

## 3.2 | Seed maturation is mostly affected by drought when applied as a single-stress

The dynamics of seed moisture content and dry weight (DW) throughout the reproductive phase combined with morphological inspections of dissected seeds revealed the three critical stages of seed development (Figure 1A): mid-to-end of morphogenesis between T1-T2 stages (280–400 GDD respectively) with a slow increase in seed DW, a strong decrease in moisture and embryos from torpedo to U shape stage; seed filling ranging from T2-T6 stages (400–860 GDD respectively) with a drastic increase of seed DW, a slower decrease in water content than the previous phase and embryos from stick cotyledons to fully mature stage; and finally the late maturation phase between T6-T9 (860–1130 GDD respectively) stages with a stable seed DW at ≈3.5 mg/seed, a continuing decrease in water content, and the progression of brown pigmentation over the seed integuments concomitant to the loss of chlorophyll in the embryo. No statistical differences for these physiological and morphological parameters were found between the different stresses (WS, Pb, and Pb + WS) compared to the control conditions (C) (Figure 1A).

To get further insight into the impact of stresses on the acquisition of seed quality, fatty acid (FA) and protein accumulation, as well as seed longevity, were evaluated throughout seed development (Figure 2). Total FAs mainly accumulated from T1 to T6 stages with the strongest increase between T1 and T3 stages (1–40% of seed DW) to reach a maximum of ≈ 45% of seed DW (Figure 2A). A slight decrease could be observed at the end of seed development (T9 stage) consistent with previous reports (Baud et al., 2002; Chia et al., 2005). During seed filling (T1-T6), all the stress conditions
<figure><img src="image_2.png" /><figcaption>**FIGURE 1** Seed development of *Brassica napus* (cv. 'Express') under various stress conditions. A. Evolution of the seed water content (dotted line) and dry weight (DW, solid line) during seed development. Immature seeds were harvested at different timepoints (T) as depicted on pictures (bar 1 mm for all). Thermal times are given in growing degree-days (GDD). Data are expressed as the mean ± SE from three biological replicates of five seeds for each treatment. B. Mean values (solid line) with standard deviations (dashed lines) of soil water potentials (mbar) in eight (C and WS) or four (Pb and Pb + WS) replicated tanks over a 20-day window that frames the application of the drought stress. C. Relative expression level of RD20 (qRT-PCR) in mature leaves of *B. napus* plants grown under the different conditions (three biological rep. each). D. Number of plants with clubroot symptoms evaluated at the onset of WS application (4 rep. × 2 plants) or at the end of the crop cycle (4 rep. × 30 plants). C, control; WS, water shortage; Pb, *P. brassicae* inoculation; Pb + WS, *P. brassicae* inoculation and water shortage; rep., biological replicate.</figcaption></figure>

induced a significant over-accumulation of total lipids, at least for one
stage, compared to seeds produced under control conditions
(Figure 2A). However, in mature seeds (T8-T9 seed stages), the FA
content was no longer significantly different between treated vs. con-
trol seeds, with the exception for seeds produced under water short-
age, for which a significant decrease was observed (−5.2%;
Figure 2A). Similar modulations were observed in the oleic acid
(C18:1) content, which is the primary constituent of total FAs (70%)
(Figure 2B). As with total FAs, over-accumulation of linoleic (C18:2)
and α-linolenic (C18:3) acid contents were measured during seed fill-
ing for at least one stress condition (Figure 2C-D). However, in mature
seeds, the linoleic acid content produced under all stress conditions
was consistent with that of control seeds, while a lower content of
alpha-linolenic acid was measured in both Pb + WS and WS treat-
ments. When applied as a single constraint, WS induced a significant
increase in seed total proteins throughout seed development (+28%
of the control). In contrast, total protein contents were almost not
modulated by Pb and Pb + WS compared to the control condition.
Indeed, slightly lower protein content was only measured in Pb seeds
at T5 (Figure 2E). Seed longevity, scored through survival curves dur-
ing artificial aging, decreased only in mature seeds (T9) produced
under drought conditions (Figures 2F, S1). These data demonstrate

that the seed maturation program is significantly modified by drought
when applied as a single constraint, resulting in significantly lower
accumulation of lipids, higher content in proteins and reduced longev-
ity of the WS-treated seeds compared to the control conditions. In
contrast, when water shortage was preceded by clubroot infection,
these differences were completely abolished, suggesting the pre-
occurring biotic stress pre-disposed the plants towards the upcoming
abiotic stress.

3.3 | Time is the main driver of gene expression throughout seed development

In parallel to the phenotypic characterization, a molecular characteri-
zation of the transcriptional pathways during seed development was
carried out. Immature seeds were collected at eight seed stages in
three replicates for C and WS conditions (T1-T6, T8, T9) and at six
seed stages in two replicates for Pb and Pb + WS conditions (T1-T3,
T5, T8-T9) (Bianchetti et al., 2021a). A global overview of the tran-
scriptome modulations during seed development is shown by principal
component analysis (PCA) (Figures 3A, S2). The first and second axes
of the PCA captured 44.2 and 14.4% of the total transcriptome
**FIGURE 2** Modulation of quality in *B. napus* seeds produced under stress conditions. Evolution of fatty acids (A-D) or protein (E) contents in developing seeds of oilseed rape. For each sampling date, data are means of three independent replicates for C and WS and two for Pb and Pb + WS, except for T1 where only one replicate was analyzed;

bars = SEs. Statistical levels of significance (treated vs control): none, non-significant; \*p* < 0.05; \*\*p* < 0.01; \*\*\*p* < 0.001.

Longevity (F) was expressed by viability (germination) after 40 days of artificial aging at 75% RH and 35°C (solid lines) compared to germination before aging (dotted lines) (three rep. for C and WS; two rep. for Pb and Pb + WS). C, control; WS, water shortage; Pb, *P. brassicae* inoculation; Pb + WS, *P. brassicae* inoculation and water shortage; rep., biological replicate.

<figure><img src="image_3.png" /></figure>

variation, respectively and were primarily related to the kinetics of seed development. Between stages T2 and T5, the WS samples were advanced along the developmental gradient compared to the control samples, whereas in the case of Pb + WS, this advance in development was only apparent for stage T2. Variance partitioning analysis showed that 69.6% of the variance of the transcript levels was explained by the seed developmental stage, while 10% was significantly attributed to the interaction between seed stage and treatment (Figure 3B), revealing that changes in the transcriptomes under our experimental conditions were mainly driven by developmental mechanisms with discrete interactions with treatments, and in particular drought.

The timing of molecular events during OSR seed development was further investigated through the expression of orthologs of master regulatory genes forming the LAFL seed network and related genes (Figure 3C-G). For each Arabidopsis target gene, several orthologs were identified in the polyploid genome of OSR, ranging from two to four copies. Most of the orthologs showed similar transcription patterns or were not expressed (Figure S3). The transcript level of

*BnLEC1* (BnaA07g10770D) was high during early embryogenesis and rapidly decreased from T1 to T3 to become barely detectable during seed filling and maturation as described for Arabidopsis seed (Meinke et al., 1994) (Figure 3C). This decrease was initiated earlier during development for seeds produced under the different stressful conditions compared to the control condition. Transcript level of *BnLEC2* (BnaA07g08500D) was lowly abundant, and both orthologs were very early expressed, concomitant with *BnLEC1* (Figure 3D). Transcript levels of *BnFUS3* (BnaC07g48930D, Figure 3E) displayed a bell shape encompassing the embryogenesis and seed filling phases, with an increase from T1, a maximum level at T2-T4 and a sharp drop afterwards. At T5, the expression level of *BnFUS3* was lower for the seeds produced under stress conditions, with the differences being only significant for WS and Pb + WS seeds. The *BnABI3* gene (BnaA03g37510D, Figure 3F), a key regulator of maturation related to ABA signalling (Nambara et al., 1992), was expressed during seed maturation between T2 to T8 and only decreased upon final maturation drying. No significant differences were detected for treated seeds compared to
<figure><img src="image_4.png" /><figcaption>FIGURE 3 Multivariate analysis of the RNA-seq dataset and gene profiles of main regulators of seed development in *B. napus*. The transcript levels of the 42467 genes with ≥1 CPM (count per million) were used for sample discrimination by PCA (A) and variance partitioning analysis using a redundancy analysis (RDA) model (B). The percentage of variance associated to each factor of the RDA model is indicated with statistical levels of significance as following: none, non-significant; \*p < 0.05; \*\*p < 0.01; \*\*\*p < 0.001; NA, not available. C-G Expression profiles of key transcription factors: LEAFY COTYLEDON 1 (BnLEC1 BnaA07g10770D) (C), LEAFY COTYLEDON 2 (BnLEC2, BnaA07g08500D) (D), FUSCA3 (BnFUS3, BnaC07g48930D) (E), ABSCISIC ACID INSENSITIVE3 (BnABI3, BnaA03g37510D) (F), and ABSCISIC ACID INSENSITIVE5 (BnABI5, BnaC04g09030D) (G). Means of gene expression levels are given in CPM with respective standard errors. Statistical levels of significance (normalized log2 transcript levels of treated vs control samples): none, non-significant; \*p < 0.05; \*\*p < 0.01; \*\*\*p < 0.001. Gene ID are given according to the reference gene annotation (Darmor-bzh-v5). C, control; WS, water shortage; Pb, *P. brassicae* inoculation; Pb + WS, *P. brassicae* inoculation and water shortage; GDD, growing degree-days. C, control; WS, water shortage; Pb, *P. brassicae* inoculation; Pb + WS, *P. brassicae* inoculation and water shortage; T1-T9, timepoints for seed sampling (refer to Figure 1).</figcaption></figure>

the control. A second gene that is also involved in the ABA signalling pathway but preferentially expressed during late maturation is *ABI5* (Bensmihen et al., 2002). In OSR seeds, *BnABI5* (BnaC04g09030D, Figure 3G) was barely expressed before T4. Then, transcript levels increased, reached a maximum at T6 and remained stable until the end of seed development. A significant over-expression of *BnABI5* was shown during seed filling in WS and Pb + WS seed samples. Transcription of the LAFL genes in OSR seeds was consistent with the temporal seed development observed in other plant species (Santos-Mendoza et al., 2008; Bizouerne et al., 2021). In addition, we noticed that drought significantly modified gene expression profiles, with a tendency to accelerate the shift between phases.

## 3.4 | Construction and characterization of a weighted gene co-expression network of seed development in *B. napus*

To characterize the temporal expression of the seed transcriptomes, genes with similar expression profiles were grouped in gene modules via the construction of a gene co-expression network by WGCNA on a subset of 15660 genes and 72 transcriptomes from the four developmental series studied (Table S3). A total of 1747 genes were removed from further analyses because they were weakly interconnected and poorly correlated with the modules. The remaining 13913 genes were clustered into 33 co-expression modules, with the size of
<figure><img src="image_5.png" /></figure>

<figure><img src="image_6.png" /><figcaption>**FIGURE 4** Legend on next page.</figcaption></figure>
modules ranging from 11 (module 33) to 3608 (module 01) genes (Figures 4A; S4B). For each sample, the gene expression of a given module was summarized by the module eigengene (ME) value. In our study, variations of ME between samples were mostly related to seed development. All gene modules contained genes that were expressed under all conditions, although clear differences in ME values were apparent between treatments (Figure 4A). All features of genes and samples from the WGCNA are provided in Tables S3 and S4, respectively.

An eigengene network was constructed by clustering ME distances to reduce the complexity of the network. As such, the modules were regrouped into 10 meta-modules, G1 to G10 (Figures 4A, B). The organization of the meta-modules followed the developmental phases described in Figure 1A, with genes of the meta-groups G1 (modules 09 and 07), G2 (module 15) and G3 (modules 05, 11, 02, 03 and 01) primarily expressed during embryogenesis (T1-T2), genes of G4 (modules 08, 10, 23 and 14) and G5 (modules 27, 04, 19 and 21) expressed during embryogenesis and early maturation (T1-T5), genes of G6 (modules 33, 18 and 30), G7 (modules 28, 13, 20, 17, 25 and 22), G8 (modules 32, 31, 16 and 24), G9 (modules 29 and 26) expressed during maturation (T3-T6), and genes of G10 (module 06 and 12) expressed during late maturation (T8-T9). Enrichment analysis of the biological processes in the meta-modules supported that functions related to translation, cell machinery assembly, photosynthesis or sucrose transport were primarily associated with embryogenesis (G1, G2 and G3) and concerned ≈73% (10094) of the genes from the WGCNA (Figure 4B; Table S5). In addition, biological processes related to metabolite biosynthesis (fatty acids, amino acids), accumulation of storage compounds (lipids, proteins) or oxidation-reduction were found during the seed filling and early maturation stages (G4, G5, G6 and G7) and involved ≈20% (2859) of the genes. Finally, mid-to-late seed maturation modules (G8, G9 and G10) targeted 7% (960) of the genes and were characterized by biological processes linked with reorganization at the cellular and molecular levels (e.g., cell wall macromolecule catabolism, chitin catabolism, RNA process and modification, protein folding, DNA topological change). Furthermore, 35% (23 of 66) genes from the G9 meta-groups encoded HEAT SHOCK PROTEINS (HSP) (Table S3), thus supporting a role for preparation to seed dehydration.

## 3.5 | Modulation of the seed development gene co-expression network by stresses

Although seed development was the main driver of gene expression in our conditions, the ME profiles indicated gene expression differences

between treatments. To mitigate potential biases, the global analysis of differentially expressed genes (DEG) (stresses vs control) targeted stages common to all samples, thereby excluding DEG unique to T4 and T6 for WS seeds. The analysis identified that 48% of the network genes (6,696 out of 13,913; Tables S6 and S7) were significantly deregulated throughout seed development in response to the stresses (Figure 5A). Most DEG (85%, 5682 genes) were due to WS whereas Pb led to the least amount of DEG, only 24% (1593). A total of 1231 DEGs were found for all three stresses. A comparison of WS with Pb + WS treatments showed that 2172 DEG were in common between both stresses (although the extent of differential expression differed between both conditions), and that 2101 genes were specific for WS, and 830 were specific for Pb + WS. To focus on pronounced DEG, a log2 fold change threshold of 1.5 was established. DEGs that meet this criterion were then mapped onto the seed development network, with their proportion per module indicated using a colour code ranging from 25 to 100% (Figure 5B). This analysis revealed that drought stress (WS) led to the highest degree of gene module deregulation, with 22 distinct modules being predominantly affected under WS conditions alone. In contrast, the Pb condition impacted only 3 modules significantly, while the combined Pb + WS condition influenced 14 modules (Figure 5B; Table S7). Four modules showed a strong reduction (>15%) in DEG in Pb + WS compared to WS (M08, M11, M23 and M30; Table S7). Strikingly, amongst the modules with more than 25% DEG, module 29 was the only one showing a higher proportion of DEG under Pb + WS compared to WS conditions (36% vs 32% respectively; Table S7 and Figure 5B). A close look revealed that the extent of differential expression for many DEG genes was also greater for Pb + WS than WS alone (Table S6). Module 29 was clustered to meta-module G9 with module 26 and expressed during mid-to-end maturation (Figure 4A).

## 3.6 | Combined stress conditions compared to single drought stress promote the expression of HSFA2, HSP and chaperone genes

Since meta-module G9 might provide a clue to explain the difference of impacts on the measured parameters between combined versus single stress, we conducted an enrichment analysis on the 65 DEG identified in modules 26 and 29. Meta-module G9 was strongly enriched in genes related to response to heat (GO:0009408; 77.4-fold enrichment, *p*-value 1.11e-22), response to temperature stimulus (GO:0009266; 30.5-fold enrichment, *p*-value 8.83e-17) or protein folding (GO:0006457; 77.15-fold enrichment, *p*-value 5.39e-15). Most of the genes encoded small heat shock proteins (sHSP), with

**FIGURE 4** A robust gene co-expression network of seed development in *B. napus*. A. Heatmap of mean module eigengenes (ME) for each 'stage × treatment' combination. ME are mean values of three or two biological replicates for C and WS or Pb and Pb + WS samples respectively. The module network dendogram above the heatmap was constructed by clustering ME distances. The red vertical line represents the threshold (0.983) used for defining the ten meta-modules (G1 to G10). B. Eigengene network visualization and annotation. Network nodes are modules identified from the WGCNA. The graphical representation (ggplot2) used the Euclidean distances between modules. Only the first quartile Euclidean distances with a maximum threshold set at 0.45 are represented. For each meta-module, the five most significant GO terms are indicated.
<figure><img src="image_7.png" /><figcaption>**FIGURE 5** Overview of the differentially expressed genes (DEG) in the seed development gene co-expression network upon stress conditions. A. Venn diagram of the DEG. The number in each circle represents the amount of DEG between the different comparisons (stress vs control) for at least one seed developmental stage. B. For each stress condition, modules are color-coded based on the proportion of differentially expressed genes (DEG) that meet specific criteria: a minimum log2 fold change of 1.5, set to filter genes initially identified as DEG using an FDR-adjusted *p*-value threshold of <0.05 (Benjamini-Hochberg correction) in the comparison of stress vs control conditions. WS, water shortage; Pb, *P. brassicae* inoculation; Pb + WS, *P. brassicae* inoculation and water shortage.</figcaption></figure>

differential expression mainly at the time point T3 of seed development (42 genes out of 65, Figure 6, Table S6). For almost all of the genes, transcript levels were higher for the Pb + WS compared to the WS (Figure 6A). Considering that the upstream regulators of HSP are heat stress transcription factors (HSF), we searched for DEG encoding HSF. In module 13, with expression profiles slightly preceding that of meta-module G7, we identified two HSFA2 orthologs (*BnHSFA2*.*A*, BnaA03g22890D and *BnHSFA2.C*, BnaC03g26940D) that showed stronger upregulation in Pb + WS compared to WS (Figures 6B, C). Several deregulated genes identified in Modules 26 and 29 with higher expression in Pb + WS compared to the WS (sHSP orthologs, *HSP101*, Bcl-2-associated athanogene 6 (BAG6), *GALACTINOL SYNTHASE 1* (GOLS1)) are known targets of HSFA2 in Arabidopsis (Nishizawa et al., 2006), suggesting that these two genes could be the underlying activators of the increased sHSP expression (Figures 6B, C).

## 3.7 | Quantification of module-trait associations identifies two gene modules linked to longevity acquisition under stress conditions

Our next objective was to incorporate traits linked to the accumulation of seed storage reserves (FA, proteins) and longevity (P50) into the gene network to identify modules and hub genes among the DEG that might govern these traits under stress conditions. To achieve this, we first determined which ME significantly correlated with each of the measured traits based on the Pearson correlation coefficient using a trait-based gene significance measure (Figure S5A). The fatty acid accumulation showed significant correlations with modules 26 (r = 0.8) and 29 (r = 0.74), which are the two modules containing the strong enrichment in heat shock proteins, and module 16 (r = 0.79). However, no genes that could explain a role in lipid accumulation were identified. Storage protein accumulation did not correlate with any specific module (Figure S5A). However, for the longevity trait, modules 06 (611 genes) and 12 (138 genes) displayed very significant correlations with the acquisition of seed longevity (r = 0.94) (Figure S5A). A GO enrichment analysis of these two modules revealed functions in ribosome biogenesis, RNA processing, fatty acid biosynthetic process, pseudouridine synthesis and chromatin silencing at telomere (Table S5).

Since genes with high connectivity are more likely to exert large effects on physiological traits (Ghazalpour et al., 2006), those genes were identified in both modules with the highest intramodular connectivity (MM) as well as with the highest correlation with the acquisition of longevity (GS). Therefore, the GS values of the DEG found in modules 06 and 12 were plotted against the MM values of each individual gene (Figure S5B). The resulting correlations allowed us to identify 46 and 2 hub genes having both high MM and GS values >0.8 for modules 06 and 12, respectively (Figure S5A; Table S8). The gene list was mostly enriched with GO terms related to ribosomal large subunit biogenesis, tRNA processing and ncRNA metabolic processes, similar functions as those observed for genes linked to the acquisition of longevity in *Medicago truncatula*, *Arabidopsis thaliana* (Righetti et al., 2015)
<figure><img src="image_8.png" /><figcaption>**FIGURE 6** Stronger upregulation of small heat shock proteins (sHSP) enriched meta-module G9 genes by the combined stress treatment compared to single stress treatment. **A.** Heatmap depicting the transcript level of the meta-module G9 genes (modules 26 and 29) that are differentially expressed by stress treatments at time T3 (520 GDD) using mean centered log2 TMM values. Also added are the two *HSFA2* genes from Module 13 that are upregulated in Pb + WS and WS. DEG with a minimum log2 fold change of 1.5 and FDR adjusted *p*-value <0.05 (Benjamini-Hochberg correction) (stress vs control). **B-D.** Expression profiles of *HSFA2A.A* (BnaA03g22890D) (B), *HSFA2.C* (BnaC03g26940D (C), and *HSP18.5* (BnaAnng26270D) (D). Means of gene expression levels are given in CPM with respective standard errors. Statistical levels of significance (normalized log2 transcript levels of treated vs control samples): none, non-significant; \* *p* < 0.05; \*\* *p* < 0.01; \*\*\* *p* < 0.001. Gene ID are given according to the reference gene annotation (Darmor-bzh-v5). C, control; WS, water shortage; Pb, *P. brassicae* inoculation; Pb + WS, *P. brassicae* inoculation and water shortage.</figcaption></figure>

and tomato (Bizouerne et al., 2021) (Table S8). Amongst the hub genes that were differentially expressed in WS (FC in WS > 3FC in Pb + WS), we identified the ortholog of DORMANCY ASSOCIATED GENE 2 (DAP2; BnaA05g09830D) that is also referred to as DRM2 in Arabidopsis (Figure S5C). Gene annotation revealed that AtDAP2 is a pathogen-inducible gene that negatively regulates disease defence against virulent pathogens and systemic acquired resistance (Roy et al., 2020).

## 3.8 | Genes encoding precursors of seed storage proteins are upregulated when drought is applied as a single stress

Although there was a significant increase in seed protein content in WS seeds compared to the other treatments (Figure 2E), no module was identified in correlation with protein content (Figure S5A). Therefore, we used a candidate gene-based strategy and projected a list of

genes related to seed storage protein (SSP) on the network. Our criteria included genes already annotated as protein synthesis precursors in the Darmor-bzh-v5 genome and those whose Arabidopsis orthologs are involved in seed protein synthesis (Verdier et al., 2008). This analysis revealed a high number of SSP genes in module 13, which contains 137 genes expressed during seed filling. Specifically, among the 72 most highly connected genes in the module (MM >0.82), 87% exhibited deregulation under WS conditions and were enriched in SSP genes. This enrichment included a total of 12 SSP genes, comprising eight encoding 2S napins and four encoding 12S cruciferins (Table 1, Table S10). SSP related genes were expressed between T1 and T6 seed stages with a peak of expression at T3 for the 2S and a bit later at T3-T4 for the 12S followed by a decrease from T4-T5 onwards (Figures 7A-D). Genes encoding 2S, such as *Bn2S.C* (BnaC01g19320D) (Figure 7A) and *Bn2S.A* (BnaA03g48490D) (Figure 7B) displayed a striking upregulation in WS seeds. Contrary to what was observed for most genes in our study, these differences did not result from a more
<table><caption>TABLE 1 Impact of treatments on gene expression during seed development: focus on genes of module M13. The 72 most correlated genes of the module M13 are presented with module membership (MM) correlation value and corresponding *p*-value, maximum fold change differences between transcript levels in seeds produced under stress and control conditions (based on log2 (CPM+1) values). Arabidopsis orthologs information (accession number and gene function), and presence of putative transcription factor binding sequences in the corresponding promoter (with respective *p*-values). Genes coding for seed storage protein precursors are written in boldface. WS, water shortage; Pb, *P. brassicaceae* inoculation; Pb+WS, *P. brassicaceae* inoculation and water shortage. The genes related to seed storage protein within the 72 most correlated genes of the module M13 are presented with module membership (MM) correlation value and corresponding *p*-value, maximum fold change differences between transcript levels in seeds produced under stress and control conditions (based on log2(CPM+1) values). Arabidopsis orthologs information (accession number and gene function), and presence of putative transcription factor binding sequences in the corresponding promoter (with respective *p*-values). Genes coding for seed storage protein precursors are written in boldface. WS, water shortage; Pb, *P. brassicaceae* inoculation; Pb+WS, *P. brassicaceae* inoculation and water shortage.</caption><thead><tr><th rowspan="2">B.napus</th><th colspan="3">Module 13</th><th colspan="3">Max FC_log2 (CPM+1)</th><th colspan="3">Transcription factor binding sequences</th></tr><tr><th>MM</th><th>p-value</th><th>WS</th><th>Pb+WS</th><th>Pb</th><th>Arabidopsis ortholog ID</th><th>Arabidopsis ortholog name / function</th><th>Gen. bHLH (7.88e-05)</th><th>FUS3 (4.26e-04)</th><th>bZIP68 (4.46e-04)</th><th>EMB1789 (1.17e-03)</th><th>Gen. GRF (1.80e-03)</th></tr></thead><tbody><tr><td>BnaC01g19320D</td><td>0.957</td><td>1.99E-39</td><td>3.91</td><td></td><td></td><td>AT4G27170 Napin (2S) precursor</td><td></td><td>X</td><td>X</td><td>X</td><td>X</td><td></td></tr><tr><td>BnaA06g05220D</td><td>0.851</td><td>2.98E-21</td><td>-1.93</td><td></td><td></td><td>AT1G08940 Phosphoglycerate mutase family protein</td><td></td><td>X</td><td></td><td>X</td><td></td><td></td></tr><tr><td>BnaA03g35700D</td><td>0.958</td><td>7.52E-40</td><td>1.88</td><td></td><td></td><td>AT3G20380 TRAF-like family protein</td><td></td><td>X</td><td>X</td><td></td><td></td><td>X</td></tr><tr><td>BnaC03g48750D</td><td>0.909</td><td>2.68E-28</td><td>1.85</td><td></td><td></td><td>AT5G65890 ACR1: amino acid binding</td><td></td><td>X</td><td></td><td></td><td></td><td></td></tr><tr><td>BnaA09g06380D</td><td>0.960</td><td>2.51E-40</td><td>-1.78</td><td></td><td></td><td>NA Thionin family protein</td><td></td><td>X</td><td></td><td></td><td></td><td></td></tr><tr><td>BnaA09g50780D</td><td>0.981</td><td>1.03E-51</td><td>1.72</td><td></td><td></td><td>AT1G03890 Cruciferin (12S) precursor</td><td></td><td></td><td></td><td>X</td><td></td><td></td></tr><tr><td>BnaC01g01280D</td><td>0.899</td><td>9.48E-27</td><td>1.67</td><td></td><td></td><td>AT4G37760 SQE3: squalene epoxidase 3</td><td></td><td>X</td><td></td><td></td><td></td><td></td></tr><tr><td>BnaCnn643010D</td><td>0.971</td><td>3.13E-45</td><td>-1.58</td><td></td><td></td><td>AT5G51500 PME59: pectinesterase/pectinesterase inhibitor 59</td><td></td><td>X</td><td></td><td></td><td></td><td>X</td></tr><tr><td>BnaA01g31980D</td><td>0.878</td><td>3.78E-24</td><td>-1.51</td><td></td><td></td><td>AT3G10080 RmlC-like cupins superfamily protein</td><td></td><td></td><td></td><td></td><td>X</td><td></td></tr><tr><td>BnaC08g30880D</td><td>0.968</td><td>8.93E-44</td><td>1.51</td><td></td><td></td><td>AT3G60730 PME36: pectinesterase/pectinesterase inhibitor 36</td><td></td><td>X</td><td></td><td></td><td>X</td><td></td></tr><tr><td>BnaA06g34020D</td><td>0.948</td><td>1.78E-36</td><td>1.46</td><td></td><td></td><td>AT2G02780 Leucine-rich repeat protein kise family protein</td><td></td><td></td><td></td><td></td><td></td><td></td></tr><tr><td>BnaC09g647920D</td><td>0.917</td><td>1.42E-29</td><td>1.45</td><td></td><td></td><td>AT5G07500 TZF6/PE1: embryo specific zinc finger TF</td><td></td><td>X</td><td>X</td><td>X</td><td>X</td><td></td></tr><tr><td>BnaA02g14080D</td><td>0.914</td><td>4.12E-29</td><td>1.45</td><td></td><td></td><td>NA Chaperonin 60 subunit beta 4</td><td></td><td></td><td></td><td></td><td></td><td></td></tr><tr><td>BnaC01g19300D</td><td>0.948</td><td>1.36E-36</td><td>1.45</td><td></td><td></td><td>AT4G27170 Napin (2S) precursor</td><td></td><td>X</td><td>X</td><td>X</td><td>X</td><td></td></tr><tr><td>BnaC08g41500D</td><td>0.856</td><td>9.68E-22</td><td>-1.39</td><td></td><td></td><td>AT1G12200 FMO: putative flavin monooxygese</td><td></td><td>X</td><td></td><td></td><td>X</td><td></td></tr><tr><td>BnaA10g23630D</td><td>0.939</td><td>3.22E-34</td><td>1.39</td><td></td><td></td><td>AT5G07190 ATS3: seed gene 3</td><td></td><td></td><td></td><td></td><td></td><td></td></tr><tr><td>BnaA01g29070D</td><td>0.969</td><td>3.83E-44</td><td>1.35</td><td></td><td></td><td>AT3G14950 TTL2: tetratricopeptide repeat thioredoxin like 2</td><td></td><td></td><td></td><td>X</td><td></td><td></td></tr><tr><td>BnaC07g21880D</td><td>0.861</td><td>2.78E-22</td><td>1.35</td><td></td><td></td><td>AT2G02780 Leucine-rich repeat protein kise family protein</td><td></td><td></td><td></td><td></td><td></td><td></td></tr><tr><td>BnaA03g22670D</td><td>0.914</td><td>4.40E-29</td><td>-1.34</td><td></td><td></td><td>AT2G26580 YAB5: transcription factor YABBY family protein</td><td></td><td>X</td><td></td><td></td><td></td><td></td></tr><tr><td>BnaA06g07870D</td><td>0.973</td><td>3.14E-46</td><td>1.32</td><td></td><td></td><td>AT1G12850 Phosphoglycerate mutase family protein</td><td></td><td></td><td></td><td></td><td>X</td><td></td></tr><tr><td>BnaC08g43570D</td><td>0.974</td><td>4.13E-47</td><td>1.31</td><td></td><td></td><td>AT1G03890 Cruciferin (12S) precursor</td><td></td><td></td><td></td><td>X</td><td>X</td><td></td></tr><tr><td>BnaC04g25160D</td><td>0.966</td><td>1.16E-42</td><td>-1.28</td><td></td><td></td><td>NA Kinase family protein</td><td></td><td></td><td></td><td></td><td></td><td></td></tr><tr><td>BnaC01g43250D</td><td>0.951</td><td>2.50E-37</td><td>1.28</td><td></td><td></td><td>NA Napin (2S) precursor</td><td></td><td>X</td><td>X</td><td>X</td><td>X</td><td></td></tr><tr><td>BnaA10g02240D</td><td>0.972</td><td>8.47E-46</td><td>1.26</td><td></td><td></td><td>AT1G03890 Cruciferin (12S) precursor</td><td></td><td>X</td><td></td><td></td><td>X</td><td></td></tr><tr><td>BnaC05g02160D</td><td>0.970</td><td>1.04E-44</td><td>1.24</td><td></td><td></td><td>AT1G03890 Cruciferin (12S) precursor</td><td></td><td>X</td><td></td><td></td><td>X</td><td></td></tr><tr><td>BnaA10g27740D</td><td>0.921</td><td>2.10E-30</td><td>1.24</td><td></td><td></td><td>AT1G04380 ACO: 1 aminocyclopropane 1 carboxylic acid oxidase</td><td></td><td></td><td></td><td></td><td></td><td></td></tr></tbody></table>

(Continues)
<table><caption>TABLE 1 (Continued)</caption><thead><tr><th rowspan="2">B.napus</th><th colspan="2">Module 13</th><th colspan="3">Max FC_log2 (CPM+1)</th><th colspan="5">Transcription factor binding sequences</th></tr><tr><th>MM</th><th>p-value</th><th>WS</th><th>Pb+WS</th><th>Pb</th><th>Arabidopsis ortholog ID</th><th>Arabidopsis ortholog name / function</th><th>Gen. bHLH (7.88e-05)</th><th>FUS3 (4.26e-04)</th><th>bZIP68 (4.46e-04)</th><th>EMB1789 (1.17e-03)</th><th>Gen. GRF (1.80e-03)</th></tr></thead><tbody><tr><td>BnaA02g36500D</td><td>0.931</td><td>2.06E-32</td><td>-1.20</td><td></td><td></td><td>AT4G02370</td><td>Probable pectinesterase</td><td>X</td><td></td><td></td><td></td><td></td></tr><tr><td>BnaA01g17200D</td><td>0.942</td><td>5.13E-35</td><td>1.20</td><td></td><td></td><td>NA</td><td>Napin (2S) precursor</td><td>X</td><td>X</td><td>X</td><td>X</td><td></td></tr><tr><td>BnaA04g00580D</td><td>0.960</td><td>1.39E-40</td><td>-1.20</td><td></td><td></td><td>AT3G61680</td><td>PLIP1: plastid localized phospholipase</td><td>X</td><td></td><td></td><td>X</td><td></td></tr><tr><td>BnaC06g19660D</td><td>0.941</td><td>9.93E-35</td><td>-1.19</td><td></td><td></td><td>AT1G79690</td><td>Nudix hydrolase 3</td><td></td><td></td><td></td><td></td><td></td></tr><tr><td>BnaC03g21710D</td><td>0.967</td><td>2.18E-43</td><td>-1.18</td><td></td><td></td><td>AT2G38900</td><td>Predicted PR-6 protease inhibitor family</td><td></td><td></td><td></td><td></td><td>X</td></tr><tr><td>BnaC02g10990D</td><td>0.857</td><td>8.56E-22</td><td></td><td>-2.47</td><td></td><td>AT5G58070</td><td>TIL: temperature-induced lipocalin</td><td>X</td><td></td><td></td><td></td><td>X</td></tr><tr><td>BnaA0ng06380D</td><td>0.867</td><td>7.76E-23</td><td>-2.00</td><td></td><td></td><td>AT5G07680</td><td>C080</td><td></td><td></td><td></td><td></td><td></td></tr><tr><td>BnaC03g70800D</td><td>0.915</td><td>2.23E-29</td><td>-1.92</td><td></td><td></td><td>NA</td><td>Unknown</td><td></td><td></td><td></td><td></td><td></td></tr><tr><td>BnaA03g22890D</td><td>0.828</td><td>3.08E-19</td><td></td><td>1.49</td><td></td><td>AT2G26150</td><td>HSFA2: heat shock TF</td><td>X</td><td></td><td></td><td></td><td>X</td></tr><tr><td>BnaC03g07860D</td><td>0.928</td><td>1.02E-31</td><td></td><td>3.23</td><td></td><td>NA</td><td>Serine carboxypeptidase like 7</td><td></td><td></td><td></td><td></td><td></td></tr><tr><td>BnaC07g39000D</td><td>0.899</td><td>7.77E-27</td><td></td><td>1.88</td><td></td><td>AT4G24550</td><td>AP4M: SSPs vacuolar sorting</td><td></td><td></td><td></td><td></td><td>X</td></tr><tr><td>BnaA09g07630D</td><td>0.944</td><td>1.79E-35</td><td></td><td>1.49</td><td></td><td>AT5G66985</td><td>HUP44</td><td>X</td><td></td><td></td><td>X</td><td></td></tr><tr><td>BnaC04g03970D</td><td>0.821</td><td>1.11E-18</td><td></td><td>1.41</td><td></td><td>AT2G45080</td><td>Cyclin P3 1</td><td>X</td><td></td><td></td><td></td><td></td></tr><tr><td>BnaC05g07030D</td><td>0.900</td><td>5.91E-27</td><td>5.28</td><td>5.02</td><td></td><td>AT1G09510</td><td>Cinamoyl-reductase 1 like</td><td>X</td><td></td><td></td><td></td><td></td></tr><tr><td>BnaA03g48490D</td><td>0.990</td><td>1.86E-61</td><td>2.45</td><td>1.47</td><td></td><td>AT4G27170</td><td>Napin (2S) precursor</td><td>X</td><td>X</td><td></td><td>X</td><td>X</td></tr><tr><td>BnaC07g40770D</td><td>0.986</td><td>4.33E-56</td><td>2.37</td><td>1.41</td><td></td><td>AT4G27170</td><td>Napin (2S) precursor</td><td>X</td><td>X</td><td>X</td><td>X</td><td>X</td></tr><tr><td>BnaA06g02910D</td><td>0.878</td><td>4.24E-24</td><td>-2.04</td><td>-2.25</td><td></td><td>NA</td><td>delta-sterol reductase</td><td>X</td><td></td><td></td><td></td><td></td></tr><tr><td>BnaC04g50170D</td><td>0.904</td><td>1.44E-27</td><td>-2.08</td><td>-2.19</td><td></td><td>AT2G45220</td><td>PME17: pectin methyltransferase</td><td>X</td><td></td><td></td><td></td><td></td></tr><tr><td>BnaC02g01470D</td><td>0.873</td><td>1.53E-23</td><td>-1.77</td><td>-1.97</td><td></td><td>AT5G07680</td><td>C080</td><td></td><td></td><td></td><td></td><td></td></tr><tr><td>BnaA05g18360D</td><td>0.847</td><td>6.52E-21</td><td>-1.77</td><td>-1.91</td><td></td><td>AT1G34630</td><td>Transmembrane protein</td><td>X</td><td></td><td></td><td>X</td><td>X</td></tr><tr><td>BnaA01g11750D</td><td>0.943</td><td>3.68E-35</td><td>-1.83</td><td>-1.54</td><td></td><td>AT4G21902</td><td>Unknown</td><td></td><td></td><td></td><td></td><td></td></tr><tr><td>BnaC01g23420D</td><td>0.976</td><td>5.29E-48</td><td>-1.70</td><td>-1.64</td><td></td><td>AT3G51450</td><td>Calcium-dependent phosphotriesterase superfamily protein</td><td></td><td></td><td></td><td></td><td></td></tr><tr><td>BnaA05g02120D</td><td>0.963</td><td>2.06E-41</td><td>-1.64</td><td>-1.42</td><td></td><td>AT2G41290</td><td>SSL2: probable strictosidine synthase</td><td></td><td></td><td></td><td></td><td></td></tr><tr><td>BnaC06g38080D</td><td>0.975</td><td>2.01E-47</td><td>-1.63</td><td>-1.50</td><td></td><td>AT1G77380</td><td>AAP3: amino acid permease</td><td>X</td><td></td><td></td><td>X</td><td></td></tr><tr><td>BnaC03g40770D</td><td>0.895</td><td>2.61E-26</td><td>-1.50</td><td>1.56</td><td></td><td>AT3G19010</td><td>2OGX: 2-Oxoglutarate Oxygeses</td><td>X</td><td></td><td>X</td><td></td><td>X</td></tr><tr><td>BnaA08g23110D</td><td>0.946</td><td>7.55E-36</td><td>1.55</td><td>1.44</td><td></td><td>AT1G17430</td><td>Alpha/Beta-hydrolase protein family</td><td></td><td></td><td></td><td></td><td></td></tr><tr><td>BnaC07g32100D</td><td>0.866</td><td>9.15E-23</td><td>-1.47</td><td>1.48</td><td></td><td>AT3G51150</td><td>Kinesin like protein</td><td></td><td></td><td></td><td></td><td></td></tr><tr><td>BnaA01g16180D</td><td>0.970</td><td>8.75E-45</td><td>1.45</td><td>1.15</td><td></td><td>AT4G27170</td><td>Napin (2S) precursor</td><td>X</td><td>X</td><td>X</td><td>X</td><td>X</td></tr><tr><td>BnaCnn67960D</td><td>0.941</td><td>9.68E-35</td><td>-1.44</td><td>-1.36</td><td></td><td>NA</td><td>GEBP: GLABROUS1 enhancer-binding protein</td><td></td><td></td><td></td><td></td><td>X</td></tr><tr><td>BnaC03g08840D</td><td>0.959</td><td>3.10E-40</td><td>1.42</td><td>1.15</td><td></td><td>NA</td><td>Napin (2S) precursor</td><td>X</td><td>X</td><td>X</td><td>X</td><td>X</td></tr><tr><td>BnaA03g55570D</td><td>0.982</td><td>5.17E-52</td><td>1.36</td><td>1.31</td><td></td><td>AT5G07010</td><td>ST2A: sulfotransferase</td><td></td><td></td><td></td><td></td><td>X</td></tr></tbody></table>
<figure><img src="image_9.png" /></figure>

<table><caption>TABLE 1 (Continued)</caption><thead><tr><th rowspan="2">B.napus</th><th colspan="2">Module 13</th><th colspan="3">Max FC_log2 (CPM+1)</th><th rowspan="2">Arabidopsis ortholog ID</th><th rowspan="2">Arabidopsis ortholog name / function</th><th colspan="5">Transcription factor binding sequences</th></tr><tr><th>MM</th><th>p-value</th><th>WS</th><th>Pb+WS</th><th>Pb</th><th>fnKC: Probable ictive serine/threonine-protein kise</th><th>FUS3 (7.88e-05)</th><th>bZIP68 (4.26e-04)</th><th>EMB1789 (1.17e-03)</th><th>Gen. bHLH (1.80e-03)</th></tr></thead><tbody><tr><td>BnaC09g03710D</td><td>0.979</td><td>2.36E-50</td><td>-1.34</td><td>-1.25</td><td></td><td>NA</td><td>fnKC: Probable ictive serine/threonine-protein kise</td><td>X</td><td></td><td>X</td><td></td><td></td></tr><tr><td>BnaA09g42760D</td><td>0.908</td><td>4.13E-28</td><td>1.18</td><td>1.28</td><td></td><td>AT2G22240</td><td>MIPS2: myo inositol phosphate synthase 2</td><td>X</td><td>X</td><td>X</td><td></td><td></td></tr><tr><td>BnaC02g01710D</td><td>0.940</td><td>1.71E-34</td><td>-1.22</td><td>-1.15</td><td></td><td>AT5G07200</td><td>GA20OX3: gibberellin 20 oxidase 3 like</td><td></td><td></td><td></td><td></td><td></td></tr><tr><td>BnaC09g07570D</td><td>0.922</td><td>1.23E-30</td><td>3.08</td><td></td><td>1.70</td><td>AT5G66985</td><td>HUP44</td><td>X</td><td></td><td></td><td>X</td><td></td></tr><tr><td>BnaA05g04170D</td><td>0.848</td><td>5.69E-21</td><td>1.90</td><td></td><td>1.99</td><td>AT2G44600</td><td>Unknown</td><td>X</td><td></td><td></td><td></td><td></td></tr><tr><td>BnaA04g09880D</td><td>0.909</td><td>2.92E-28</td><td>1.55</td><td></td><td>1.49</td><td>AT2G04690</td><td>PNOX: pyridoxamine phosphate oxidase family protein</td><td></td><td>X</td><td></td><td></td><td></td></tr><tr><td>BnaA02g24250D</td><td>0.969</td><td>3.80E-44</td><td>1.14</td><td></td><td>1.33</td><td>AT5G46050</td><td>PTR3: peptide transporter</td><td>X</td><td></td><td></td><td></td><td></td></tr><tr><td>BnaC03g26940D</td><td>0.831</td><td>1.55E-19</td><td></td><td>1.40</td><td>1.45</td><td>AT2G26150</td><td>HSFA2: heat shock TF</td><td>X</td><td>X</td><td>X</td><td></td><td></td></tr><tr><td>BnaA08g31350D</td><td>0.905</td><td>1.27E-27</td><td>5.05</td><td>1.70</td><td>1.74</td><td>NA</td><td>TGG2: myrosise 2 precursor</td><td></td><td></td><td></td><td></td><td></td></tr><tr><td>BnaC02g11470D</td><td>0.964</td><td>4.62E-42</td><td>2.07</td><td>-1.30</td><td>1.35</td><td>NA</td><td>XTH protein family</td><td>X</td><td></td><td>X</td><td></td><td></td></tr><tr><td>BnaA02g35320D</td><td>0.976</td><td>5.75E-48</td><td>-1.78</td><td>-1.43</td><td>1.46</td><td>NA</td><td>XTH protein family</td><td>X</td><td></td><td>X</td><td></td><td></td></tr><tr><td>BnaA07g24400D</td><td>0.858</td><td>6.05E-22</td><td>1.75</td><td>1.67</td><td>1.64</td><td>NA</td><td>TGG2: myrosise 2 precursor</td><td></td><td></td><td>X</td><td></td><td></td></tr><tr><td>BnaA09g26590D</td><td>0.846</td><td>8.49E-21</td><td>1.58</td><td>1.50</td><td>1.39</td><td>NA</td><td>Myrosine</td><td></td><td></td><td>X</td><td></td><td></td></tr><tr><td>BnaA01g02380D</td><td>0.966</td><td>5.95E-43</td><td>1.51</td><td>1.44</td><td>1.38</td><td>AT4G34810</td><td>SAUR5: small auxin upregulated R 5</td><td></td><td></td><td></td><td></td><td></td></tr><tr><td>BnaA04g25160D</td><td>0.981</td><td>7.49E-52</td><td>1.28</td><td>1.27</td><td>1.34</td><td>AT2G43535</td><td>DEFL: defensin like protein</td><td></td><td></td><td></td><td>X</td><td></td></tr></tbody></table>

<figure><img src="image_10.png" /></figure>
precocious transcript accumulation in WS seeds compared to the other treatment but from an increase in overall transcript level. In contrast, the differences observed in 12S gene expression levels between WS and control conditions (Figures 7C, D) were due to precocious transcription of the 12S genes upon drought. In addition, deregulation of the SSP genes was not observed for the Pb + WS samples, suggesting that precocious inoculation with *P. brassicae* abolishes the impact of WS on the gene transcription, consistent with the difference in protein content between WS and Pb + WS (Figure 2E).

To get further insight into the upstream regulators of the SSP genes, we analysed the potential cis-element enrichment in promoter sequences of the most highly connected genes of module 13 using the SeqEnrich software (Becker et al., 2017). Promoter sequences of 2S and 12S genes contained binding sequences for transcription factors such as FUS3 and bZIP68 that have been shown to regulate seed development (Table 1, Table S10). bZIP68 is a redox-sensitive transcription factor that normally suppresses the expression of stress tolerance genes and promotes the expression of growth-related genes, and its inactivation by oxidative modification enhances stress tolerance but suppresses growth (Li et al., 2019). bZIP68 might balance stress tolerance with growth through the extent of its oxidative inactivation according to the environment.

The only gene of module 13 with similar promoter binding sequences to the SSP genes was BnaC09g47920D (Table 1, Table S10), encoding an ortholog of Arabidopsis PEI1/TZF6, an embryo-specific transcription factor required at the globular to heart stage transition in embryo development (Li and Thomas, 1998; Bogamuwa and Jang, 2016). As shown for the SSP genes, *BnPEI1/TZF6* was significantly upregulated when WS was applied alone. Therefore, we concluded that module 13 contained a number of critical target genes whose modulations could explain the increase in seed protein content upon drought. However, it appears that these genes were not deregulated under combined stresses.

## 3.9 | Impact of water deficit on transcription timing of genes involved in fatty acid synthesis

Since drought stress also induced a reduction of oil accumulation in seeds (Figure 2A), we searched whether expression patterns of genes related to fatty acid (FA) biosynthesis were modulated under the applied stresses. For that purpose, we selected genes directly involved in FA biosynthesis based on previous knowledge (Baud et al., 2008) and looked for the OSR orthologs throughout the gene co-expression network. In total, we found 47 genes from modules 01, 05, 06, 07 and 10 that covered the entire process of FA biosynthesis and whose transcript levels were significantly impacted by the stress, especially WS, during seed filling (Table S9). Indeed, expression levels of most FA-related genes peaked during early seed development (T2-T3) and decreased afterwards. Expression profiles were similar under stress conditions but generally decreased earlier, especially under WS as observed for genes encoding proteins involved in the first step of FA biosynthesis (*BnPKPα*, BnaA01g24280D and *BnBCCP2*, BnaA03g05470D), oleic acid

saturation (*BnFAB2*, BnaA03g20420D) or triacylglycerol assembly (*BnDAGAT*, BnaA01g19390D) (Figures 7E-H). In addition, 13 of the 47 genes selected within the co-expression network encoded oil body proteins such as *BnHSD1* (BnaCnng57830D) or *BnCLO1* (BnaC09g31660D) (Figures 7I-J, Table S9). In contrast with the other FA genes, transcript levels of the oil body related genes increased during seed filling and reached a peak at T5-T6, with an earlier accumulation under WS conditions (Table S9). Finally, we identified four genes coding for AWPM19-like (ABA-induced Wheat Plasma Membrane Polypeptide-19 like) proteins (Table S9) within modules 05 and 07 that showed similar expression profiles to the oil body genes (Figure 7K; BnaC08g50170D). Interestingly, predicted functional partners of the Arabidopsis AWPM19 gene (At1g04560) included oleosin (OLE1),caleosin (CLO1) and cysteine peroxiredoxin (PER1), are all present in module 07 (Table S9) and display comparable gene expression profiles under stress or control conditions (Figure 7L). In conclusion, expression profiles of FA-related genes displayed a shift with a faster decrease of transcript levels under WS conditions, which was not observed when plants were precociously inoculated with *P. brassicae*. On the other hand, oil body related genes were found to be expressed earlier in stressed seeds, with no significant changes in the maximum expression levels. Moreover, we identified the BnAWPM19-like genes as potential regulators of the oil body related genes in OSR seeds.

# 4 | DISCUSSION

## 4.1 | Seed development is resilient to moderate environmental stresses experienced by the mother plants

Under the stress conditions applied in our study, we showed that physiological and transcriptomic modulations were mostly related to the seed developmental stages rather than to the constraints themselves. These findings could be explained by a combination of factors regarding the stress modalities and plant responses.

Concerning the biotic stress, we chose to infect our plants with a dose of inoculum ten times lower than that usually used in clubroot assays of Brassicaceae species ($10^6$ instead of $10^7$ resting spores/ml; Manzanares-Dauleux et al., 2000b; Strelkov et al., 2021). This choice was made in order not to compromise the survival of the plants and to secure seed production, knowing that at this concentration of inoculum, the infection can also cause significant symptoms (Strehlow et al., 2015). For the same reasons, we used the OSR Express genotype, known to be moderately susceptible to the *P. brassicae* eH isolate (Aigu et al., 2020). This may explain why only 40% of the plants analyzed showed clear root galls. Therefore, the “Pb” seed pool originated from at least 40% of infected plants, potentially more, and allowed us to examine the pathogen’s impact on seed development.

Regarding the drought condition, the major threat to seed development came from the reduction of the mother plant photosynthesis production (Sehgal et al., 2018). In our experiment, water shortage
<figure><img src="image_11.png" /><figcaption>**FIGURE 7** Effect of single and combined stress on expression profiles of genes encoding seed storage compounds. **A-D.** Seed storage protein precursors and co-expressed genes: *Bn2S.C* (BnaC01g19320D) (**A**), *Bn2S.A* (BnaA03g48490D) (**B**), *Bn12S* (BnaC08g43570D) (**C**), and *BnPEI1/TZF6* (BnaC09g47920D) (**D**). **E-L.** Fatty acid biosynthesis and co-expressed genes: *BnPKPa* (BnaA01g24280D) (**E**), *BnBCCP2* (BnaA03g05470D) (**F**), *BnFAB2* (BnaA03g20420D) (**G**), *BnDAGAT* (BnaA01g19390D) (**H**), *BnHSD1* (BnaCnng57830D) (**I**), *BnCLO1* (BnaC09g31660D) (**J**), *BnAWPM19-like* (BnaC08g50170D) (**K**), and *BnPER1* (BnaCnng77870D) (**L**). Means of gene expression levels are given in CPM with respective standard errors. Statistical levels of significance (normalized log2 transcript levels of treated vs control samples): none, non-significant; \*p < 0.05; \*\*p < 0.01; \*\*\*p < 0.001. Gene ID are given according to the reference gene annotation (Darmor-bzh-v5). C, control; WS, water shortage; Pb, *P. brassicae* inoculation; Pb + WS, *P. brassicae* inoculation and water shortage; GDD, growing degree-days. C, control; WS, water shortage; Pb, *P. brassicae* inoculation; Pb + WS, *P. brassicae* inoculation and water shortage; T1-T9, timepoints for seed sampling (refer to Figure 1).</figcaption></figure>

reduced the total plant biomass and seed yield by around 35% (Bianchetti et al. in prep.), which was consistent with previous reports on OSR response to drought (Champolivier and Merrien, 1996; Hatzig et al., 2018). In addition, studies in a number of plant species have demonstrated that water deficits occurring at flowering led preferentially to seed abortion rather than an incomplete seed filling. Indeed, a reduction in seed number is a way to ensure a sufficient amount of nutrients is supplied to the remaining sinks (Sun et al., 2004; Liu et al., 2013). We thus hypothesized that seed development was resilient to the drought conditions because of a sufficient nutrient allocation ensured by a reduction of seed number more precisely on the branches (Bianchetti et al. in prep.).

## 4.2 | Seedling infection by *P. brassicae* could prime the plants towards subsequent water shortage

Stress monitoring in our experimental design proved that plants grown under abiotic or combined stress conditions were subjected to a similar water deficit. However, stronger modulations of plant and seed development were observed when drought was applied alone compared to when it occurred in combination with *P. brassicae* infection. These findings suggested that *P. brassicae* infection somewhat dampens the effects of drought. Although severe clubroot symptoms were known to damage the plant root system and increase its sensibility to water deficit, the effects of low pathogen pressure on plant drought tolerance remained less
studied (Dixon, 2009; Ramegowda and Senthil-Kumar, 2015). In accordance with our observations, studies demonstrated that a precocious pathogen infection could increase the resistance of several plant species to subsequent drought stress (Aguilar et al., 2017; Xu et al., 2008). In these experiments, plant drought resistance was induced by viral infections promoting the accumulation of osmoprotectants and/or antioxidant metabolites such as raffinose, trehalose, putrescine or proline. Thus, a deeper metabolic analysis of plant vegetative organs would pave the way to new hypotheses to support the possible beneficial effects of *P. brassicae* precocious infection in plant drought responses.

## 4.3 | Accumulation of heat shock proteins is known to participate in plant acclimation to multiple stresses

Our transcriptomic data pointed out several heat shock proteins (HSP), namely *BnHsp17.6* and *BnHsp22*, that were up-regulated during seed development under stress combination. HSP were reported as critical determinants of plant response to drought and/or pathogen infections, notably by preventing the misfolding of proteins and participating in plant innate immunity (Park and Seo, 2015; ul Haq et al., 2019). In particular, Hsp17.6 was positively identified with drought or pathogen resistance in several studies (Barcala et al., 2008). Moreover, it has been demonstrated that genes coding for HSP were highly expressed in Brassicaceae roots infected by *P. brassicae* (Ciaghi et al., 2019). To add to this, a recent transcriptomic study demonstrated that a *P. brassicae* resistant variety of OSR was characterized by higher HSP gene expression in roots during the infection (Summanwar et al., 2021). These findings led us to hypothesize that infection by *P. brassicae* could have triggered the HSP synthesis in plant roots, thus leading to enhanced water stress tolerance. Another hypothesis is that the infection might have triggered a priming stress, eliciting a cellular machinery that allows plants to faithfully evade a later-arriving, often more severe, so-called triggering stress (Hilker et al., 2016). Two *HSFA2* genes were found to be expressed in the Pb and WS seeds (Figure 6, Table S6). HSFA2 has been associated with transcriptional memory of heat stress, which describes a modified transcriptional response following a second stress that depends on the first stress exposure after an intervening period of inactivity, as described for drought (Ding et al. 2012). It would be tempting to suggest that the genes detected in modules 26 and 29 present a subset of memory genes that show a stronger re-activation upon drought after the infection of *P. brassicae*. In Arabidopsis, genes that show this behaviour are highly enriched in H3K4me2 and H3K4me3, making it interesting to look if *P. brassica* infection might modulate chromatin modifications.

## 4.4 | Plant water limitation modulates the acquisition of seed quality traits

Modulations of seed quality traits and related gene expression levels were identified during seed development under stress conditions, more specifically water limitation, in line with previously published studies.

First, water stress usually lowers seed oil content by inducing early embryo development and shortening the seed-filling phase (Baldini et al., 2002; Farooq et al., 2014). In the present study, we showed that expression levels of 29 genes related to fatty acid metabolism decreased faster under WS compared to well-watered conditions, supporting a shortened seed-filling process promoted by drought. At the same time, it was described that seed oil decreased at the very end of seed maturation in Brassicaceae species (Baud et al., 2002), likely due to some enzymatic activities such as the malate synthase (Chia et al., 2005). Interestingly, malate content and malate synthase expression levels were higher in WS mature seeds than in the controls (Bianchetti et al., in prep.). This suggested that a prolonged maturation phase could interact with a faster seed-filling phase, thus resulting in reduced oil accumulation in WS seeds due to limited oil biosynthesis and activated oil degradation process.

Regarding triglyceride packaging in the cell, we identified a huge number of genes related to oil body biogenesis that were deregulated in WS seeds. In addition, four genes encoding *BnAWPM19-like* proteins were strongly associated with oil body-related genes, a result that corroborated the fact that AWPM19-like protein played an important role in drought response and ABA signalling in rice (Yao et al., 2018). In addition, it was recently reported that AWPM19-like functioned as an osmosensor during germination, regulating water uptake of seeds under osmotic stress (Alexander et al., 2020). In our experiment, seeds produced under drought conditions germinated faster, probably because of a better water uptake capacity (Bianchetti et al., in prep.), which could be linked to a critical role of BnAWPM19-like in regulating germination abilities.

WS seeds showed protein enrichment, as has already been described in other studies (Bouchereau et al., 1996). Considering that seed protein content mostly relies on N remobilization from vegetative tissues (Havé et al., 2017), it can be hypothesized that the drought-driven reduction in seed number has increased the amount of N available and allocated per remaining seed, thus at the same time promoting the transcription of seed storage precursors genes. Finally, gene co-expression network led us to point out a specific module that gathered potential regulators of seed storage proteins in addition to the genes encoding 2S and 12S. For instance, module 13 contained the seed-specific *PEI1/TZF6* gene known for its role during embryogenesis (Li and Thomas, 1998; Bogamuwa and Jang, 2016). Altogether, these genes displayed cis-regulatory sequences related to *FUS3* and *bZIP68* binding sites in their promoter, suggesting a common regulatory network. Both transcription factors have been described in the seed storage protein regulatory pathway (Chern et al., 1996) and were shown to be modulated under drought conditions. However, their precise roles in promoting seed protein content under drought remain to be clarified.

## 4.5 | Hub genes correlate with seed longevity acquisition and suggest the potential implication of *BnDAP2*

Seed acquisition of P50 strongly correlated with the eigengene values of modules 06 and 12. Genes of these late maturation modules were mostly involved in the RNA process and ribosome synthesis.
During maturation, seeds accumulate mRNA and ribosome units, both allowing a faster restart of cellular machinery upon imbibition (Basbouss-Serhal et al., 2015). Protective molecules such as LEA are positively involved in desiccation tolerance and, thus, longevity of orthodox seeds (Rajjou and Debeaujon, 2008). In these two modules, 48 hub genes were identified for high connectivity within their respective module and strong positive correlation with P50 evolution. Based on the 'guilt-by-association' principle, these genes represent prime targets for the identification of molecular regulators of seed longevity acquisition (Wolfe et al., 2005). Among these hub genes, BnDAP2 stood out for the striking similarity of his transcript level kinetics with the P50 evolution of all seeds. In Arabidopsis, DAP2 (also referred to as DRM2) is notably involved in the negative regulation of local and systemic acquired resistance against virulent bacterial pathogens (Roy et al., 2020). During the maturation of M. truncatula and Arabidopsis seeds, many genes associated with longevity acquisition are involved in plant defence response (Righetti et al., 2015). Moreover, through a co-expression network analysis, the DAP2/DRM2 splicing variant, DRM1 (Dormancy-Associated Gene-1) is a member of a plant-specific gene family that has been associated with plant dormancy and was similarly pinpointed for its strong correlation with tomato seed longevity acquisition (Rae et al., 2014; Bizouerne et al., 2021). These results suggest that the reduction of DAP2 transcript levels during late maturation of seeds produced under drought conditions could be related to the reduction of their longevity.

## AUTHOR CONTRIBUTIONS

J.B. and N.N. designed the study and obtained funding. G.B., C.B., S.B., and A.T. carried out the experiments. G.B., V.C., F.L., and K.G. analyzed the RNA-seq data. G.B. and V.C. deposited the RNA-seq data in the NCBI Sequence Read Archive. M.J.M.D. brought expertise in clubroot infection. G.B., J.B., and N.N. interpreted the results and wrote the manuscript. All authors read and approved the final version.

## ACKNOWLEDGMENTS

The authors kindly acknowledge colleagues from INRAE (Le Rheu, France), especially Y. Dellero, F. Le Cahérec, and L. Leport for helpful discussions, M. Rousseau-Gueutin and C. Falentin for help with the genomic resources of Brassica napus, L. Charlon, J-P. Constantin, F. Le Tertre and P. Rolland for supplying greenhouse services, A. Carrillo for running the experimental design and supporting the seed sampling, P. Glory and C. Lariagon for assistance with *P. brassicae* infection. We thank the BrA-CySol biological resource center (INRAE Ploudaniel, France) for providing seeds of Express. The BGI sequencing platform (https://www.bgi.com) and the GenOuest bioinformatics platform (https://www.genouest.org/) are acknowledged for production, storage and analysis of the RNAseq data. The work described in this manuscript was performed within the framework of the transnational FACCE-JPI-ERA-NET+ collaborative project SYBRACLIM, with funding from ANR (France) to N.N. (grant ANR-14-JFAC-007-01). J.B. and N.N. were also supported by INRAE (grant INRAE-BAP-IB2017-SQUAL). G.B. was the recipient of a 3 years PhD fellowship from INRAE and the Région de Bretagne.

## CONFLICT OF INTEREST STATEMENT

The authors declare no conflicts of interest.

## DATA AVAILABILITY STATEMENT

Datasets supporting the findings of this study are openly available. The transcriptomic dataset can be accessed at https://doi.org/10.1016/j.dib.2021.107392 and the metabolic/physiologic datasets are available at https://doi.org/10.1016/j.dib.2021.107247.

## ORCID

Grégoire Bianchetti https://orcid.org/0000-0003-1192-982X

Vanessa Clouet https://orcid.org/0000-0001-8541-7291

Fabrice Legeai https://orcid.org/0000-0002-6472-4839

Cécile Baron https://orcid.org/0000-0002-3713-370X

Kévín Gazengel https://orcid.org/0000-0001-9054-2281

Sébastien Baud https://orcid.org/0000-0002-5507-327X

Maria J. Manzanares-Dauleux https://orcid.org/0000-0002-6452-0393

Julia Buitink https://orcid.org/0000-0002-1457-764X

Nathalie Nesi https://orcid.org/0000-0002-2552-472X

## REFERENCES

Aguilar, E., Cutrona, C., Del Toro, F.J., Vallarino, J.G., Osorio, S., Pérez-Bueno, M.L., Barón, M., Chung, B.-N., Canto, T. and Tenllado, F. (2017) Virulence determines beneficial trade-offs in the response of virus-infected plants to drought via induction of salicylic acid. *Plant Cell Environment*, 40, 2909–2930. https://doi.org/10.1111/pce.13028

Aigu, Y., Laperche, A., Mendes, J., Lariagon, C., Guichard, S., Gravot, A. and Manzanares-Dauleux, M.J. (2018) Nitrogen supply exerts a major/minor switch between two QTLs controlling *Plasmodiophora brassicae* spore content in rapeseed. *Plant Pathol.*, 67, 1574–1581. https://doi.org/10.1111/ppa.12867

Aigu, Y., Cao, T., Strelkov, I.S., Manolii, V.P., Lemoine, J., Manzanares-Dauleux, M.J., Strelkov, S.E. and Gravot, A. (2020) Identification of winter and spring *Brassica napus* genotypes with partial resistance to Canadian isolates of *Plasmodiophora brassicae*. *Can. J. Plant Pathol.*, 42, 538–546. https://doi.org/10.1080/07060661.2020.1723870

Aigu, Y., Daval, S., Gazengel, K., Marnet, N., Lariagon, C., Laperche, A., Legeai, F., Manzanares-Dauleux, M.J. and Gravot, A. (2022) Multi-omic investigation of low-nitrogen conditional resistance to clubroot reveals *Brassica napus* genes involved in nitrate assimilation. *Front. Plant Sci.*, 11. https://doi.org/10.3389/fpls.2022.790563

Alexa, A., Rahnenführer, J. and Lengauer, T. (2006) Improved scoring of functional groups from gene expression data by decorrelating GO graph structure. *Bioinformatics*, 22, 1600–1607. https://doi.org/10.1093/bioinformatics/btl140

Alexander, R.D., Castillejo-Pons, P., Alsaif, O., Stahl, Y., Seale, M. and Morris, P.C. (2020) The conserved plant PM19 protein functions as an osmosensor and regulator of germination. *bioRxiv*, https://doi.org/10.1101/2020.08.10.244889.

Atkinson, N.J. and Urwin, P.E. (2012) The interaction of plant biotic and abiotic stresses: from genes to the field. *J. Exp. Bot.*, 63, 3523–3543. https://doi.org/10.1093/jxb/ers100

Baldini, M., Giovanardi, R., Tahmasebi Enferadi, S. and Vannozzi, G. (2002) Effects of water regime on fatty acid accumulation and final fatty acid composition in the oil of standard and high oleic sunflower hybrids. *Ital. J. Agron.*, 6, 119–126.

Barcala, M., García, A., Cubas, P., Almoguera, C., Jordano, J., Fenoll, C. and Escobar, C. (2008) Distinct heat-shock element arrangements that mediate the heat shock, but not the late-embryogenesis induction of
small heat-shock proteins, correlate with promoter activation in root-knot nematode feeding cells. *Plant Mol. Biol.*, 66, 151-164. https://doi.org/10.1007/s11103-007-9259-3

Barnabas, B., Jager, K. and Feher, A. (2008) The effect of drought and heat stress on reproductive processes in cereals. *Plant Cell Env.*, 31, 11-38. https://doi.org/10.1111/j.1365-3040.2007.01727.x

Barthet, V.J. and Daun, J.K. (2011) Seed morphology, composition, and quality. In *Canola: chemistry, production, processing, and utilization* (Daun, J.K., Michael Eskin, N.A. and Hickling D., eds). Canola, AOCS Press, pp. 119-162. https://doi.org/10.1016/B978-0-9818936-5-5.50009-7

Basbouss-Serhal, I., Soubigou-Taconnat, L., Bailly, C. and Leymarie, J. (2015) Germination potential of dormant and nondormant Arabidopsis seeds is driven by distinct recruitment of messenger RNAs to poly-somes. *Plant Physiol.*, 168, 1049-1065. https://doi.org/10.1104/pp.15.00510

Baud, S., Boutin, J-P., Miquel, M., Lepiniec, L. and Rochat, C. (2002) An integrated overview of seed development in *Arabidopsis thaliana* ecotype WS. *Plant Physiol. Biochem.*, 40, 151-160. https://doi.org/10.1016/S0981-9428(01)01350-X

Baud, S., Dubreucq, B., Miquel, M., Rochat, C. and Lepiniec, L. (2008) Storage reserve accumulation in Arabidopsis: metabolic and developmental control of seed filling. In *The arabidopsis book*, 6, e0113. https://doi.org/10.1199/tab.0113

Becker, M.G., Walker, P.L., Pulgar-Vidal, N.C. and Belmonte, M.F. (2017) SeqEnrich: A tool to predict transcription factor networks from co-expressed Arabidopsis and *Brassica napus* gene sets. *PLoS One*, 12, e0178256. https://doi.org/10.1371/journal.pone.0178256

Bensmihen, S., Rippa, S., Lambert, G., Jublot, D., Pautot, V., Granier, F., Giraudat, J. and Parcy, F. (2002) The homologous ABI5 and EEL transcription factors function antagonistically to fine-tune gene expression during late embryogenesis. *Plant Cell*, 14, 1391-1403. https://doi.org/10.1105/tpc.000869

Bewley, J.D., Bradford, K.J., Hilhorst, H.W.M. and Nonogaki, H. (2013) *Seeds: Physiology of development, germination and dormancy*. 3rd Edition, Springer, New York. https://doi.org/10.1007/978-1-4614-4693-4

Bianchetti, G., Clouet, V., Legeai, F., Gazengel, K., Baron, C., Carrillo, A., Manzanares-Dauleux, M.J., Buitink, J. and Nesi, N. (2021a) RNA sequencing data for responses to drought stress and/or clubroot infection in developing seeds of *Brassica napus*. *Data Br.*, 38, 107392. https://doi.org/10.1016/j.dib.2021.107392

Bianchetti, G., Baron, C., Carrillo, A., Berardocco, S., Marnet, N., Wagner, M-H., Demilly, D., Ducournau, S., Manzanares-Dauleux, M.J., Le Cahérec, F., Buitink, J. and Nesi, N. (2021b) Dataset for the metabolic and physiological characterization of seeds from oilseed rape (*Brassica napus* L.) plants grown under single or combined effects of drought and clubroot pathogen *Plasmodiophora brassicae*. *Data Br.*, 37, 107247. https://doi.org/10.1016/j.dib.2021.107247

Bizouerne, E., Buitink, J., Ly Vu, B., Ly Vu, J., Esteban, E., Pasha, A., Provart, N., Verdier, J. and Leprince, O. (2021) Gene co-expression analysis of tomato seed maturation reveals tissue-specific regulatory networks and hubs associated with the acquisition of desiccation tolerance and seed vigour. BMC Plant Biol., 21, 124. https://doi.org/10.1186/s12870-021-02889-8

Bogamuwa, S. and Jang, J.C. (2016) Plant tandem CCCH zinc finger proteins interact with ABA, drought, and stress response regulators in processing-bodies and stress granules. *PLoS One*, 11, e0151574. https://doi.org/10.1371/journal.pone.0151574

Borisjuk, L., Neuberger, T., Schwender, J., Heinzel, N., Sunderhaus, S., Fuchs, J., Hay, J.O., Tschiersch, H., Braun, H.-P., Denolf, P., et al. (2013) Seed architecture shapes embryo metabolism in oilseed rape. *Plant Cell*, 25, 1625-1640. https://doi.org/10.1105/tpc.113.111740

Bouchereau, A., Clossais-Besnard, N., Bensaoud, A., Leport, L. and Renard, M. (1996) Water stress effects on rapeseed quality. *Eur. J. Agron.*, 5, 19-30. https://doi.org/10.1016/S1161-0301(96)02005-9

Brunel-Muguet, S., D'Hooghe, P., Bataillé, M-P., Larré, C., Kim, T-H., Trouverie, J., Avice, J-C., Etienne, P. and Dürr, C. (2015) Heat stress during seed filling interferes with sulfur restriction on grain composition and seed germination in oilseed rape (*Brassica napus* L.). *Front. Plant Sci.*, 6, 213. https://doi.org/10.3389/fpls.2015.01236

Chalhoub, B., Denoeud, F., Liu, S., Parkin, I.A.P., Tang, H., et al. (2014) Early allopolyploid evolution in the post-Neolithic Brassica napus oilseed genome. *Science*, 345, 950-953. https://doi.org/10.1126/science.1253435

Champolivier, L. and Merrien, A. (1996) Effects of water stress applied at different growth stages to *Brassica napus* L. var. oleifera on yield, yield components and seed quality. *Eur. J. Agron.*, 5, 153-160. https://doi.org/10.1016/S1161-0301(96)02004-7

Chern, M.S., Eiben, H.G. and Bustos, M.M. (1996) The developmentally regulated bZIP factor ROM1 modulates transcription from lectin and storage protein genes in bean embryos. *Plant J.*, 10, 135-148. https://doi.org/10.1046/j.1365-313X.1996.10010135.x

Chia, T.Y.P., Pike, M.J. and Rawsthorne, S. (2005) Storage oil breakdown during embryo development of Brassica napus (L.). *J. Exp. Bot.*, 56, 1285-1296. https://doi.org/10.1093/jxb/eri129

Ciaghi, S., Schwelm, A. and Neuhauser, S. (2019) Transcriptomic response in symptomless roots of clubroot infected kohlrabi (*Brassica oleracea* var. gongyloides) mirrors resistant plants. BMC Plant Biol., 19, 288. https://doi.org/10.1186/s12870-019-1902-z

Creus, C.M., Sueldo, R.J. and Barassi, C.A. (2004) Water relations and yield in Azospirillum-inoculated wheat exposed to drought in the field. *Can. J. Bot.*, 82, 273-281. https://doi.org/10.1139/b03-119

Ding, Y., Fromm, M. and Avramova, Z. (2012) Multiple exposures to drought "train" transcriptional responses in Arabidopsis. *Nat. Commun.*, 3, 740. https://doi.org/10.1038/ncomms1732

Dixon, G.R. (2009) The occurrence and economic impact of *Plasmodiophora brassicae* and clubroot disease. *J. Plant Growth Regul.*, 28, 194-202. https://doi.org/10.1007/s00344-009-9090-y

Elferjani, R. and Soolanayakanahally, R. (2018) Canola responses to drought, heat, and combined stress: shared and specific effects on carbon assimilation, seed yield, and oil composition. *Front. Plant Sci.*, 9, 1224. https://doi.org/10.3389/fpls.2018.01224

Ewels, P.A., Peltzer, A., Fillinger, S., Patel, H., Alneberg, J., Wilm, A., Garcia, M.U., Di Tommaso, P. and Nahnsen, S. (2020) The nf-core framework for community-curated bioinformatics pipelines. *Nat. Biotechnol.*, 38, 276-278. https://doi.org/10.1038/s41587-020-0439-x

Faroog, M., Hussain, M. and Siddique, K.H.M. (2014) Drought stress in wheat during flowering and grain-filling periods. *Critic. Rev. Plant Sci.*, 33, 331-349. https://doi.org/10.1080/07352689.2014.875291

Fatihi, A., Boulard, C., Bouyer, D. Baud, S., Dubreucq, B. and Lepiniec, L. (2016) Deciphering and modifying LAFL transcriptional regulatory network in seed for improving yield and quality of storage compounds. *Plant Sci.*, 250, 198-204. https://doi.org/10.1016/j.plantsci.2016.06.013

Ghazalpour, A., Doss, S., Zhang, B., Wang, S., Plaisier, C., Castellanos, R., et al. (2006) Integrating genetic and network analysis to characterize genes related to mouse weight. *PLoS Genet.*, 2, 1182-1192. https://doi.org/10.1371/journal.pgen.0020130

Goldberg, R.B., de Paiva, G. and Yadegari, R. (1994) Plant embryogenesis: Zygote to seed. *Science*, 266, 605-614. DOI: https://doi.org/10.1126/science.266.5185.605

Gravot, A., Lemarié, S., Richard, G., Lime, T., Lariagon, C. and Manzanares-Dauleux, M.J. (2016) Flooding affects the post-invasive development of *Plasmodiophora brassicae* in *Arabidopsis* roots during the secondary phase of infection. *Plant Pathol.*, 65, 1153-1160. https://doi.org/10.1111/ppa.12487

Hatzig, S.V., Nuppenau, J-N., Snowdon, R.J. and Schießl, S.V. (2018) Drought stress has transgenerational effects on seeds and seedlings in winter oilseed rape (*Brassica napus* L.). BMC Plant Biol., 18, 297. https://doi.org/10.1186/s12870-018-1531-y
Havé, M., Marmagne, A., Chardon, F. and Masclaux-Daubresse, C. (2017) Nitrogen remobilization during leaf senescence: lessons from Arabidopsis to crops. *J. Exp. Bot.*, 68, 2513-2529. https://doi.org/10.1093/jxb/erw365

Hilker, M., Schwachtje, J., Baier, M., Balazadeh, S., Bäurle, I., Geiselhardt, S., et al. (2016) Priming and memory of stress responses in organisms lacking a nervous system. *Biol. Rev. Camb. Phil. Soc.*, 91, 1118-1133. https://doi.org/10.1111/brv.12215

Kim, D., Paggi, J.M., Park, C., Bennett, C. and Salzberg S.L. (2019) Graph-based genome alignment and genotyping with HISAT2 and HISAT-genotype Graph-based genome alignment and genotyping with HISAT2 and HISAT-genotype. *Nat. Biotechnol.*, 37, 907-915. https://doi.org/10.1038/s41587-019-0201-4

Kotak, S., Vierling, E., Bäumlein, H. and von Koskull-Döring, P. (2007) A novel transcriptional cascade regulating expression of heat stress proteins during seed development of Arabidopsis. *Plant Cell*, 19, 182-195. https://doi.org/10.1105/tpc.106.048165

Lancashire, P.D., Bleiholder, H., Boom, T.V.D., Langelüddeke, P., Stauss, R., Weber, E. and Witzenberger, A. (1991) A uniform decimal code for growth stages of crops and weeds. *Ann. Appl. Biol.*, 119, 561-601. https://doi.org/10.1111/j.1744-7348.1991.tb04895.x

Langfelder, P. and Horvath, S. (2008) WGCNA: An R package for weighted correlation network analysis. *BMC Bioinformatics*, 9, 559. https://doi.org/10.1186/1471-2105-9-559

Laperche, A., Aigu, Y., Jubault, M., Ollier, M., Guichard, S., Glory, P., Strelkov, S.E., Gravot, A. and Manzanares-Dauleux, M.J. (2017) Clubroot resistance QTL are modulated by nitrogen input in *Brassica napus*. *Theor. Appl. Genet.*, 130, 669-684. https://doi.org/10.1007/s00122-016-2842-8

Le, B.H., Cheng, C., Bui, A.Q., Wagmaister, J.A., Henry, K.F., Pelletier, J., Kwong, L., Belmonte, M., Kirkbride, R., Horvath, S., Drews, G.N., Fischer, R.L., Okamuro, J.K., Harada, J.J. and Goldberg, R.B. (2010) Global analysis of gene activity during Arabidopsis seed development and identification of seed-specific transcription factors. *Proc. Natl. Acad. Sci. USA*, 107, 8063-8070. https://doi.org/10.1073/pnas.1003530107

Leprince, O., Pellizzaro, A., Berrir, S. and Buitink, J. (2017) Late seed maturation: drying without dying. *J. Exp. Bot.*, 68, 827-841. https://doi.org/10.1093/jxb/erw363

Lê, S., Josse, J., and Husson, F. (2008). FactoMineR: an R package for multivariate analysis. *Journal of statistical software*, 25, 1-18. Li, Z. and Thomas, T.L. (1998) PEI1, an embryo-specific zinc finger protein gene required for heart-stage embryo formation in Arabidopsis. *Plant Cell*, 10, 383-398. https://doi.org/10.1105/tpc.10.3.383

Li, Y., Beisson, F., Pollard, M. and Ohlrogge, J. (2006) Oil content of Arabidopsis seeds: the influence of seed anatomy, light and plant-to-plant variation. *Phytochemistry*, 67, 904-915. https://doi.org/10.1016/j.phytochem.2006.02.015

Li, Y., Ye, W., Wang, M. and Yan, X. (2009) Climate change and drought. A risk assessment of crop-yield impacts. *Clim. Res.*, 39, 31-46. https://doi.org/10.3354/cr00797

Li, Y., Liu, W., Zhong, H., Zhang, H.-L. and Xia, Y. (2019) Redox-sensitive bZIP68 plays a role in balancing stress tolerance with growth in Arabidopsis. *Plant J.*, 100, 768-783. https://doi.org/10.1111/tpj.14476

Liao, Y., Smyth, G.K. and Shi, W. (2014) featureCounts: an efficient general purpose program for assigning sequence reads to genomic features. *Bioinformatics*, 30, 923-930. https://doi.org/10.1093/bioinformatics/btt656

Liégard, B., Baillet, V., Etcheverry, M., Joseph, E., Lariagon, C., Lemoine, J., Evrard, A., Colot, V., Gravot, A., Manzanares-Dauleux, M.J. and Jubault, M. (2019) Quantitative resistance to clubroot infection mediated by transgenerational epigenetic variation in Arabidopsis. *New Phytol.*, 222, 468-479. https://doi.org/10.1111/nph.15579

Liu, Y.-H., Offler, C.E. and Ruan, Y.-L. (2013) Regulation of fruit and seed response to heat and drought by sugars as nutrients and signals. *Front. Plant Sci.*, 4, 282. https://doi.org/10.3389/fpls.2013.00282

Locascio, A., Roig-Villanova, I., Bernardi, J. and Varotto, S. (2014) Current perspectives on the hormonal control of seed development in Arabidopsis and maize: a focus on auxin. *Front. Plant Sci.*, 5, 412. https://doi.org/10.3389/fpls.2014.00412

Ludwig-Müller, J., Prinsen, E., Rolfe, S.A. and Scholes, J.D. (2009) Metabolism and plant hormone action during clubroot disease. *J. Plant Growth Regul.*, 28, 229-244. https://doi.org/10.1007/s00344-009-9089-4

Mittler, R. (2006) Abiotic stress, the field environment and stress combination. *Trends Plant Sci.*, 11, 15-19. https://doi.org/10.1016/j.tplants.2005.11.002

Manzanares-Dauleux, M.J., Divaret, I., Baron, F. and Thomas G. (2000a) Evaluation of French *Brassica oleracea* landraces for resistance to *Plasmodiophora brassicae*. *Euphytica*, 113, 211-218. https://doi.org/10.1023/A:1003997421340

Manzanares-Dauleux, M.J., Delourme, R., Baron, F. and Thomas, G. (2000b) Mapping of one major gene and of QTLs involved in resistance to clubroot in *Brassica napus*. *Theor. Appl. Genet.*, 101, 885-891. https://doi.org/10.1007/s001220051557

Meinke, D.W., Franzmann, L.H., Nickle, T.C. and Yeung, E.C. (1994) Leafy cotyledon mutants of Arabidopsis. *Plant Cell*, 6, 1049-1064. https://doi.org/10.1105/tpc.6.8.1049

Nambara, E., Naito, S. and McCourt, P. (1992) A mutant of Arabidopsis which is defective in seed development and storage protein accumulation is a new abi3 allele. *Plant J.*, 2, 435-441. https://doi.org/10.1111/j.1365-313X.1992.00435.x

Nesi, N., Delourme, R., Brégeon, M., Falentin, C. and Renard, M. (2008) Genetic and molecular approaches to improve nutritional value of *Brassica napus* L. seed. *C. R. Biol.*, 331, 763-771. https://doi.org/10.1016/j.crvi.2008.07.018

Nishizawa, A., Yabuta, Y., Yoshida, E., Maruta, T., Yoshimura, K. and Shigeoka, S. (2006) Arabidopsis heat shock transcription factor A2 as a key regulator in response to several types of environmental stress. *Plant J.*, 48, 535-547. https://doi.org/10.1111/j.1365-313X.2006.02889.x

Oksanen, J., Blanchet, F. G., Friendly, M., Kindt, R., Legendre, P., McGlinn, D., Minchin, P. R., O'Hara, R. B., Simpson, G. L., Solymos, P., Stevens, M. H. H., Szoecs, E., & Wagner, H. (2015) Vegan: community ecology package. R package vegan, vers. 2.5-7. Retrieved on September 21th, 2021, from http://CRAN.Rproject.org/package=vegan2015

Pageau, D., Lajeunesse, J. and Lafond, J. (2006) Impact of clubroot (*Plasmodiophora brassicae*) on the yield and quality of canola. *Can. J. Plant Pathol.*, 28, 137-143.

Pandey, P., Irulappan, V., Bagavathiannan, M.V. and Senthil-Kumar, M. (2017) Impact of combined abiotic and biotic stresses on plant growth and avenues for crop improvement by exploiting physio-morphological traits. *Front. Plant Sci.*, 8, 537. https://doi.org/10.3389/fpls.2017.00537

Park, C.J. and Seo, Y.S. (2015) Heat shock proteins: a review of the molecular chaperones for plant immunity. *Plant Pathol. J.*, 31, 323-333. https://doi.org/10.5423/PPJ.RW.08.2015.0150

Rae, G.M., Uversky, V.N., David, K. and Wood, M. (2014) DRM1 and DRM2 expression regulation: potential role of splice variants in response to stress and environmental factors in Arabidopsis. *Mol. Gen. Genomics*, 289, 317-332. https://doi.org/10.1007/s00438-013-0804-2

Rajjou, L. and Debeaujon, I. (2008) Seed longevity: Survival and maintenance of high germination ability of dry seeds. *C. R. Biol.*, 331, 796-805. https://doi.org/10.1016/j.crvi.2008.07.021

Ramegowda, V. and Senthil-Kumar, M. (2015) The interactive effects of simultaneous biotic and abiotic stresses on plants: mechanistic
understanding from drought and pathogen combination. *J. Plant Physiol.*, 176, 47–54. https://doi.org/10.1016/j.jplph.2014.11.008

Rasmussen, S., Barah, P., Suarez-Rodriguez, M.C., Bressendorff, S., Friis, P., Costantino, P., Bones, A.M., Nielsen, H.B. and Mundy, J. (2013) Transcriptome responses to combinations of stresses in Arabidopsis, *Plant Physiol.*, 161, 1783–1794. https://doi.org/10.1104/pp.112.210773

Righetti, K., Vu, J.L., Pelletier, S., Vu, B.L., Glaab, E., Lalanne, D., Pasha, A., Patel, R.V., Provart, N.J., Verdier, J., et al. (2015) Inference of longevity-related genes from a robust coexpression network of seed maturation identifies regulators linking seed storability to biotic defense-related pathways. *Plant Cell*, 27, 2692–2708. https://doi.org/10.1105/tpc.15.00632

Rousseau-Gueutin, M., Belser, C., Da Silva, C., Richard, G., Istace, B., Cruaud, C., Falentin, C. et al. (2020) Long-read assembly of the *Brassica napus* reference genome Darmor-bzh. *GigaScience*, 9, giaa137.

Roy, S., Saxena, S., Sinha, A. and Nandi, A.K. (2020) DORMANCY/AUXIN ASSOCIATED FAMILY PROTEIN 2 of Arabidopsis thaliana is a negative regulator of local and systemic acquired resistance. *J. Plant Res.*, 133, 409–417. https://doi.org/10.1007/s10265-020-01183-2

Santos-Mendoza, M., Dubreucq, B., Baud, S., Parcy, F., Caboche, M. and Lepiniec, L. (2008) Deciphering gene regulatory networks that control seed development and maturation in Arabidopsis. *Plant J.*, 54, 608–620. https://doi.org/10.1111/j.1365-313X.2008.03461.x

Sehgal, A., Sita, K., Siddique, K.H.M., Kumar, R., Bhogireddy, S., Varshney, R.K., HanumanthaRao, B., Nair, R.M., Prasad, P.V.V. and Nayyar, H. (2018) Drought or/and heat-stress effects on seed filling in food crops: impacts on functional biochemistry, seed yields, and nutritional quality. *Front. Plant Sci.*, 9, 1705. https://doi.org/10.3389/fpls.2018.01705

Stagnari, F., Galieni, A. and Pisante, M. (2016) Drought stress effects on crop quality. In *Water Stress and Crop Plants: A Sustainable Approach* (Ahmad, P., eds). Chichester, UK: Wiley, pp. 375–392.

Strehlow, B., de Mol, F. and Struck, C. (2015) Risk potential of Clubroot disease on winter oilseed rape. *Plant Dis.*, 99, 667–675. https://doi.org/10.1094/PDIS-05-14-0482-RE

Strelkov, S.E., Hwang, S.F., Manolii, V.P., Turnbull, G., Fredua-Agyeman, R., Hollman, K. and Kaus, S. (2021) Characterization of clubroot (*Plasmodiophora brassicae*) from canola (*Brassica napus*) in the Peace Country of Alberta, Canada. *Can. J. Plant Pathol.*, 43, 155–161. https://doi.org/10.1080/07060661.2020.1776931

Summanwar, A., Farid, M., Basu, U., Kav, N. and Rahman, H. (2021) Comparative transcriptome analysis of canola carrying clubroot resistance from ‘Mendel’ or Rutabaga and the development of molecular markers. *Physiol. Mol. Plant Pathol.*, 114, 101640. https://doi.org/10.1016/j.pmpp.2021.101640

Sun, K., Hunt, K. and Hauser, B.A. (2004) Ovule abortion in Arabidopsis triggered by stress. *Plant Physiol.*, 135, 2358–2367. https://doi.org/10.1104/pp.104.043091

Suzuki, N., Rivero, R.M., Shulaev, V., Blumwald, E. and Mittler, R. (2014) Abiotic and biotic stress combinations. *New Phytol.*, 203, 32–43. https://doi.org/10.1111/nph.12797

Tellmann G. (2006) The E-Method: a highly accurate technique for gene-expression analysis. *Nat. Methods*, 3. https://doi.org/10.1038/nmeth894

ul Haq, S., Khan, A., Ali, M., Khattak, A.M., Gai, W.X., Zhang, H.X., Wei, A.M. and Gong, Z.H. (2019) Heat shock proteins: dynamic bio-molecules to counter plant biotic and abiotic stresses. *Int. J. Mol. Sci.*, 20, 5321. https://doi.org/10.3390/ijms20215321

Vandesompele, J., De Preter, K., Pattyn, F., Poppe, B., Van Roy, N., De Paepe, A. and Speleman, F. (2002) Accurate normalization of real-time quantitative RT-PCR data by geometric averaging of multiple internal control genes. *Genome Biol.*, 3, research0034.1. https://doi.org/10.1186/gb-2002-3-7-research0034

Verdier, J. and Thompson, R.D. (2008) Transcriptional regulation of storage protein synthesis during dicotyledon seed filling. *Plant Cell Physiol.*, 49, 1263–1271. https://doi.org/10.1093/pcp/pcn116

Wickham, H. (2009) *Ggplot2* elegant graphics for data analysis 2nd edn. Springer, New York.

Wolfe, C.J., Kohane, I.S. and Butte, A.J. (2005) Systematic survey reveals general applicability of “guilt-by-association” within gene coexpression networks. *BMC Bioinformatics*, 6, 227. https://doi.org/10.1186/1471-2105-6-227

Xu, P., Chen, F., Mannas, J.P., Feldman, T., Sumner, L.W. and Roossinck, M.J. (2008) Virus infection improves drought tolerance. *New Phytol.*, 180, 911–921. https://doi.org/10.1111/j.1469-8137.2008.02627.x

Yao, L., Cheng, X., Gu, Z., Huang, W., Li, S., Wang, L., Wang, Y.-F., Xu, P., Ma, H. and Ge, X. (2018) The AWPM-19 family protein OsPM1 mediates abscisic acid influx and drought response in rice. *Plant Cell*, 30, 1258–1276. https://doi.org/10.1105/tpc.17.00770

Zhang, B. and Horvath, S. (2005) A general framework for weighted gene co-expression network analysis. *Stat. Appl. Genet. Mol. Biol.*, 4, Article17. https://doi.org/10.2202/1544-6115.1128

## SUPPORTING INFORMATION

Additional supporting information can be found online in the Supporting Information section at the end of this article.

**How to cite this article:** Bianchetti, G., Clouet, V., Legeai, F., Baron, C., Gazengel, K., Vu, B.L. et al. (2024) Identification of transcriptional modules linked to the drought response of *Brassica napus* during seed development and their mitigation by early biotic stress. *Physiologia Plantarum*, 176(1), e14130. Available from: https://doi.org/10.1111/ppl.14130