Genome-Wide Identiﬁcation and Characterization of Oil-Body-Membrane Proteins in Polyploid Crop Brassica napus

: Oil-body-membrane proteins (OBMPs) are essential structural molecules of oil bodies and also versatile metabolic enzymes involved in multiple cellular processes such as lipid metabolism, hormone signaling and stress responses. However, the global landscape for OBMP genes in oil crops is still lacking. Here, we performed genome-wide identiﬁcation and characterization of OBMP genes in polyploid crop Brassica napus . B. napus contains up to 88 BnaOBMP genes including 53 oleosins, 20 caleosins and 15 steroleosins. Both whole-genome and tandem duplications have contributed to the expansion of the BnaOBMP gene family. These BnaOBMP genes have extensive sequence polymorphisms, and some harbor strong selection signatures. Various cis-acting regulatory elements involved in plant growth, phytohormones and abiotic and biotic stress responses are detected in their promoters. BnaOBMPs exhibit differential expression at various developmental stages from diverse tissues. Importantly, some BnaOBMP genes display spatiotemporal patterns of seed-speciﬁc expression, which could be orchestrated by transcriptional factors such as EEL, GATA3, HAT2, SMZ, DOF5.6 and APL. Altogether, our data lay the foundations for studying the regulatory mechanism of the seed oil storage process and provide candidate genes and alleles for the genetic improvement and breeding of rapeseed with high seed oil content.


Introduction
Brassica napus L. (known as rapeseed, oilseed rape and rapa), a relatively recent allotetraploid formed from hybridization between Brassica rapa and Brassica oleracea, is an important oilseed crop and widely cultivated in the world [1,2]. Rapeseed provides more than 13% of the world's supply of vegetable oil [3]. Thus, increasing seed oil content (SOC) is one of the most important breeding objectives for rapeseed with the greatest economic significance [4]. Over the past decades, many QTLs (quantitative trait loci) controlling SOC have been identified in B. napus, and significant progress has been achieved in breeding rapeseed varieties with high oil content [5,6]. However, the molecular mechanism underlying the regulation of seed oil metabolism, particularly the storage process, is still poorly understood in B. napus.
Lipids, mainly triacylglycerols (TAGs), are the main reserves in B. napus seeds and stored in oil bodies (OBs; known as lipid droplets, oil droplets or oleosomes) [7]. OBs are specialized organelles acting as neutral lipid storage compartments, which are widespread in both prokaryotic and eukaryotic cells [8][9][10][11]. The formation of OBs begins on the endoplasmic reticulum (ER), where TAGs are synthesized in maturing seeds, and then bud off from the ER membrane as droplets consisting of an oil core and a hydrophilic surface with a phospholipid monolayer [12]. The structure and size of OBs are highly dynamic during seed development, and various oil-body-membrane proteins (OBMPs) wrap around the surface to stabilize the OBs and exert metabolic functions [13][14][15]. Previous studies showed that overexpression of some OBMP genes could regulate oil body size, as well as seed size and weight, to influence oil accumulation in seeds [16][17][18].
The OBMP gene family is composed of three main classes including oleosin, caleosin and steroleosin [19]. Oleosin is a relatively small protein of 15~26 kilodaltons (kDa), with a conserved central domain of approximately 70 uninterrupted nonpolar residues, which could form a hydrophobic hairpin [20,21]. The two arms of this hairpin are linked by a proline (Pro) loop (called Pro knot) with three conserved Pro residues and one Ser residue. The N-and C-terminal peptides of oleosin could form amphipathic α-helical structures interacting horizontally with the phospholipid monolayer, acting as a receptor for the binding of metabolic enzymes or regulatory proteins. Oleosin genes are widespread from green algae to higher plants [20]. Six oleosin lineages have been recognized including primitive (P), universal (U), tapetum (T), mesocarp (M), and seed low-molecular-weight (SL) and high-molecular-weight (SH) oleosins [22][23][24]. The Arabidopsis thaliana genome has 17 oleosin genes of which OLE1 and OLE2 are the most abundant in seeds. Double oleosin mutant (i.e., ole1 and ole2) seeds have enlarged OBs and exhibit severely defective germination [25]. In addition, the seeds of the oleosin mutants are sensitive to freezing stress [26], indicating that oleosins function in enhancing plant survival during winter.
Like oleosins, caleosins are also structural proteins of OBs, consisting of a conserved hydrophobic central domain with a proline knot motif. Particularly, they contain an N-terminal hydrophilic domain with a single Ca 2+ -binding site, EF-hand motif, and possess an enzyme activity as peroxygenase [27]. Caleosin genes have been found ubiquitously in plants [28]. Among the eight caleosins identified in A. thaliana, CLO1 is preferentially expressed in developing seeds, and clo1 mutant seeds exhibit distorted vacuole morphology and a significant delay in the storage lipid degradation [29]. In addition, caleosins display calcium-binding properties and play crucial roles in many aspects of growth and development such as cell division, photosynthesis, polarity formation, apoptosis and stress resistance [30,31]. Steroleosins are known as sterol dehydrogenases due to the fact of their sterol-binding and sterol-coupling dehydrogenase activity [9,32,33], and they are reported to have functions related to sterol-regulated signal transduction, seed maturation and germination [33]. Unlike oleosin and caleosin, steroleosin has two main structural domains including an N-terminal hydrophobic domain, which contains a conserved proline knob instead of the proline knot, and a C-terminal domain exhibiting an NADP(H)-binding subdomain and a sterol-binding subdomain [9,32]. Various steroleosin genes have been discovered in A. thaliana [34], Pinus massoniana [35] and Sesamum indicum [36].
Although OBMP family genes have been studied in several plant species, the comprehensive identification and characterization of this gene family from important oil seed corps is lacking. In the present study, we performed a genome-wide identification of OBMP family genes in B. napus and 54 other genomes from green algae to angiosperms. The gene structures, phylogenetic relationships, cis-acting regulatory elements in promoters, sequence polymorphisms and expression patterns of BnaOBMPs were analyzed. We also identified transcriptional factors that potentially regulate the expression of BnaOBMP genes during seed development. Our results help to further elucidate the molecular mechanism of the seed oil storage process and provide new targets for the selection of rapeseed with high SOC.

Identification, Properties and Genomic Distributions of BnaOBMP Genes
To identify all the members of the OBMP gene family in rapeseed, we used protein sequences of 33 A. thaliana OBMP genes, including 17 oleosins, 8 caleosins and 8 steroleosins (Table S1), as queries to search against the protein dataset of B. napus var. ZS11 with BLASTP setting the E-value at 1 × 10 −5 . The peptides of putative BnaOBMPs with the best hits on A. thaliana oleosin, caleosin or steroleosin proteins were further used to predict functional domains in the Pfam database to confirm the presence of the oleosin domain (PF01277), the caleosin domain (PF05042) or the steroleosin domain (PF00106). Finally, a total of 88 BnaOBMP genes were identified in rapeseed of which 43 and 45 genes originated from the A and C subgenomes, respectively (Table 1). Based on the functional domains of each BnaOBMP contained, these genes could be divided into 53 oleosins, 20 caleosins and 15 steroleosins ( Table 1). The 53 oleosin genes could be further classified into four subfamilies of 27 T, 8 SL, 9 SH and 9 U oleosins in rapeseed. Meanwhile, we also found 44 and 47 OBMP genes in its two ancestors of B. rapa var. Z1 and B. oleracea var. HDEM, respectively (Table S1). Additionally, we identified OBMP genes in fifty-one other genomes covering green algae, mosses, gymnosperms and angiosperms and found that B. napus contained the greatest number of OBMP family genes ( Figure S1).

Phylogeny, Gene Structures and Crucial Motifs of BnaOBMPs
To explore the molecular evolution of the OBMP gene family in B. napus, a total of 212 OBMP genes from B. napus, B. rapa, B. oleracea and A. thaliana were used to construct an unrooted phylogenetic tree using the maximum likelihood (ML) method. According to the phylogenetic relationships of these OBMP genes, they could be divided into three independent classes, which is consistent with classification by the functional domains (i.e., oleosin, caleosin and steroleosin) contained ( Figure 2a). Moreover, the oleosin class could be further categorized into four subclasses, corresponding to the SH, SL, T and U subfamilies (Figure 2a). Based on the gene information of the genome available in the BnPIR database (http: //cbi.hzau.edu.cn/bnapus/ (accessed on 8 January 2022), we performed a gene structure analysis on the BnaOBMPs. The adjacent BnaOBMP genes in the phylogenetic tree, which were derived from the homologous A and C subgenomes of rapeseed, respectively, exhibited similar gene structures (Figure 2b). For the BnaOBMP genes of the caleosin or steroleosin subfamilies, their exon features of order, length and number were largely conserved among different members (Figure 2b). On the contrary, the BnaOBMP genes from the oleosin subclasses showed varied gene structures. For examples, the SH and SL oleosins and most of the T oleosins contained two exons and one intron, while all of the U oleosins and three T oleosin genes (i.e., BnaOBMP.T6, BnaOBMP.T12 and BnaOBMP.T25) had only one exon, and the BnaOBMP.T8 and BnaOBMP.T21 oleosin genes had four and three exons, respectively (Figure 2b). Additionally, the conserved motifs of BnaOBMP proteins were identified using MEME, and three distinct conserved motifs (i.e., Motifs 1, 2 and 3) were found (Figure 2b). Motif 1 was a proline knot with four invariable residues of the proline knot sequence (-PX 5 SPX 3 P-), Motif 2 was a Ca 2+ -binding site, and Motif 3 was an NADPH-binding domain (Figure 2c). The sequence homology of BnaOBMP proteins was determined through multiple sequence alignment. In B. napus, the protein sequences of caleosins consisted of an H-form insertion, a Ca 2+ -binding motif and a proline knot ( Figure S2). All the oleosins comprised an N-arm of a Pro loop, a proline knot and a C-arm of a Pro loop ( Figure S3), and the steroleosins harbored one proline knob, one NADPH-binding motif and a sterol-binding domain in the C-terminal ( Figure S4).

Synteny and Gene Duplication Modes of BnaOBMP Genes
Gene duplications are considered to be one of the major driving forces in the evolution of genomes and expansions of gene families [37]. Whole-genome duplication, segmental duplication and tandem duplication are the major causes of gene family expansion in plants [38]. Based on the synteny and collinearity of B. napus genome from all-against-all pairwise alignment of all the proteins, we identified the duplicated events for BnaOBMP genes. Approximately 82.9% (73/88) of the BnaOBMP genes were found associated with at least one collinear gene pair (Figure 3a). A total of 165 gene pairs with 72 (accounting for 81.82%) genes were identified as WGD/segmental duplication (Figure 3a,b). In addition, 15 (accounting for 17.04%) tandem duplicated genes were identified within 6 tandem duplicated gene clusters, and they all belonged to the T oleosin subfamily (Figure 3a,b). To better understand different selective constraints on the BnaOBMP gene family, the nonsynonymous substitution/synonymous substitution (Ka/Ks) ratios of the OBMP gene pairs between B. napus and A. thaliana were calculated ( Figure S5). Except for BnaOBMP.T5, BnaOBMP.T10, BnaOBMP.T15, BnaOBMP.T18, BnaOBMP.T24, BnaOBMP.SH9, BnaOBMP.SL7 and BnaOBMP.S13, all other orthologous OBMP gene pairs had Ka/Ks < 1. In addition, the Ka/Ks ratios of T oleosin genes were substantially higher than the other OBMP genes ( Figure S5).

Cis-Acting Regulatory Elements in the Promoters of BnaOBMP Genes
The cis-acting elements in the promoter region play vital roles in the regulation of gene expression and also function to coordinate responses to developmental and environmental cues in plants [39]. To identify the putative cis-acting elements of BnaOBMP genes, the genomic DNA sequences from the transcription initiation site (+1) to 2 kilo bases (Kb) upstream of the BnaOBMP genes were extracted for cis-acting elements profiling using PlantCARE. All the BnaOBMP gene promoters possessed various cis-acting regulatory elements ( Figure 4). Based on the functional annotation of these cis-acting elements, they could be classed into three main groups including plant growth and development, phytohormone responsive and abiotic and biotic stress responsive. Specifically, the BnaOBMP genes contained multiple phytohormone responsive elements, such as ABRE (abscisic acid-responsive element), AuxRE (auxin-responsive element), ERE (ethylene-responsive element), GARE (gibberellin-responsive element), MeJARE (MeJA-responsive element) and SARE (salicylic acid-responsive element), which suggests that the expression of BnaOBMP genes might be induced by different phytohormones. Among these phytohormone-responsive elements, the ABRE and MeJARE elements showed higher frequency, whereas AuxRE and GARE had a lower frequency in the promoters of BnaOBMP genes (Figure 4). Among the cisacting regulatory elements, LTR (light-responsive element) exhibited the most abundance in the promoters of BnaOBMP genes. In addition, diverse cis-acting regulatory elements involved in abiotic and biotic stress-responsive elements, such as ARE (anoxic-responsive element), DSRE (defense-and stress-responsive element), DIRE (drought-responsive element), DRE (damage-responsive element), LTRE (low-temperature-responsive element) and WRE (wound-responsive element), were also identified in the promoter regions of BnaOBMP genes.

Extensive Sequence Polymorphisms of the BnaOBMP Genes
DNA sequence polymorphisms in the gene can provide insights into the evolutionary forces acting on populations and species adapting to different environments [40]. Based on the genomic resequencing dataset of worldwide rapeseed accessions [41], we performed a variation analysis on the BnaOBMP genes to assess their sequence polymorphisms. The polymorphism sites of exons and introns in the BnaOBMP genes ranged from 2 (BnaOBMP.T15) to 86 (BnaOBMP.T26) and 2 (BnaOBMP.SH8) to 951 (BnaOBMP.S1), respectively (Table S2). Among the BnaOBMP genes, twelve genes contained 1 to 5 splicing variants, twenty-seven genes had one to three stop-gain variants, and four genes harbored one stop-loss variant, which indicates that these genes exhibited loss of function (Table S2). Among the BnaOBMP genes, BnaOBMP.U9 exhibited the highest coding sequence variation ratio of 0.167, while BnaOBMP.T15 had the lowest variation ratio of 0.003 (Table S2). The pi (π) values of the nucleotide diversity parameters extended from 0.000087 (BnaOBMP.T16) to 0.01528 (BnaOBMP.T9). The BnaOBMP genes of BnaOBMP.T9 (π = 0.01528), BnaOBMP.U9 (π = 0.00986), BnaOBMP.T14 (π = 0.00871), BnaOBMP.U7 (π = 0.00834), BnaOBMP.C4 (π = 0.00772) and BnaOBMP.SL3 (π = 0.00764) showed relatively higher nucleotide diversity (Figure 5a). In addition, the π values of BnaOBMP genes derived from the A and C subgenomes showed no significant difference (Student's t-test, p-value = 0.055). To study the population selection pressures for BnaOBMP genes, we conducted neutrality test statistics using Tajima's D-test method. Of all the BnaOBMP genes, approximately 27.3% (24/88) had positive Tajima's D values, while the others were negative (Figure 5b). Notably, Tajima's D values for four BnaOBMP genes, including BnaOBMP.S13 BnaOBMP.T25 (−2.247), were less than −2, suggesting that these genes were under strong selection (positive or negative). The selection pressures on the twelve BnaOBMPs in rapeseed were also explored based on the Ka/Ks ratios. A Ka/Ks ratio greater than 1, equal to 1 and less than 1 indicated positive selection, neutral evolution and purifying selection at a low evolutionary rate, respectively. Among these BnaOBMP genes, except for BnaOBMP.S13, the Ka/Ks ratio values of the others were all less than 1 ( Figure S6).

Quantitative Real-Time PCR (qRT-PCR) Analysis of the BnaOBMP Genes
To further validate the functional roles of BnaOBMPs under abiotic stresses and during seed development, six of the BnaOBMP genes from different subfamilies were selected for examining expression levels using qRT-PCR in the cultivated species of rapeseed (ZS11). These genes included two caleosins (i.e., BnaOBMP.C19 and BnaOBMP.C7), three oleosins (i.e., BnaOBMP.SH8, BnaOBMP.SL4 and BnaOBMP.U9) and one steroleosin (i.e., BnaOBMP.S13). After a 3 h heat treatment, the expression level of BnaOBMP.C7 was significantly downregulated by approximately five-fold compared with the control, while BnaOBMP.SH8 was significantly upregulated in the leaves of rapeseed (Figure 7). After a 3 day drought treatment, the expression levels of BnaOBMP.C7, BnaOBMP.U9 and BnaOBMP.SH8 were reduced (Figure 7). In contrast, BnaOBMP.C19, BnaOBMP.SL4 and BnaOBMP.S13 showed no detectable expression in rapeseed leaves. Compared to seeds of 15 DAF (early seed development), the expression levels of BnaOBMP.C19, BnaOBMP.SH8, BnaOBMP.S13 and BnaOBMP.U9 increased in seeds of 25 DAF (seed filling stage) and reached a peak in seeds of 35 DAF (seed filling stage) and then decreased in seeds of 50 DAF (seed maturation stage). Both BnaOBMP.C7 and BnaOBMP.SL4 exhibited significantly induced expression during the development of seed (Figure 7). qRT-PCR analysis validated the gene expression patterns identified by RNA-seq. Figure 7. qRT-PCR verification that the expression of some representative BnaOBMP genes responded to heat and drought stresses (leaves) as well as during seed development. Statistically significant differences (Student's t-test) are indicated as followed: * p < 0.05, ** p < 0.01, *** p < 0.001.

Co-Expression and Gene Regulatory Network of BnaOBMPs during Seed Development
In order to comprehensively analyze the gene interactions and regulatory relationships of OBMP genes during seed development, the previous seed RNA-seq data were employed to build a co-expression and gene regulatory network. We calculated the Pearson correlation coefficients (PCCs) of the expression levels between BnaOBMP gene pairs ( Figure S8 To explore transcription factors involved in the regulation of OBMP genes during seed development, we further constructed a seed-specific gene regulatory network (GRN) using GENIE3 algorithm in B. napus. By filtering for the genes expressed > 1 FPKM in at least one of the seed samples, we captured 61,907 genes, including 4643 transcription factors (TFs) and 64 BnaOBMP genes, for GRN construction. The top 5 TFs of each BnaOBMP connected were extracted from the GRN and used to construct a subnetwork of BnaOBMP-GRN (Figure 8b). The classification of these TFs in BnaOBMP-GRN revealed that the expression of BnaOBMP genes were potentially regulated by various types of TFs such as ZIP, C2C2-GATA, HB-HD-ZIP, HB-WOX, C2H2, NF-YB, C3H, bHLH, AP2/ERF-ERF, NAC, GARP-G2-like, MYB and C2C2-dof ( Figure 8b). Notably, a total of 159 hub TFs were identified in BnaOBMP-GRN (Table S3)

Discussion
OBs share common features in all kingdoms of life [42], consisting of a densely packed hydrophobic core of neutral lipids enclosed by a phospholipid monolayer decorated by three main classes of oil-body-membrane proteins (OBMP) including oleosins, caleosins and steroleosins [13,42]. Oleosins are the most abundant protein constituents and are sufficient to cover the whole surface of a seed OB [43]. Oleosins play important roles in the formation and stabilization of OBs during seed and pollen grain development as well as in OB turnover [26]. Caleosin is a common OB surface protein found in a wide range of plant species [44]. Caleosin has a Ca 2+ -binding motif, which has the ability to bind calcium. Recent studies have suggested that caleosins also possess peroxygenase activities that convert hydroperoxides of α-linolenic acid to various oxylipins as phytoalexins [27,45]. Steroleosins belong to hydroxysteroid dehydrogenases (HSDs), consisting of a sterolbinding site and an NADPH-binding site involved in some biological functions related to membrane remodeling and lipid signaling [32,46,47]. Thus, OBMP is not only the key structural molecule for the formation and stabilization of OBs but may also exert a myriad of cellular functions related to carbon, energy and lipid metabolisms; stress responses; hormone signaling pathways; be involved in various aspects of plant growth and development.
As the second-largest source of vegetable oil, rapeseed is an important worldwide oilseed crop [2]. Their seeds contain lipids as major storage reserves, which is up to 50% of the dry weight, and the main component of lipids is triacylglycerol stored in oil bodies [48]. To investigate the function of oil-body-membrane proteins from important oil crops, we performed identification and characterization of the OBMP gene family in the polyploid crop B. napus in the present study. A total of 88 OBMP genes were found in all of the nineteen chromosomes of rapeseed, which were classed into 53 oleosins (27 T, 8 SL, 9 SH and 9 U), 20 caleosins and 15 steroleosins based on their functional domains and phylogenetic relationships (Table 1). Compared to the number of OBMP genes identified in the species from green algae to higher plants, B. napus contained the most abundant OBMP genes, suggesting the massive expansion of this gene family in rapeseed ( Figure S1). The A and C subgenomes of B. napus contained 43 and 45 BnaOBMP genes, respectively, which is comparable to the 44 and 47 genes identified in B. rapa and B. oleracea, respectively (Tables 1 and S1). Moreover, the BnaOBMP genes from the A and C subgenomes showed conserved synteny and gene order (Figure 3). These results indicate that B. napus retained the vast majority of OBMP genes from its two ancestors during the allopolyploidization.
Gene and genome duplications gave rise to the number of genes, resulting in functional redundancy and differentiation of genes, enabling genome-wide adaptation to various environments during evolution [37]. The previous study revealed that the crucifer (Brassicaceae) lineage experienced two whole-genome duplications (WGDs) and one whole-genome triplication event (WGT), shared by most dicots [49]. Moreover, the Brassica species experienced an extra WGT event approximately 15.9 million years ago compared with A. thaliana [50]. B. napus is a relatively new species of the Brassica genus, with a short history of post-Neolithic speciation (~7500 years) and domestication (~700 years). Consistent with these WGD and WGT events, the OBMP genes experienced gene duplication events leading to an expanded OBMP gene family in rapeseed. The WGD or segmentally duplicated genes accounted for the majority of BnaOBMP gene family. Particularly, all of the caleosins and steroleosins resulted from WGD or segmental duplication in rapeseed. The tandemly duplicated genes were also detected in the BnaOBMP gene family, and all fifteen tandemly duplicated genes belonged to T oleosins, which might be the result of the aggregated distribution of T oleosins in the chromosomes of the rapeseed genome. Furthermore, the BnaOBMP genes showed different levels of polymorphism. Most of the BnaOBMP genes had Ka/Ks ratios less than 1 and Tajima's D values less than 0, suggesting that the BnaOBMP gene family experienced strong purifying (stabilizing) selection rather than positive selection during the evolution. In addition, the Ka/Ks ratios were substantially highest among the T oleosins than other BnaOBMP genes, implying that T oleosins evolved faster than the other BnaOBMP genes in rapeseed. Altogether, these results revealed that the BnaOBMP gene family expanded primarily by gene duplications with WGD/segmental duplication being the major driving force in B. napus.
Oil bodies, as essential lipid storage organelles in the seeds of plants, play important roles in seed germination and the postgerminative growth of seedlings, as well as many other cellular processes such as stress responses, lipid metabolism, organ development, and hormone signaling. These biological functions of seed OBs depend on OBMP proteins, which are embedded in the OB phospholipid monolayer. Our results revealed that various cis-acting regulatory elements exist in the promoters of BnaOBMP genes. The cis-acting regulatory elements exert various functions associated with plant growth and development, phytohormone responsive, and abiotic and biotic stress responsive such as ABRE (abscisic acid-responsive element), AuxRE (auxin-responsive element), ERE (ethyleneresponsive element), ARE (anoxic-responsive element), DIRE (drought-responsive element) and LTRE (low-temperature-responsive element). This suggests that BnaOBMP genes could be induced by different phytohormone and stress signals so as to adjust OBs to different environmental conditions. Unraveling the expression pattern of different OB proteins throughout seed development is crucial for improving our understanding of OB formation. A previous study revealed that accumulation of oleosins S1-S5 and caleosin CLO1 began at approximately 12 days after pollination (DAP), while steroleosin SLO1 accumulated later at approximately 25 DAP, and then they all increased rapidly and reached a peak at 55-60 DAP in A. thaliana, as analyzed by immunoblot [51]. We found that the expression of most BnaOBMP genes were upregulated along the development of the seed and showed the highest expression levels at the late stages, which is consistent with the previous study. Although OBs occur minimally in nonstorage vegetative organs, we also observed that some BnaOBMP genes exhibited high expression in flower bud, pistils, leaves, stems and roots, indicating that OBs are present in these tissues. For example, T oleosins, except for BnaOBMP.T1 and BnaOBMP.T8, were preferentially highly expressed in flower bud. Some caleosins showed higher expression levels in flower bud, stamen, pistil, leaf and root. For steroleosins, BnaOBMP.S1 showed the highest expression level in pistil, while BnaOBMP.S15 was highly expressed in root. These results suggested that OBs could be not only present in seed as storage warehouses but also exist in nonstorage vegetative organs as detoxification refuges.
Genes involved in the same process usually have similar expression patterns, and they typically grouped into the same module in co-expression analysis. Here, based on the co-expression networks of BnaOBMPs, the BnaOBMP family genes showed a similar tendency and gathered closely together in one cluster during seed development in B. napus. Strong positive correlations were observed between the members of the BnaOBMP gene family. Meanwhile, only a few significant negative correlations appeared between BnaOBMP.C15 and the others. The results suggested that the BnaOBMP genes might share functional redundancy in B. napus. To discern gene transcriptional regulatory mechanisms of BnaOBMP genes, we constructed a regulatory network incorporating TF information by GENIE3. Previous studies demonstrated that a GENIE3 network could provide biologically relevant transcription factor-target relationships in wheat [52,53]. Our BnaOBMP-GRN network revealed that the BnaOBMP genes could be regulated by various transcript factors. After prioritization of the candidate regulatory genes, the top hub transcription factor was the bZIP gene of BnaA04G0262900ZS encoded ABSCISIC ACID-INSENSITIVE 5 (ABI5)-like protein 3 (EEL). The Arabidopsis EEL (known as AtbZIP12) is transcription factor homologous to ABI5, which is a key player in light-, abscisic acid-, and gibberellin-signaling pathways to precisely control seed maturation and germination [54][55][56]. In addition, some other hub TFs, such as GATA3, HAT2, SMZ, DOF5.6 and APL, might also play important roles in the regulation of BnaOBMPs during OB formation in rapeseed, which needs further investigation.
In conclusion, our results reveal that B. napus had an expansion of the OBMP gene family due to the fact of WGD and tandem duplications. These BnaOBMP genes contain extensive sequence polymorphisms, and some members may have experienced strong selection. Various cis-acting regulatory elements involved in plant growth, phytohormone and abiotic and biotic stress responses were found in their promoter regions. In addition, both transcriptomic and qRT-PCR analyses corroborated that BnaOBMPs exhibited spatiotemporal expression patterns and are preferentially expressed in seeds. The genetic variations (i.e., SNPs or InDels) of BnaOBMP genes can be used as molecular markers to select rapeseed cultivars with high seed oil content (SOC). Moreover, further manipulating the expression patterns of some candidate BnaOBMPs during seed development using genetic engineering techniques, such as transgenic technology and CRISPR/Cas9 tools, would contribute to the increase in the SOC in rapeseed.

Identification and Property Analysis of BnaOBMP Genes
The OBMP genes in A. thaliana were retrieved from the TAIR (http://www.arabidopsis.

Gene Structure and Chromosomal Localization Analysis
Based on the genome annotation of B. napus var. ZS11, available from the BnPIR database (http://cbi.hzau.edu.cn/bnapus/ (accessed on 8 January 2022)), a graphical representation of the exon-intron structure of each BnaOBMP gene was drawn using Gene Structure Display Server 2.0 (http://gsds.cbi.pku.edu.cn/ (accessed on 12 March 2022)). The schematic map of the BnaOBMP genes on chromosomes was plotted using R software with RIdeogram package (version 0.2.2) according to their physical chromosomal locations on B. napus.

Cis-Acting Regulatory Elements and Motif Analysis
The 2 Kb sequence upstream from the transcription start site of each BnaOBMP gene was defined as the promoter region and was extracted from the B. napus genome sequence. The cis-acting regulatory elements within the promoters were analyzed by PlantCARE (https: //bioinformatics.psb.ugent.be/webtools/plantcare/html/ (accessed on 23 March 2022)).

Multiple Sequence Alignments and Conserved Motif Analysis
MUSCLE (version 3.8) was used to align the OBMP protein sequences of B. napus, B. rapa, B. oleracea and A. thaliana with default parameters. The MEME suite (https://memesuite.org/meme/ (accessed on 14 April 2022)) was employed to identify the motifs with conserved amino acids in the protein sequences of BnaOBMPs with default parameters.

Phylogenetic Analysis
Based on the multiple sequence alignments of OBMP proteins from B. napus, B. rapa, B. oleracea and A. thaliana, IQ-TREE (version 2.0.3) was employed to reconstruct a maximum likelihood (ML) gene tree with 1000 replicates. The VT + F + R4 model was the best-fit evolutionary model selected by ModelFinder implemented in IQ-TREE. The obtained ML gene trees were visualized using iTOL (https://itol.embl.de/ (accessed on 24 March 2022)).

Synteny and Duplicate Gene Analysis
The syntenic blocks and gene duplications were identified within the B. napus, B. rapa, B. oleracea and A. thaliana genomes using MCScanX (version 3.8.31) with the default parameters. The circos graph with the genomic collinearity of BnaOBMP genes was plotted using the Rcircos package in R software (version 3.6.1).

Sequence Variations and Polymorphism Analysis
The publicly available genomic resequencing dataset of different B. napus accessions from around the world (database accession: SRP155312) was collected from the National Center of Biotechnology Information (NCBI) database [41]. The alignment and variant calling were performed using Sentieon DNASeq Variant Calling Workflow [57]. Variant annotation was achieved by ANNOVAR (version 20220320) based on the annotation of the B. napus var. ZS11 genome. The sequence diversity metrics, including pairwise nucleotide variation as a measure of variability (π) and selection statistics (Tajima's D) of BnaOBMP genes, were calculated by vcftools (version 0.1.13) based on SNP distribution.

Expression Analysis of BnaOBMP Genes
The publicly available RNA-seq datasets (database accession: PRJNA394926 and PRJNA311067) of eight tissues (i.e., root, stem, leaf, flower bud, sepal, stamen, pistil and seed) collected from different developmental stages [58,59] and an RNA-Seq dataset (database accession: CRA003544) of 20 DAF and 40 DAF seeds of 280 B. napus germplasm accessions [60] were downloaded from the NCBI and NGDC databases, and they were used for gene expression profiling. Mapping of these RNA-Seq reads against the B. napus var. ZS11 reference genome using HISAT2 (version 2.1.0) with the default settings. The reads count per gene was calculated with HTSeq (version 0.9.172) and was further used to calculate the FPKM values for the quantification of gene expression. The heatmaps were visualized using the pheatmap package (version 1.0.12) in R software (version 3.6.1). The principal component analysis was performed using the function prcomp() in the R software.

Plant Materials and Treatments
The seeds of the rapeseed cultivar ZS11 were germinated on a filter paper saturated with distilled water in darkness at 22 • C for 3 days. The seedling plants were transplanted to soil culture pots in a greenhouse to grow for six weeks under well-controlled conditions as follows: a temperature of 25 • C, light intensity of 150 µmol m −2 s −1 provided by a highpressure sodium lamp, and a humidity of 50-60%. Then, the rapeseed plants before the bolting stage were chosen for analysis. The top third of fully expanded leaves from rapeseed plants were sampled as a control before being stressed for the following experiments. The drought and heat stress treatments were conducted in plant growth chambers with wellcontrolled temperature and humidity. For the drought treatment, water was withdrawn for 7 days, and then the plants were rewatered to recover from the stress. The growth chamber was programmed as follows: 40% humidity in 16 h light at a temperature of 25 • C; 45% humidity in 8 h dark at a lower temperature of 22 • C. Leaf samples were collected 3 days after drought treatment. For heat treatment, the growth chamber was set at 60% humidity in 16 h light at a high temperature of 40 • C, followed with 55% humidity in 8 h dark at a temperature of 35 • C. Heat-treated leaf samples were collected at a time point of 3 h during the stress treatment. Rapeseed seeds at 15, 25, 35 and 50 days after flowering of the ZS11 plants cultivated in the experimental field were sampled for analysis. All samples were immediately frozen using liquid nitrogen after being detached, and they were stored at −80 • C for further assay.

RNA Isolation and qRT-PCR Analysis
Total RNA was extracted from each sample using an RNA extraction kit (Takara, Dalian, China) following the manufacturer's procedure. Two micrograms of total RNA were used to synthesize the first-strand cDNA using the Prime Script RT reagent Kit (Takara, Dalian, China) according to the manufacturer's protocol. Quantitative real-time PCR was performed using 2 µL of cDNA in a 20 µL reaction volume with SYBR Premix Ex Taq (Takara) on a 7500-Fast real-time PCR System (Applied Biosystems). Gene-specific primers were designed and are listed in the Supplementary Materials (Table S4). The thermal cycler was set as follows: an initial incubation at 50 • C for 2 min and 95 • C for 5 min, followed by 40 cycles at 95 • C for 30 s, 55 • C for 30 s and 72 • C for 30 s. All qRT-PCR reactions were assayed in triplicates. The relative quantification of transcription level was determined by the methods described previously [61].

Pearson Correlation and Gene Regulatory Network Analysis
Based on the previous seed RNA-seq (database accession: CRA003544) analysis results, the PCC values and significant of correlations among BnaOBMP genes were calculated using the function rcorr() of R software (version 3.6.1). The genes with FPKM >1 in at least one of the 560 seed samples of B. napus were selected and used to construct the gene regulatory network (GRN) using GENIE3 (version 1.19.0). The TFs of B. napus identified by iTAK (version 1.7) were used as candidate regulators for the GRN construction. Gephi (version 0.9.2) was used to calculate network metrics and create visualization graphs.
Supplementary Materials: The following supporting information can be downloaded at: https: //www.mdpi.com/article/10.3390/plants11172241/s1, Figure S1: Number of OBMP genes identified from each of the plant genomes; Figure S2: Pileup of the sequences of all caleosins in B. napus aligned with the MUSCLE method in MEGA software; Figure S3: Pileup of the sequences of all U, SL, SH and T oleosins in B. napus aligned with the MUSCLE method in MEGA software; Figure S4: Pileup of the sequences of all steroleosins in B. napus aligned with the MUSCLE method in MEGA software; Figure S5: Selective pressures for OBMP genes in B. napus as revealed in nonsynonymous and synonymous substitution (Ka/Ks) ratios; Figure S6: Ka, Ks and Ka/Ks ratios of 12 OBMP genes under selection in B. napus; Figure S7: PCA plot of all the seed samples of 280 B. napus accessions; Figure S8: Correlation matrix of gene expression levels among BnaOBMP genes; Table S1: Summary information of the OBMP gene family in A. thaliana, B. rapa and B. oleracea; Table S2: Statistic of the variations of BnaOBMP genes; Table S3: Hub genes of BnaOBMP-GRN network in B. napus; Table S4. List of primers for qRT-PCR.

Conflicts of Interest:
The authors declare no conflict of interest.