*Article* **Exogenous Phytohormones and Fertilizers Enhance** *Jatropha curcas* **L. Growth through the Regulation of Physiological, Morphological, and Biochemical Parameters**

**Rahmatullah Jan 1,2,†, Murtaza Khan 3,†, Muhammad Adnan 4, Sajjad Asaf 5, Saleem Asif 1, Kyung-Min Kim 1,2,\* and Waheed Murad 6,\***


**Abstract:** *Jatropha curcas* L. is a perennial plant, that emerged as a biodiesel crop attracting the great interest of researchers. However, it is considered a semi-wild plant and needed to apply cropimproving practices to enhance its full yield potential. This study was conducted to improve the growth and development of the *J. curcas* plant by exogenous application of Gibberellic acid (GA), indole acetic acid (IAA), and fertilizer (nitrogen, phosphorus, potassium (NPK)). The experiment was conducted in pots in triplicate and 100 ppm and 250 ppm of GA and IAA were applied separately while NPK was applied in two levels (30 and 60 g/pot). The results revealed a significant difference in growth parameters with the application of hormones and fertilizer. The highest shoot length (47%), root length (63%), root fresh weight (72%), and root dry weight (172%) were shown by plants treated with GA 250 ppm. While plants treated with NPK 60 g showed the highest increases in shoot fresh weight and shoot dry weight compared to control plants. The highest increase in leaves number (274%) and branches number (266%) were shown by the plants treated with GA 100 ppm and GA 250 ppm, respectively, while GA 250 ppm and IAA 250 ppm highly increased stem diameter (123%) and stem diameter was also shown by GA 250 ppm-treated plants. NPK 60 g highly increased proximate composition (protein content, carbohydrate, fat, moisture content, and ash content) compare with hormones and control plants. Our results concluded the optimized concentration of IAA, GA, and NPK significantly increases *J. curcas* growth vigor.

**Keywords:** biodiesel; gibberellin; growth parameters; indole acetic acid; *Jatropha curcas*; proximate composition

#### **1. Introduction**

*Jatropha curcas* L. is a perennial shrub, belonging to the family Euphorbiaceae. It is about 5 m tall and has smooth grey bark, leaves are large and usually pale green. Fruits are produced in winter or throughout the year depending on temperature and soil humidity. The indications show that it originated from South and Central America and some other parts of the tropical and subtropical regions of Africa and Asia [1]. Specifically, it is grown in Benin, Brazil, China, Egypt, Ethiopia, Ghana, Guinea, India, Madagascar, Mali, Mexico, Mozambique, Namibia, Senegal, South Africa, Sudan, Tanzania, Uganda, Zambia, and Zimbabwe [2]. Its name indicates that it is used as a medicinal plant in the Portuguese in the 16th century, as its name is derived from the Greek word "iatros" which means doctor,

**Citation:** Jan, R.; Khan, M.; Adnan, M.; Asaf, S.; Asif, S.; Kim, K.-M.; Murad, W. Exogenous Phytohormones and Fertilizers Enhance *Jatropha curcas* L. Growth through the Regulation of Physiological, Morphological, and Biochemical Parameters. *Plants* **2022**, *11*, 3584. https://doi.org/10.3390/ plants11243584

Academic Editors: Andrzej Bajguz and Othmane Merah

Received: 27 October 2022 Accepted: 13 December 2022 Published: 19 December 2022

**Publisher's Note:** MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

**Copyright:** © 2022 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (https:// creativecommons.org/licenses/by/ 4.0/).

and "trophe" which means food [3]. Literature shows that the Portuguese established the plantation of *J. curcas* for the first time to make soap, lamp oil, and other medicine [1].

However, *J. curcas* is exotic to Pakistan and introduced in 2007 to Karachi. It is an oilseed plant that does not compete with food crops and its seed contains 34–60% oil contents [3,4]. It has the ability to adapt to a high range of agro-climatic conditions. It grows in gravelly, sandy, degraded, acidic, and poor stony soils. It is easy to establish, grows quickly, drought tolerant, grow well in low rainfall condition, and can be used to reclaim eroded areas [5]. Therefore, it is suggested for cultivation on poor degraded soil. However, it shows stunt growth under heavy metals stress [6]. It can survive in long dry periods and it is considered to be well adapted to arid and semi-arid conditions. In the initial stage, it requires a large amount of water but after maturity, it may survive without water for sixty days. In the rainfed area, it required 250 up to 3000 mm/year of rainfall for best growth [7]. Gadallah and Sayed (2001) reported that the exogenous application of hormones increases *J. curcas* resistance against environmental constraints and promotes plant growth and development, however, the effectiveness of growth regulators depends on their concentration and the method of application [8].

To promote plant growth and development in non-suitable conditions, phytohormones or macronutrients are important candidates. Trail studies indicated that phytohormones and fertilizers increase plant growth and development under stress conditions. Gibberellin (GA) increases plant growth, leaves number, bud formation, cell division and elongation, and flowering [9,10]. Indol acetic acid (IAA) is also involved in several developmental processes such as the development of vascular tissue, cell elongation, and apical dominance [11]. The growth regulator enhances the growth of *J. curcas* by regulating its morphophysiological and biochemical processes [8]. Besides hormones, many other parameters also affect the growth of *J. curcas*. Fertilizer is one of the essential parameters which promote plant development in a harsh environment. Nitrogen, phosphorous, and potassium (NPK) are the main nutrients needed for plant growth and development. Without providing adequate NPK, plants cannot reach maximum growth [12]. The NPK supplementation increases the biomass and oil contents in *J. curcas,* which shows the importance of fertilizers for the growth and development of *J. curcas*. Being a biofuel plant, *J. curcas* attracted the researcher's attention. As it grows well in drought and high temperatures; therefore, it is suitable for arid regions where it is not competing with the crops for water and land.

Due to the potential demand and better opportunities, cultivation of *J. curcas* appears viable. It adapts well to marginal lands and the large-scale cultivation on wasteland with low water and rainfall could generate employment and increase the income of the locality. However, the ability of *J. curcas* to grow in marginal land and dry soil has not been properly explored [13]. Several studies have been conducted to evaluate the *J. curcas* performance under low water availability in marginal land [13]. In order to understand the plant growth performance in the marginal land of Kohat, Khyber Pakhtun Khwa, Pakistan, region. Kohat is an arid region containing a large quantity of marginal land, which is suitable for *J. curcas* cultivation. Based on *J. curcas* plant adaptation to marginal land, we conducted a trail base studyin the Kohat region (arid region) district in Khyber Pakhtunkhwa, Pakistan. However, the result of plant growth was not satisfactory under normal growth conditions (data not published). To further investigate the root cause of our previous project failure, we performed the current study to investigate the possible way to promote the growth of *J. curcas* in the same region. In the current study, we focus on exogenous hormonal application and fertilizer supplementation. We hypothesized that GA, IAA and NPK promote *J. curcas* growth and development by regulation of morphological, physiochemical and biochemical characteristics of *J. curcas*. Therefore, the current study aimed to investigate that the exogenous hormones and fertilizer induces morphological, physiological, and biochemical parameters of the *J. curcas* plant and enable it to grow and develop in arid regions. Based on our study, we can further suggest the cultivation of the *J. curcas* plant in the barren land after phytohormones and fertilizers supplementation. This study is of great importance and thereby suggests that *J. curcas* growth and development increases with

the exogenous application of phytohormones and fertilizers. *J. curcas* is an economically important plant and unlike crop plants, it does not need significant attention and does not compete with other crops for land and water. The current study added a new insight to the field that application of phytohormones and fertilizers can promote *J. curcas* growth and development in the marginal land under limited water availability.

#### **2. Results**

#### *2.1. Phytohormones and Fertilizer Enhances J. curcas Growth Vigor*

In the present study, it is investigated that IAA, GA, and NPK increased *J. curcas* shoot length (Figure 1). IAA 100 ppm and IAA 250 ppm increased significantly the shoot length up to 25% and 29%, respectively, as compared to control plants after four months of growth. The GA 100 ppm, GA 250 ppm, and NPK 30 g significantly increased the shoot length up to 17, 47, and 41%, respectively, after four months of growth, compared with control plants (Figure 1A). The higher concentration of IAA, GA, and NPK also increased shoot fresh and dry weight (Figure 1B,C). The results indicated that IAA 250 ppm significantly increased fresh and dry weight up to 141 and 183%, respectively. GA 250 ppm increased shoot fresh and dry weight up to 141 and 179%, respectively, while NPK 60 g increased by 185 and 267%, respectively, compared to control plants. NPK-applied plants showed higher fresh and dry weight than the GA and IAA. On the other hand, hormones and fertilizer also induced root biomass. The IAA 250 ppm, GA 100 ppm, GA 250 ppm, and NPK 60 g significantly increased root length up to 31, 50, 63, and 40%, respectively (Figure 1D). Root fresh weight was significantly increased by 21% by IAA 250 ppm, 72% by GA 250 ppm, and 30% by NPK 60 g compared with control plants (Figure 1E). The root dry weight followed the same pattern of the root fresh weight. The IAA 250 ppm increased the root dry weight to 81%, GA 250 ppm increased to 172%, and NPK 30 g and 60 g increased to 63 and 118%, respectively (Figure 1F). We investigated that, among the hormones and fertilizer, the highest shoot length was shown by the GA 250 ppm-treated plants, while the highest shoot fresh and dry weight was shown by NPK 60 g-treated plants. However, the GA 250 ppm-treated plant showed the highest increase in root length, root fresh, and dry weight which shows that GA can induce *J. curcas* root growth and development.

**Figure 1.** IAA, GA, and NPK regulate *J. curcas* plant soot, root length, shoot, root fresh and dry weight. (**A**) represents shoot length, (**B**) represent shoot fresh weight, (**C**) represent shoot dry weight. (**D**–**F**) represent foot length, root fresh weight, and root dry weight, respectively. Graphs show mean ± standard deviation, and asterisks show significant differences (\* *p* ≤ 0.05, \*\* *p* ≤ 0.01, and \*\*\* *p* ≤ 0.001) according to two-way ANOVA and Bonferroni post hoc tests.

#### *2.2. Phytohormones and Fertilizer Treatment Promotes Leaf and Branch Number and Stem Diameter*

In the current study, we found that IAA, GA, and NPK efficiently stimulated *J. curcas* number of leaves and branches, and stem diameter in four months (Figure 2). Both low and high concentrations of hormones and NPK application significantly increased the number of leaves, branches, and stem diameter as compared to control plants. IAA 100 ppm and 250 ppm increased leaves number by 262 and 211%, respectively, while branches number increased by 200 and 160%, respectively, and stem diameter increased by 100 and 123%, respectively. GA 100 ppm and 250 ppm increased leaves number up to 274 and 211%, branches number 229 and 260%, and stem diameter 109 and 123%, respectively, as compared to control plants after four months. The 30 g and 60 g of NPK application increased leaves number by 199 and 211%, branches number by 129 and 229%, and stem diameter by 76 and 80%, respectively, after four months of treatment. Among the IAA, GA, and NPK, IAA 100 ppm and GA 250 ppm showed the highest number of leaves increase while GA 250 and NPK 60 g showed the highest induction of branches. Whereas the highest increase in stem diameter was found in IAA 250 ppm and GA 250 ppm treated plants followed by GA 100 ppm. However, the lowest stem diameter and branches number were found in the NPK 30 g treated plants, whereas the lowest leaves number were found in the IAA 100 ppm, GA 250 ppm and NPK 60 g treated plants (equally reduced 211%). These results show that the optimized concentration of IAA, GA, and NPK is needed to increase *J. curcas* leaves, branches, and stem diameter.

**Figure 2.** *Cont*.

#### *2.3. Phytohormones and Fertilizers Promote J. curcas Plant Growth by Regulating Proximate Compositions*

In the current study, we determined the protein, carbohydrate, and fat contents after four months of treatment of IAA, GA, and NPK in different concentrations. The higher concentration of hormones and NPK such as IAA 250 ppm, GA 250 ppm, and NPK 60 g showed a significant increase of protein, carbohydrates, and fat contents compared with control plants (Figure 3). The NPK 60 g treated plants showed the highest increase (34%) in protein contents followed by IAA 250 ppm (27%), and GA 250 ppm (25%) compared with control plants (Figure 3A). IAA 100 ppm also increased significantly the protein content up to 20%, which indicates that increasing the IAA concentration can increase protein contents. Carbohydrates were significantly increased by 61% by NPK 60 g followed by GA 250 ppm 42% and IAA 250 ppm 30%, compared with control plants (Figure 3B). Similarly, the highest increase (46%) in fat contents was found in NPK 60 g treated plants followed by IAA 250 ppm (39%) and GA 250 ppm (32%), compared with control plants (Figure 3C). The highest moisture and ash contents were found in NPK 60 g treated plants followed by IAA 250 ppm and GA 250 ppm treated plants compared with control plants (Figure 3E,F). These results suggested that phytohormones and NPK treatment enhance the growth of *J. curcas* by induction of proximate compositions. Among the phytohormones and NPK, it was found that NPK application highly increased the overall proximate composition compared with IAA and GA. These results show that fertilizers can better increase the *J. curcas* succession in temperate and barren land compared with hormones.

**Figure 3.** IAA, GA, and NPK induces biosynthesis of *J. curcas* plant proximate composition contents. (**A**) represents protein contents, (**B**) represents carbohydrate contents, (**C**) represents fat contents, (**D**) represents moisture contents, and (**E**) represents ash contents. Graphs show mean ± standard deviation, and asterisks show significant differences (\* *p* ≤ 0.05, \*\* *p* ≤ 0.01, and \*\*\* *p* ≤ 0.001) according to two-way ANOVA and Bonferroni post hoc tests.

#### **3. Discussion**

In the present study, we reported that IAA, GA, nitrogen, phosphorus, and potassium (NPK) act as growth enhancers of the *J. curcas* plant. We provided physiological, morphological, and biochemical evidence that *J. curcas* growth and development were enhanced by the application of hormones and NPK. Generally, the *J. curcas* plant grows in almost all kinds of soil and drought conditions; however, hormones and fertilizers increase their growth and development. The applied hormones (GA, IAA) and NPK greatly influenced the plant's morphological, physiological, and biochemical traits.

IAA is an essential auxin with significant in vivo roles such as stem growth, root growth, stem cambium cell activation, and lateral bud formation [14,15]. In most plants, IAA is the supreme active form of auxin. Recent studies revealed that IAA application has positive effects on growth parameters such as plant root and shoot length, fresh and dry weight, number of leaves, chlorophyll contents, carbohydrates, amino acids, and phenolic contents [16,17]. Our results also showed that IAA increased plant growth parameters (shoot length, root length, fresh and dry weight, branches and leaves number, and stem diameter) in *J. curcas* plants compared with control plants. In a recent study, it is predicted that phytohormones including IAA alter the sugar metabolism which is responsible for the modulation of biological processes that are involved in plant growth promotion [18,19]. IAA

interacts with sucrose and alters plant morphogenesis which regulates leaf morphology [20]. In addition, IAA increases endogenous GA accumulation, which is a prominent growth regulator [18]. In our study, we found that IAA increased leaves and branches number. The possible reason for the increased number of leaves and branches number might be the suppression of the ABA hormone. ABA generally inhibits growth or keeps the apical tissue dormant while IAA reduces ABA and increases GA which breaks dormancy and results in plant growth and increasing branches number [18]. Besides other growth parameters, we also found that the IAA-treated plant's stem diameter was increased significantly than the control plants, which predicts that IAA is efficiently involved in vascular cell division and differentiation. It is reported that IAA increases *Persea americana* cell differentiation which increases vascular vessel density [21]. In *Glycin max*, IAA stimulated plant height, leaves number and area, number of branches, and seed per plant [22]. However, different plants show different reactions to IAA in various concentrations. In the *J. curcas* plant, between the IAA 100 ppm and 250 ppm, the IAA 250 ppm showed a higher increase in morphological, physiological, and biochemical parameters. However, IAA 100 ppm increased the *G. max* height more than that of 200 pmm [22]. Foliar application of IAA on cowpea plants also increased plant height, fresh and dry weight, number of branches, number of leaves, and yield components [23]. As far as IAA is concerned to plant height, some researchers provided evidence that IAA promoted GA synthesis. For instance, it is reported that exogenous IAA increased GA1 and GA3 biosynthesis through the activation of GA1 and GA3 synthesis enzymes [24,25]. IAA not only induced morphological and physiological parameters of the *J. curcas* plant, but also induced their proximate composition (protein, carbohydrate, fat, moisture, and ash contents) compared to control plants. To the best of our knowledge, there is no data published on the effect of IAA and GA on the proximate composition of *J. curcas*; however, researchers investigated the effects of IAA and GA on the proximate composition of various plants. A recent study investigated that, exogenous application of IAA enhanced the photosynthesis rate in *Gossypium hirsutum* which resulted in increased biomass including fresh and dry weight [18]. In the *Balanites aegyptiaca* plant, different concentrations of IAA increased total protein and carbohydrates when compared to non-treated plants [26]. IAA also increased fresh and dry weight, relative water contents (moisture contents), and chlorophyll contents in white clover plants [24]. These results are in line with our findings, which show that IAA is involved in the biosynthesis of protein and carbohydrates and increases moisture contents in *J. curcas*.

GA is a key hormone that induces many parameters of plant growth and development such as seed germination, shoot elongation, leaf expansion, and flower and fruit development [27]. In our study, we found that GA exogenous application increased shoot and root length, fresh and dry weight, and proximate composition percentage in *J. curcas* plants compared with control plants. Plant biomass accumulation is associated with endogenous GA accumulation. GA enhances carbonic anhydrase activity which promotes CO2 fixation in photosynthesis as this enzyme takes part in the hydration of CO2 and is strictly related to chloroplast [28]. This phenomenon provides enough CO2 to the site of fixation and increases the photosynthesis rate which as a result increases biomass accumulation [29,30]. In comparison to IAA and NPK, GA greatly influenced the shoot length, root length, and root fresh and dry weight (Figure 1). Researchers evaluated the basic mechanism of GA that regulates plant growth and development. Previous studies investigated that GA induces the transcription of genes that are involved in cell division and cell elongation occurring during plant growth [31]. GA also stimulates the induction of hydrolytic enzymes which are involved in the conversion of starch to sugar and the controlling of starch and sugar accumulation by GA can significantly influence plant growth and development [32]. Mostly the GA induces plant growth by alteration of certain genes expression. GA also plays a key role in several metabolic pathways that effecting plant growth such as chlorophyll biosynthesis, nitrogen metabolism and redistribution, and translocation of assimilates [28]. GA influences the differentiation of phloem fiber and enhances the length of bast fibers by inducing internode length. Reports show that a high level of GA results in primary phloem

fiber elongation in the *Coleus blumei* plant and the length of differentiating internode is associated with the length of the primary phloem [33]. The increase in phloem fiber is associated with an increase in plant height and an increase in plant intermodal length [34]. Previous studies show that GA efficiently stimulates the lateral branches outgrowth in the *J. curcas* plant and it was investigated that increasing concentration increases shot branching [35]. In *Arabidopsis thaliana*, GA insensitive mutant shows a reduction in apical dominance and an increased number of axillary shoots [36]. GA sometimes shows an inhibitory role in lateral branching depending on plant species. In *Pisum sativum* plants, GA plays an inhibitory role in lateral bud formation [37]. The overexpression of the GA biosynthesis gene increased tillers or branches in *Paspalum distichum* and *Populus tremula*, suggesting that GA is significantly involved in the branching of these species [38–40]. However, GA-induced bud formation in citrus, sweet cherry, and rose [41–43]. Our results were consistent with the previously reported research; therefore, it is evident that GA is a positive regulator of stem elongation and buds and branches development. Our study further showed that GA-treated plants showed increased stem diameter, which suggested that GA is significantly involved in the secondary growth of *J. curcas*. Melanie Mauriat et al., (2011) reported that the expression of the GA biosynthesis gene in *Populus tremula* and *Nicotiana tabacum* enhanced the internode length and stem diameter [39]. Reports show that GA enhances the elongation and division of xylem and fiber cells in the vascular bundle region and increases cambium activity which increases stem diameter [44–46]. Compared to IAA and NPK, GA showed an increased number of leaves, branches, and stem diameter (Figure 2). Leaf number and growth are major determining factors contributing to shoot biomass and yield production. Leaves and branches number are regulated and controlled genetically and depend on developmental stage and species [47]. However, studies reported that growth regulators regulate leaves and branches development. In *Gladiolus grandifloras,* GA enhanced the leaves number and sprouting emergence compared with control plants [48]. Although there is a lack of data about GA induction of proximate composition in the *J. curcas* plant; however, it is reported that GA application increases carbohydrate contents in *Phalaenopsis* apex [49]. Reports show that GA accumulation enhances certain enzyme activity and increases cell wall plasticity which enhances membrane permeability and facilitates the uptake of mineral nutrients and transport of photosynthates [50–53]. The GA-facilitated uptake of nutrients and photosynthates transportation promotes plant growth and development. GA accelerates the cell cycle and starch hydrolysis to provide energy and a carbon skeleton for the synthesis of soluble sugar (carbohydrates) and other metabolites [54]. Our results showed that GA application increased protein contents in *J. curcas* which is in line with the results of Satendra Singh et al., (2014), they determined that GA3 significantly increased the total amino acid (protein contents) of *Phaseolus vulgaris* L. plant [55]. Taken together with these results, it is revealed that the exogenous application of GA enhances *J. curcas* growth and development mediated through the regulation of morphological, physiological, and biochemical parameters.

*J. curcas* is a biodiesel-producing plant and possesses great characteristics of growing on barren land, low water, and harsh climatic condition. It is a nutrient-reactive plant and its requirement for nutrients varies with the soil fertility [56]. N, P, K, Zn, and B are the main fertilizers needed for the full-size growth of the plant to produce seeds [57,58]. N and P are the main nutrients that affect significantly the seed yield of *J. curcas* [59]. Researchers have suggested that, in the developmental stage of *J. curcas*, NPK is needed to build up the plant architecture such as root, stem, leaves, flowers, and seeds [57]. Analysis of our results revealed that NPK application to *J. curcas* showed significant differences in root shoot length, fresh dry weight, number of leaves, and branches. Nitrogen has a key role in chlorophyll biosynthesis and increases the rate of photosynthesis which results in increases in the dry biomass of the plant [60]. Similarly, phosphorus is involved in energy metabolism and photosynthesis while potassium plays a key role in carbohydrates and protein metabolism during plant growth [61,62]. In general, NPK fertilizer is easily absorbed by plants which play an important role in growth by supporting vegetative development including leaf, stem, and root development. Nitrogen is a building factor of protein, which greatly affect plant biomass [63]. Phosphorus is also used in protein and fat biosynthesis and transform adenosine diphosphate into adenosine triphosphate to generate energy [64]. Potassium enable CO2 during photosynthesis to enter through stomata, photosinthate transport, water, sugar, and protein and sugar synthesis [63]. Potassium availability increases energy that results into growth and development regulation [65]. Compared to GA and IAA, NPK 60 g showed the highest significant increase in total protein, carbohydrate, fat, moisture, and ash content percentage (Figure 3), which indicates that NPK higher concentration increases the proximate composition of *J. curcas* plant. Ali Sher Chandio et al., (2016) investigated that the highest amount of NPK application per hector showed maximum growth, fruit, and seed oil yield in the *J. curcas* plant [56]. Another research also revealed that the highest value of dry matter and seed fatty acid were obtained from NPK-treated *J. curcas* plants compared with non-treated plants [12]. Nitrogen has a key importance in plant growth than the other material as it plays a central role in many physiological and biochemical processes in plants [66]. It is a basic part of the structure of chlorophyll, protein, fats, and nucleic acid [66]. After nitrogen, phosphorus is an important macronutrient that affects plant growth. Without enough phosphorus, it is difficult for a plant to attain maximum growth and development as it has a key role in the storage and transfer of energy in plants [67,68]. Whereas potassium also plays a key role in metabolism and its application influence *J. curcas* seed oil and fatty acids [12]. All three macronutrients (NPK) are important growth inducers of plants but instead of using them separately, they are more efficient when used in combination [34]. Among the NPK, N plays a more important role in enhancing agriculture production by promoting chlorophyll, soluble protein, proline contents, and promoting fiber yield [69,70]. Compared to previous research, our study also determined that NPK application promoted *J. curcas* growth and development by enhancing morphological, physiological, and biochemical regime.

#### **4. Materials and Methods**

#### *4.1. Experimental Design and Material Used*

In the present study, six months old seedlings of equal size were collected from a common dealer in Multan, Pakistan. The seedling weight, stem width, and root length were different; therefore, the initial weight, stem width, and root length were measured and recorded. The whole experiment was conducted in pots and single seedlings were grown in each pot filled with an equal amount of soil and treated with hormones and fertilizer separately. The hormones Gibberellic acid GA3 (GA) and Indole acetic acid (IAA) were applied as 100 ppm and 250 ppm each to every single pot and NPK was applied 30 g and 60 g to each pot separately (mixed in the soil). Control plants were treated with only water. The plants were watered after two weeks with an equal amount of water. The data were collected after each week for four months and the whole experiment was conducted in triplicate. In fertilizer, nitrogen was 28%, phosphorus was 18% and potassium was 16%.

#### *4.2. Parameters Studied*

In the present study, three main parameters were studied, i.e., physiological, morphological, and biochemical parameters. In physiological parameters, root, shoot length, and fresh and dry weight of root and stem were measured. In morphological parameters, the number of leaves, number of branches, and stem diameter were studied. In biochemical parameters, carbohydrates, proteins, lipids, moisture, and ash contents were studied. The data were taken every week until four months.

#### *4.3. Proximate Composition Analysis*

To determine the effect of hormones and fertilizers on leaf moisture contents, the fresh leaves were randomly collected in triplicate from each treatment. The collected leaves were washed, weighed, and recorded as fresh weight (FW). The selected leaves were at 105 ◦C in the oven for 3 h. The fully dried leaves were again weighed and recorded as dried weight (Dw). The moisture content was calculated as the following:

$$\text{Moisture content} \left(\% \right) = \left(\text{Fw} - \text{Dw}\right) \times 100$$

To measure the ash content percentage, we collected 1 g of leaves and washed them to remove the contamination and the tissue dried. The leaves were burned in the muffle furnace at 450–550 ◦C to remove water and other volatile substances. The leaves samples were weighed before and after burning in the muffle furnace. The ash content was calculated as follows:

$$\text{Ash} \left( \% \right) = \left( \text{M}\_{\text{ash}} / \text{M}\_{\text{dry}} \right) \times 100 \text{ J}$$

where Mash is the mass of fresh ash and Mdry is the mass of dry ash.

To determine fat contents, we followed the method used by Aurea M. Almazan and Samuel O. Adeyeye (1998) [71]. About 1 g of fresh leaves were collected from *J. curcas* and dried in the oven at 60 ◦C for 15 h and then ground into fine powder. The powder was homogenized with hexane and extracted fat by using the Soxhlet method (AOAC, 1990). The extract was treated with methanolic sodium and methanolic boron trifluoride to convert the fatty acid into methyl ester, using the method followed by Paquot and Hautfenne in 1987 [72]. The fatty acid methyl ester in the hexane layer was dried at 90 ◦C by passing the nitrogen gas on its surface. The solution of hexane and methyl ester was filtered by a 45 μm filter and injected 1 μL into a Shimadzu gas chromatograph with a Perkin-Elmer PE-WAX capillary column (30 m × 025 mm) and flame ionization detector. The initial column temperature was 80 ◦C and increased to 260 ◦C kept for 7 min with the flow rate of gas (H2) being 30 cm/s. The percentage of fat was calculated based on the total area of fatty acids. All these determinations were carried out according to the Association of Official Analytical Chemists (AOAC, 1990). For crude protein determination, the Micro Kjeldahl method was followed according to the AOAC international [73]. About 1 g of the leaf samples was collected randomly then ground into a fine powder and digested with 15 mL H2SO4 by heating in the presence of K2SO4 and selenium using a heating block at 420 ◦C for 2 h. The digested sample was then neutralized by adding NaOH, to convert ammonium sulfate into ammonia, which is further distilled off and collected in a flask, and added boric acid was to form ammonium borate. The residual boric acid was further titrated with H2SO4 with the use of an endpoint indicator to determine the total nitrogen contents. The amount of total nitrogen in the raw material was multiplied by the traditional conversion factor of 6.25 and the specific conversion factor [74]. Carbohydrate content was determined by calculating the difference between the sum of all the proximate compositions from 100% [75].

#### *4.4. Soil Analysis*

To find the nutrient deficiency in the soil, the soil was analyzed at the Barani Agriculture Research Center (BARS), Kohat. Two soil samples were randomly selected from the different sites of the field. The BARS report of soil is presented in Table 1.


**Table 1.** Parameters studied in soil analysis.

#### *4.5. Statistical Analysis*

All experiments were performed in triplicate, and the data from each replicate were pooled. Data were analyzed using one-way ANOVA with Bonferroni post hoc tests (\* shows *p* < 0.05, \*\* shows *p* < 0.01, and \*\*\* shows *p* < 0.001 significant difference). A completely randomized design was used to compare the mean values of different treatments. Data were graphically plotted, and statistical analyses were performed using the GraphPad Prism software (version 5.01, GraphPad, San Diego, CA, USA).

#### **5. Conclusions**

This study demonstrated that the exogenous application of phytohormones and fertilization of the *J. curcas* plant promoted growth and development. The results confirmed that different concentrations of GA, IAA, and NPK induced various parameters differentially. GA 250 ppm increased shoot root length, root fresh and dry weight, branches number, and stem diameter while GA 100 ppm increased leaves number as compared to IAA and NPK. While IAA 250 ppm increased stem diameter, and NPK increased proximate compositions compared to hormones. Our study concluded that optimum hormones and NPK level is essential for the efficient promotion of morphological, physiological, and biochemical aspects of the *J. curcas* plant. Furthermore, the hormones were used for scientific validation, however on commercial basis application of phytohormones are less economical than the fertilizers. Our current study opens an important research area for the future study to investigate, how to improve the indigenous phytohormones of *J. curcas* to improve its growth and development in the barren land.

**Author Contributions:** Conceptualization, K.-M.K., W.M. and R.J.; methodology, R.J., S.A. (Sajjad Asaf) and M.A.; Formal analysis and software, R.J., S.A. (Saleem Asif) and M.K.; writing—original draft, R.J., K.-M.K. and W.M.; review, visualization and supervision, K.-M.K. and W.M. All authors have read and agreed to the published version of the manuscript.

**Funding:** This research received no external funding.

**Data Availability Statement:** The data presented in this study are available on request from the corresponding author.

**Acknowledgments:** This work was supported by a grant from the New Breeding Technologies Development Program (Project No. PJ016531012022), Rural Development Administration, Korea.

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

#### **References**


## *Article Brassica napus* **Transcription Factor** *Bna.A07.WRKY70* **Negatively Regulates Leaf Senescence in** *Arabidopsis thaliana*

**Tiantian Liu, Yuxin Li, Chang Wang, Da Zhang, Jiajia Liu, Mingyuan He, Mingxun Chen \* and Yuan Guo \***

National Yangling Agricultural Biotechnology & Breeding Center, Shaanxi Key Laboratory of Crop Heterosis, and College of Agronomy, Northwest A&F University, Yangling 712100, China

**\*** Correspondence: cmx786@nwafu.edu.cn (M.C.); guoyuan2109@163.com (Y.G.)

**Abstract:** Leaf senescence is the final stage of leaf development and is essential for storage properties and crop productivity. WRKY transcription factors have been revealed to play crucial roles in several biological processes during plant growth and development, especially in leaf senescence. However, the functions of *Brassica napus* WRKY transcription factors in leaf senescence remain unclear. In the present study, *Bna.A07.WRKY70*, one paralogue of *Brassica napus WRKY70*, was cloned from the *B. napus* cultivar "Zhongshuang11 (ZS11)". We found that Bna.A07.WRKY70 contains a highly conserved WRKY domain and is most closely related to *Arabidopsis thaliana* WRKY70. The subcellular localization and transcriptional self-activation assays indicated that Bna.A07.WRKY70 functions as a transcription factor. Meanwhile, RT-qPCR and promoter-GUS analysis showed that *Bna.A07.WRKY70* is predominantly expressed in the leaves of *B. napus* and rosette leaves of *A. thaliana*. In addition, our results demonstrated that ectopic expression of *Bna.A07.WRKY70* in *A. thaliana wrky70* mutants could restore the senescence phenotypes to wild-type levels. Consistently, the expression levels of three senescence-related marker genes of *wrky70* mutants were restored to wild-type levels by ectopic expression of *Bna.A07.WRKY70*. These findings improve our understanding of the function of *Bna.A07.WRKY70* in *B. napus* and provide a novel strategy for breeding the new stay-green cultivars in rapeseed through genetic manipulation.

**Keywords:** *Bna.A07.WRKY70*; leaf senescence; *Arabidopsis thaliana*; *Brassica napus*

#### **1. Introduction**

Rapeseed (*Brassica napus* L., AACC, 2n = 38), an allotetraploid species, originated from spontaneous hybridization between two diploid *Brassica* species: *Brassica rapa* (AA, 2n = 20) and *Brassica oleracea* (CC, 2n = 18) [1]. It is a major oilseed crop grown worldwide for the production of edible oil in the human diet, livestock feed, and industrial materials [2]. Therefore, there is important social and economic significance for studying its associated biological processes, including leaf development. The leaf is the primary organ of photosynthesis and can produce nutrition and gather energy during plants' growth and maturation stages. Leaf senescence, as a type of programmed cell death (PCD), is the terminal stage of leaf development. During leaf senescence, the chloroplast first starts disassembling, and is followed by a loss of chlorophyll together with the catabolism of macromolecules such as protein, lipids, nucleic acids, and RNA [3]. By general catabolism, cellular materials are converted into easily exportable nutrients, which from senescing leaves were subsequently transported to reproductive and developing structures [4]. Consequently, leaf senescence is a critical process for crop fitness and is particularly essential for the optimization of crop productivity. Generally, leaf senescence is influenced by various external environmental and endogenous factors. The environmental cues that affect leaf senescence include high temperature, light signals, drought, and biotic stress [5–7]. The endogenous factors include the accumulation of reactive oxygen species (ROS), variation of plant hormones, and, most importantly, regulation of multiple senescence-associated genes [8–11]. Therefore, mining the key genes regulating the leaf senescence process is of great importance in rapeseed.

**Citation:** Liu, T.; Li, Y.; Wang, C.; Zhang, D.; Liu, J.; He, M.; Chen, M.; Guo, Y. *Brassica napus* Transcription Factor *Bna.A07.WRKY70* Negatively Regulates Leaf Senescence in *Arabidopsis thaliana*. *Plants* **2023**, *12*, 347. https://doi.org/10.3390/ plants12020347

Academic Editor: Vagner A. Benedito

Received: 14 December 2022 Revised: 3 January 2023 Accepted: 5 January 2023 Published: 11 January 2023

**Copyright:** © 2023 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (https:// creativecommons.org/licenses/by/ 4.0/).

The WRKY proteins are one of the largest and most important superfamilies of transcription factors (TFs) in plants. WRKY transcription factors encompass a core motif WRKYGQK (a highly conserved WRKY domain) at the N-terminus and an atypical Zinc finger motif at the C-terminus [12]. On the basis of both the number of WRKY domains and the features of the Zn-finger motif in their evolutionary history, the WRKY TFs can be divided into three different groups (I, II, and III). Group I contains two WRKY domains and a finger motif whose pattern is conserved zinc ligands (C–X4-5–C–X22-23–H–X1–H), which is the same as the Zinc finger motif of group II, but there is only one WRKY domain in group II. Instead of a C2–H2 pattern in group I and II, group III contains a pattern of C2–HC (C–X7–C–X23–H–X1–C) zinc finger-like motif and have one WRKY domain [13,14]. All three groups' members of WRKY TFs have been demonstrated to interact with the specific DNA *cis*-acting element W-box (C/TTGACT/C) in the promoter regions of downstream genes and further regulate their expression [15]. In recent decades, experimental evidence has shown that WRKY proteins act as key regulators widely involved in various plant growth and development processes, such as leaf senescence, growth of roots [16], stem elongation [17], and multiple biotic and abiotic stressors [18,19]. In *Arabidopsis thaliana*, *WRKY53* acts in a complex transcription factor signaling network regulating leaf senescencespecific gene expression [20], WRKY45 interacts with the DELLA protein RGL1 to positively regulate age-triggered leaf senescence [21], and *WRKY71* mediates ethylene (ET) signaling and synthesis to hasten leaf senescence [22].

*WRKY70* belongs to WRKY III TFs and has been reported in response to several developmental and physiological processes in diversified species [23]. In *Arabidopsis*, *WRKY70* acts as a negative regulator of leaf senescence, with gradually increasing expression during leaf development [3], and *WRKY70* is also crucial in plant defense against pathogens, controlling the cross-talk of salicylic acid (SA) and jasmonic acid (JA) signaling in plant defense [23,24]. Moreover, *WRKY70* is an important signaling component that is positively involved in brassinolide (BR)-regulated growth and negatively involved in drought responses by inhibiting drought-responsive genes [25]. In chickpeas, *WRKY70* was reported to regulate the expression of a chickpea HD-Zip transcription factor *CaHDZ12*, which improved tolerance to osmotic stresses under drought and salinity stress, and increased sensitivity to abscisic acid (ABA) in transgenic tobacco and chickpea [26]. In addition, it was suggested that GhWRKY70D13 negatively regulates cotton's resistance to *Verticillium dahliae* mainly through its effect on ET and JA biosynthesis and signaling pathways [27]. A recent study demonstrated that *TaWRKY70* positively regulates *TaCAT5* by directly binding to the *TaCAT5* promoter to enhance Cd tolerance in transgenic *Arabidopsis* [28]. In *B. napus*, the *BnWRKY70* knockout plants by CRISPR/Cas9 system enhanced *Sclerotinia sclerotiorum* resistances, while overexpression of *BnWRKY70* reduced resistance to *S. sclerotiorum* [29]. However, the roles of WRKY proteins in *B. napus* in the regulation of leaf senescence remain unclear.

In the current study, *Bna.A07.WRKY70*, one of the *AtWRKY70* orthologues in *B. napus*, was isolated and functionally characterized. We found that *Bna.A07.WRKY70* functioned as a TF and was specifically expressed in the leaves in *A. thaliana* and *B. napus*. We also demonstrated that ectopic expression of *Bna.A07.WRKY70* in the *A. thaliana wrky70* mutant restored the leaf senescence rate and chlorophyll content and greatly altered the expression of three senescence-related genes in this mutant. Our results may indicate that *Bna.A07.WRKY70* functions as a negative regulator of leaf senescence in *Arabidopsis*, which might reveal a conserved role of WRKY70 proteins in regulating leaf senescence between *A. thaliana* and *B. napus*.

#### **2. Results**

#### *2.1. Sequence Analysis of BnaWRKY70 Paralogs*

In the *B. napus* cultivar "Zhongshuang11 (ZS11)", six paralogs of BnaWRKY70 were predicted in BnPIR (http://cbi.hzau.edu.cn/bnapus/index.php, accessed on 9 September 2022) and were designated Bna.A07.WRKY70 (BnaA07G0195100ZS), Bna.C06.WRKY70 (BnaC06G0198900ZS), Bna.A04.WRKY70 (BnaA04G0035900ZS), Bna.C08.WRKY70 (BnaC08G0362900ZS), Bna.A09.WRKY70 (BnaA09G0519800ZS), and Bna.C04.WRKY70 (BnaC04G0308100ZS). With the multiple sequence alignment, we found that the WRKY70 protein from *B. napus* and *A. thaliana* possessed highly conserved WRKY domains, including WRKYGQ/KK core motif and a pattern of C2–HC zinc finger-like motif at the C-terminus (Figure 1A). Among them, Bna.A07.WRKY70 was predicted to share the highest identity in the amino acid sequence with the AtWRKY70 protein (66.01%) (Figure S1). A phylogenetic analysis was performed to investigate the evolutionary relationships between Bna.A07.WRKY70 and 20 WRKY70 proteins from seven plant species, including *A*. *thaliana*, *B*. *napus*, *B*. *rapa*, *Glycine max*, *Zea mays*, *Oryza sativa*, and *Setaria italic*. As illustrated in Figure 1B, Bna.A07.WRKY70 is most closely related to the WRKY70 protein from *B. rapa* (NP\_001288847.1) and *A. thaliana* (AtWRKY70). These results suggested preliminarily that Bna.A07.WRKY70 might have similar functions as AtWRKY70.

#### *2.2. Subcellular Localization and Transcriptional Activity of Bna.A07.WRKY70*

For subcellular localization, Bna.A07.WRKY70 was expressed in tobacco (*Nicotiana benthamiana*) leaf cells as a recombinant protein fused to a green fluorescent protein marker. The fluorescence signal was detected in the nucleus by laser scanning confocal microscopy (Figure 2A), suggesting that Bna.A07.WRKY70 might function as a transcription factor.

To further characterize Bna.A07.WRKY70 function, we investigated whether Bna.A07. WRKY70 has transcription activation activity in yeast cells. The empty vector pGBKT7 as negative control and fusion construct (*pBD-Bna.A07.WRKY70*) were transformed separately into Y2HGold yeast cells, which were cultured on SDO (SD/-Trp) and TDO (SD/-Trp/- His/-Ade) medium. As shown in Figure 2B, on SDO (SD/-Trp) medium, all yeast transformants could grow normally, indicating that the constructs were transformed successfully into the Y2HGold yeast cells. Instead, on TDO (SD/-Trp/-His/-Ade) medium, the empty vector pGBKT7 did not grow, but yeast cells with Bna.A07.WRKY70 fusion constructs grew well, which demonstrated that Bna.A07.WRKY70 could activate the expression of the reporter genes. Given these findings, the Bna.A07.WRKY70 was testified to function as a transcription activator.

#### *2.3. Analysis of Bna.A07.WRKY70 Expression Pattern*

We further investigated the spatiotemporal expression pattern of *Bna.A07.WRKY70* by analyzing the relative abundance of the mRNA in various tissues of *B. napus* cultivar "ZS11" using quantitative reverse transcription PCR (RT-qPCR). The results showed that *Bna.A07.WRKY70* was widely expressed in different organs of *B. napus*, with higher expression in leaves, moderate in stems, flowers, and roots but low in developing seeds (Figure 3A). To comprehensively investigate the spatiotemporal expression pattern of *Bna.A07.WRKY70* in *A. thaliana*, we obtained 16 *pBna.A07.WRKY70:GUS* in wild-type background independent lines of *A. thaliana* and the one representative line was used for promoter-GUS analysis because of similar GUS staining patterns in most lines. Consistent with the RT-qPCR data in *B. napus*, promoter-GUS activity staining was predominantly detected in rosette leaves of *A. thaliana* (Figure 3D) and then was also slightly detected in other organs of *A. thaliana*, including stems (Figure 3C), roots (Figure 3B), and flower abscission zones (Figure 3E). Conversely, it was not detected in the embryo (Figure 3G) and siliques (Figure 3F). In summary, these observations suggested that *Bna.A07.WRKY70* might regulate a significant function during leaf development.


**Figure 1.** Protein sequence and phylogenetic analyses of WRKY70 proteins. (**A**) Protein sequence

alignment of WRKY70 from *A. thaliana* and *B. napus* was carried out using the MUSCLE program (http://www.ebi.ac.uk/Tools/msa/muscle/, accessed on 12 September 2022). Asterisks indicate non-conservative differences. The WRKY domain 125–185, as indicated by the Conserved Domain Search program (http://www.ncbi.nlm.nih.gov/Structure/cdd/wrpsb.cgi, accessed on 12 September 2022), is underlined. The highly conserved core sequence WRKYGQK in the WRKY domain is represented by a red box, together with the C and H residues in the CCHC zinc-finger-like motif indicated by a downward black triangle. (**B**) Phylogenetic analysis of Bna.A07.WRKY70 with 20 other WRKY70 proteins from seven plant species, including AtWRKY70 (*Arabidopsis thaliana*); Bna.A07.WRKY70 (BnaA07G0195100ZS), Bna.C06.WRKY70 (BnaC06G0198900ZS), Bna.A04.WRKY70 (BnaA04G0035900ZS), Bna.C08.WRKY70 (BnaC08G0362900ZS), Bna.A09.WRKY70 (BnaA09G0519800ZS), and Bna.C04.WRKY70 (BnaC04G0308100ZS (*Brassica napus*); NP\_001288821.1, XP\_033146547.1, XP\_009104000.1, and NP\_001288847.1 (*Brassica rapa*); XP\_015637594.1 (*Oryza sativa*); XP\_035823367.1, XP\_035823368.1, and XP\_035818584.1 (*Zea mays*); XP\_012703156.1, XP\_022682393.1, XP\_022682403.1, and XP\_022682410.1 (*Setaria italica*) and XP\_006595434.1 and XP\_014621824.1 (*Glycine max*). A neighbor-joining tree (Jones–Taylor–Thornton model) with 1000 replicates of bootstrap analysis was generated by MEGA7. Bootstrap values are indicated at the nodes, and the accession numbers of the species are labeled on the phylogenetic tree.

**Figure 2.** Transcription factor characterization of Bna.A07.WRKY70. (**A**) Subcellular localization of Bna.A07.WRKY70 protein fused with GFP (*35S:GFP-Bna.A07.WRKY70*) in tobacco leaves (*Nicotiana benthamiana*). mCherry, a nuclear-localized protein fused with a red fluorescent protein; merge, merge of mCherry, GFP, and bright field images. (**B**) Transcriptional activation assays of Bna.A07.WRKY70 in yeast. BD: empty vector that contains GAL4 DNA-binding domain, *BD-Bna.A07.WRKY70*: cDNAs encoding of *Bna.A07.WRKY70* transcripts were separately cloned into the pGBKT7/BD vector containing the GAL4 DNA binding domain, which transformed into the yeast strain Y2HGold, SDO: ability of yeast transformants to grow on medium lacking Trp, TDO: ability of yeast transformants to grow on medium lacking Trp, His and Ade indicates transcriptional activation. 1, 10−1, 10−2: the transformed strains were spotted on plates by diluting to different concentrations. The images show representative results from more than four independent yeast transformants.

#### *2.4. Bna.A07.WRKY70 Negatively Regulates Leaf Senescence in A. thaliana*

To further explore the function of *Bna.A07.WRKY70* on leaf development, we introduced the construct *35S:Bna.A07.WRKY70-GFP* into *A. thaliana wrky70* mutant (Figure 4A). Twelve independent T1 transgenic plants were generated using hygromycin selection, and five independent T3 homozygous transgenic lines *wrky70 35S:Bna.A07.WRKY70-GFP* were selected and confirmed by PCR amplification with the specific primers 35S-F/Bna.A07.WRKY70- GFP-BamHI-R (Figure 4B; Supplementary Table S1). Of these lines, three representative ones, *wrky70 35S:Bna.A07.WRKY70-GFP #4*, *#6*, and *#12*, with a relatively high expression

level (Figure 4C), were selected for the follow-up experiment. As illustrated in Figure 5A,B, the loss-of-function mutants *wrky70* exhibited markedly yellowing of leaves at 35 DAG (days after germination) and indicated earlier senescence compared to wild-type plants, which is in line with previous findings [3]. Interestingly, we found that the *A. thaliana wrky70* mutant leaves were smaller than wild-type plants. Ectopic expression of *Bna.A07.WRKY70* fully restored the rate of leaf senescence to wild-type levels in *Arabidopsis wrky70* mutants (Figure 5A). Furthermore, by arranging the rosette leaves of 35-day-old Col-0, *wrky70* mutant, and transgenic plants (*#4*, *#6*, *#12*) according to their age from older to younger, we found that three *Bna.A07.WRKY70* transgenic lines in the *wrky70* background delayed the premature senescence of leaves and rescued the phenotype of leaves relatively smaller in size compared to wild-type plants (Figure 5B). The results of the chlorophyll content indicated that the chlorophyll content of the *wrky70* mutant intensified degradation from the fifth week, but the chlorophyll content of *Bna.A07.WRKY70* transgenic lines were in keeping with that of Col-0 and clearly higher than that of the *wrky70* mutant (Figure 5C).

**Figure 3.** Analysis of the *Bna.A07.WRKY70* expression pattern. (**A**) RT-qPCR analysis of the *Bna.A07.WRKY70* expression in various tissues of *B. napus* cultivar "ZS11". The RT-qPCR result was normalized against the expression of *BnACTIN7* as an internal control. Values are means ± SD (n = 3). Error bars denote standard deviations. (**B**) to (**G**), Histochemical GUS staining in 35-day-old *ProBna.A07.WRKY70:GUS* transgenic *Arabidopsis* plants. (**B**) Roots (bar = 2 mM); (**C**) Stems and leaves (bar = 2 mM); (**D**) rosette leaves (bar = 2 mM); (**E**) Flowers (bar = 2 mM); (**F**) siliques 12 days after pollination (bar = 2 mM); (**G**) Developing seeds 12 days after pollination (bar = 200 μM).

In order to further confirm that the *Bna.A07.WRKY70* regulates the progress of leaf senescence in *A. thaliana*, we assessed the transcript levels of representative genes relating to senescence in the fifth and sixth rosette leaves of *A. thaliana* wild-type, the *wrky70* mutant, *wrky70 35S:Bna.A07.WRKY70-GFP* plants at 35 DAG. Compared to the wild type, the results showed that the expression of the senescence-related gene *AtSAG13* (*senescence-associated gene 13*) and *AtSEN1* (*senescence-associated gene 1*) were significantly increased, while the expression of photosynthesis-related *AtCAB1* gene (*chlorophyll a/b-binding protein*) was significantly decreased in 35-day-old *wrky70* mutant plants (Figure 6). However, when the *35S:Bna.A07.WRKY70-GFP* was introduced into the *wrky70* mutant, we found that the transcript abundance of these three senescence-related marker genes, including *AtSEN1*, *AtCAB1*, and *AtSAG13*, was restored to wild type levels. In brief, all results containing the premature senescence phenotype, chlorophyll content, and the expression of senescenceassociated marker genes together revealed that *Bna.A07.WRKY70* may negatively regulate the leaf senescence by adjusting the expression of senescence genes in *A. thaliana* and play a similar role with *AtWRKY70* in *A. thaliana*.

**Figure 4.** Molecular characterization of *wrky70 35S:Bna.A07.WRKY70*-*GFP* transgenic plants. (**A**) Schematic diagram of constitutive expression cassette of the *Bna.A07.WRKY70* gene in the binary vector pCAMBIA-1300 used for plant transformation. RB, right border; LB, left border; 35S Pro, CaMV 35S promoter; Nos, nopaline synthase terminator; CaMV35S, CaMV 35S promoter; *hptII*, hygromycin resistance gene. (**B**) PCR-based DNA genotyping of *wrky70 35S:Bna.A07.WRKY70*- *GFP* transgenic plants using specific primers of 35S\_P/Bna.A07.WRKY70-GFP-BamHI-R. Const, *35S:Bna.A07.WRKY70*-*GFP* construct. Col-0 and *wrky70* indicate *A. thaliana* wild type and mutant plants, respectively. (**C**) Expression analysis of *Bna.A07.WRKY70* in *wrky70 35S:Bna.A07.WRKY70*-*GFP* transgenic plants using RT-qPCR. The expression level was normalized against the expression of *AtACTIN7*, which was used as an internal control. Values are the means ± SD (n = 3). Error bars indicate standard deviation. # indicates the transgenic lines.

**Figure 5.** Effects of *Bna.A07.WRKY70* overexpression in the *wrky70* mutant background on leaf senescence in *A. thaliana*. (**A**) The whole plant phenotypes of leaf senescence in the wild type (Col-0), *wrky70* mutant, and *wrky70 35S:Bna.A07.WRKY70-GFP* transgenic plants. The images were taken 35 days after germination (DAG). Bar = 1 cM. (**B**) Phenotype of rosette leaves in 35-day-old plants, excised leaves are arranged according to age, from older to younger. Bar = 1 cM. (**C**) Comparisons of chlorophyll content of the fifth to sixth rosette leaves among wild-type (Col-0), *wrky70* mutant, and *wrky70 35S:Bna.A07.WRKY70* transgenic plants at the indicated ages. Values are means ± SD (n = 3). Asterisks indicate significant differences from wild-type (two-tailed paired Student's *t*-test, *p* ≤ 0.05). Error bars indicate standard deviation. # indicates the transgenic lines.

**Figure 6.** Expression analysis of leaf senescence marker genes in the rosette leaves among wild-type (Col-0), *wrky70* mutant, and *wrky70 35S:Bna.A07.WRKY70-GFP* transgenic plants at the 35 DAG, as measured by RT-qPCR. Expression levels were normalized to the expression of the internal reference gene, *AtACTIN7*. Values are means ± SD (n = 3). Error bars indicate standard deviations. Asterisks indicate statistically significant differences from wild type plants (two-tailed paired Student's *t*-test, *p* ≤ 0.05). # indicates the transgenic lines.

#### **3. Discussion**

Leaf senescence is an indispensable portion and spans the latter half of leaf development. It is a highly intricate process regulated by multiple pathways [30]. As previously reported, the three largest groups of transcription factors, WRKY, NAC, and MYB superfamilies, are responsible for modulating transcriptional changes during leaf senescence [31], in which the *AtWRKY70* has already been confirmed with a high level of expression in the late stage of leaf development and functions as an essential repressor during leaf senescence in *A. thaliana* [32]. However, the roles of *WRKY70* transcription factors during leaf development in *B. napus* were lacking.

It has been widely known that *B. napus* was formed 7500 years ago by natural hybridization between *B. rapa* and *B. oleracea* [33]. *B. napus*, and diploid parental species *B. rapa* and *B. oleracea*, are believed to share a common ancestor with *A. thaliana*, a fact that has favored the transfer of knowledge from *Arabidopsis* to *B*. *napus*. As an allopolyploid plant, a large number and a high frequency of chromosome variation activities were identified, such as duplication, rearrangement, fusion, and deletion in the evolution processes of *B. napus*, which makes the genomics of *B*. *napus* more complicated. Generally, a single *Arabidopsis* gene is represented by two to eight paralogs in the *B. napus* genome [34,35]. Accordingly, six paralogues (Bna.A07.WRKY70, Bna.C06.WRKY70, Bna.A04.WRKY70, Bna.C08.WRKY70, Bna.A09.WRKY70, and Bna.C04.WRKY70) were found in the *B. napus* genome (Figure 1). In the WRKY transcription factor family, the WRKY domain is the major determinant of DNA-binding and specifically binds DNA *cis*-acting element W-box (C/TTGACT/C). Our results showed that all six Bna.WRKY70 had the WRKY protein domain containing the WRKYGQ/KK core motif and a pattern of C2–HC zinc finger-like motif at the C-terminus (Figure 1A). In the present study, among these BnaWRKY70 paralogues, Bna.A07.WRKY70, which was predicted to have the highest identity of protein sequence and the WRKY central conserved domains with AtWRKY70 (Figure 1), was cloned from the *B. napus* cultivar "ZS11" and functionally characterized. Bna.A07.WRKY70 was located in the nucleus of tobacco leaf cells, and we further demonstrated that Bna.A07.WRKY70 could activate the expression of the reporter genes in yeast cells (Figure 2). These results suggested that Bna.A07.WRKY70 functions as a transcription activator. Additionally, the *Bna.A07.WRKY70* transcript was broadly present in different vegetative tissues, with the highest levels observed in leaves (Figure 3), suggesting that *Bna.A07.WRKY70* might regulate a significant function during leaf development. Ectopic expression of *Bna.A07.WRKY70* in the background of *A. thaliana wrky70* mutants significantly delayed the senescence of leaves and restored the chlorophyll content to the wild type level (Figure 5). Moreover, the

expression of senescence-associated genes (*AtSEN1*, *AtSAG13*, and *AtCAB1*) was clearly regulated by *Bna.A07.WRKY70* during leaf senescence. Thus, these results may indicate that *Bna.A07.WRKY70* functions as a negative factor in leaf senescence as the *AtWRKY70*.

During leaf senescence, the leaves turned yellow, resulting in photosynthesis deficiency and beginning with chloroplast dismantling, followed by degradation of chlorophyll and chlorophyll-protein complexes. Meanwhile, leaf senescence is accompanied by decreased expression of genes related to photosynthesis and protein synthesis and increased expression of senescence-associated genes (SAGs) [36]. Consistently, our results demonstrated that compared to the wild type, the expression of the photosynthesis-related *AtCAB1* gene was decreased, and the expression of senescence-related gene *AtSAG13* and *AtSEN1* were increased in *A. thaliana wrky70* mutant plants. The expression of these three marker genes was rescued to wild-type in *wrky70 35S:Bna.A07.WRKY70-GFP* transgenic plants, which proved that *Bna.A07.WRKY70* indeed delayed the leaf senescence during plant senescence by affecting the expression of these three senescence genes in *A. thaliana*. Leaf senescence was widely influenced by a variety of external and internal factors, including environmental stresses and phytohormones. Recently, key gene regulatory networks comprising these TFs have been identified, indicating that leaf senescence is controlled by multiple cross-linking pathways, many of which are associated with stress response signaling [37–39]. *Arabidopsis WRKY71* was reported that it is able to directly upregulate the ethylene signaling pathway genes *EIN2* (*ethylene insensitive2*) and *ORE1* (*oresara1*) and promote ethylene synthesis by directly activating the *ACS2* gene to accelerate leaf senescence in *Arabidopsis* [22]. The cotton (*Gossypium hirsutum* L.) *GhWRKY91* directly targets *GhWRKY17*, a gene associated with ABA signals and reactive oxygen species (ROS) production to negatively mediate leaf senescence and provide a foundation for further functional studies on natural and stress-induced leaf senescence [40]. *OsWRKY53* of rice, as a positive regulator, repressed the transcript of ABA catabolic genes (*OsABA8ox1* and *OsABA8ox2*) by directly binding to their promoters to promote ABA accumulation, and modulated ABA-induced leaf senescence [41]. In *Arabidopsis*, *AtWRKY70* transcript levels were more strongly reduced in *npr1* (*non-expressor of PR 1*) and *pad4* (*phytoalexin-deficient 4*) and completely abolished in *NahG* (salicylate hydroxylase gene) plants compared to wild-type at 40 days post germination, among which, the *NahG*, *pad4*, and *npr1* belonged to SA mutants and exhibited a delayed senescence phenotype [3]. These findings support the role of *AtWRKY70* as a senescence-associated gene and indicate a functional requirement of SA for its normal expression. Besides, the preceding research illustrated that the pathway of plant hormones could respond to numerous abiotic stresses; for instance, *GhWRKY17* from upland cotton modulated the increased sensitivity of plants to drought by reducing the level of ABA, and repressed transcript levels of ABA-inducible genes, including *AREB* (*ABA-responsive element binding*), *DREB* (*dehydration-responsive element binding*), *NCED* (*9-cis-epoxycarotenoid dioxygenase*), *ERD* (*early responsive to dehydration*) and *LEA* (*late embryogenesis-abundant protein*) under drought and salt stress conditions, indicating that *GhWRKY17* responds to drought and salt stress through ABA signaling and the regulation of cellular ROS production in plants [42]. With the above findings in mind, whether *Bna.A07.WRKY70* of *B. napus* adjusts the signaling pathways of phytohormone by combining with some key genes during the regulation of leaf senescence and responds to plant stress resistance mediated by the signaling pathways, can be explored further.

Interestingly, it has been reported that the *Arabidopsis wrky70* knockout mutants were slightly reduced in size compared to wild-type plants during the entire period of development in *A. thaliana* [3]. However, from another report, no obvious growth phenotype was observed in a single knockout mutant of *wrky70* compared with the wild-type *A. thaliana* [25]. In this study, our results found that the *A. thaliana wrky70* mutant leaves were smaller than wild-type plants, and the leaf size of the *wrky70* mutant was restored to wild type by the ectopic expression of *Bna.A07.WRKY70*. Whether or not *Bna.A07.WRKY70* plays a role in controlling the size of leaves requires further investigation. Overall, based on the above, the multiple functions and regulation network of *Bna.A07.WRKY70* still has great research potential.

#### **4. Materials and Methods**

#### *4.1. Plant Materials and Growth Conditions*

The *A. thaliana* wild-type ecotype Columbia (Col-0), the T-DNA mutant *wrky70* (SALK\_025198) in the Col-0 background obtained from Arashare (https://www.arashare. cn/index/, accessed on 20 October 2020), and *Brassica napus* L. cultivar "Zhongshuang 11 (ZS11)", were used in this study. The *A. thaliana* plants were grown in a growth chamber at 22 ◦C under a long day duration (LD, 16 h light/8 h dark) with moderate light intensity (160 μmol m−<sup>2</sup> s<sup>−</sup>1). The *B. napus* cultivar "ZS11" was first grown in the greenhouse at 22 ◦C with a long day duration for six weeks. For vernalization, the plants were transferred to a cold chamber at 4 ◦C under LD conditions. After vernalization, the plants were returned to the initial greenhouse conditions for 10 weeks until harvest.

#### *4.2. Protein Sequence and Phylogenetic Analysis*

The protein sequences of WRKY70 were obtained from the National Center for Biotechnology Information (NCBI) database (https://www.ncbi.nlm.nih.gov/, accessed on 9 September 2022) and the *B. napus* pan-genome information resource (BnPIR) database (http://cbi.hzau.edu.cn/bnapus/index.php, accessed on 9 September 2022). Protein sequence alignment was carried out using MUSCLE (http://www.ebi.ac.uk/Tools/msa/ muscle/, accessed on 12 September 2022). The conserved WRKY domain of Bna.A07.WRKY70 was indicated using the conserved domain search program in the National Center for Biotechnology Information (http://www.ncbi.nlm.nih.gov/Structure/cdd/wrpsb.cgi, accessed on 12 September 2022). The phylogenetic tree was conducted using the neighborjoining tree (Jones–Taylor–Thornton model) by MEGA7. Bootstrap analysis with 1000 replicates was performed to assess the statistical reliability of the tree topology.

#### *4.3. Gene Cloning and Plasmid Construction*

The full-length coding domain sequence (CDS) of *Bna.A07.WRKY70* (XP\_ 013648025.1) without stop codon was amplified by the specific primer designed in NCBI (https:// www.ncbi.nlm.nih.gov/tools/primerblast/, accessed on 20 October 2020). The total RNA was extracted from young leaves of the *B. napus* cultivar "ZS11" by the SteadyPure Plant RNA Extraction Kit (Accurate Biology, Changsha, China), and the RNA concentration was determined by spectrometry (Nano Drop; Thermo Scientific, Wilmington, MA, USA) (Supplementary Table S3) and quality was checked by 1% agarose gel electrophoresis. For cloning, first-stand cDNA was synthesized from total RNA using EasyScript One-Step gDNA Removal and cDNA Synthesis SuperMix (TransGen, Beijing, China). The CDS of *Bna.A07.WRKY70* was isolated from cDNA through PCR (Thermal Cycler Block, Thermo Fisher Scientific) using the high-fidelity thermostable DNA polymerase KOD-Plus-Neo (Toyobo Co., Ltd., Osaka, Japan). The PCR conditions were as follows: predenaturation at 94 ◦C for 2 min, followed by 35 cycles of 98 ◦C for 10 s, 55 ◦C for 30 s, and 68 ◦C for 1 min, and final extension at 68 ◦C for 7 min. Cloning primers are listed in Supplementary Table S1.

To construct the plasmid *35S:Bna.A07.WRKY70-GFP*, the CDS of *Bna.A07.WRKY70* without stop code was digested with the restriction endonucleases *XbaI* and *BamHI* and cloned into P1300-35S-green fluorescent protein (GFP) vector, which was driven by the CaMV35S (35S) promoter. Similarly, the digested PCR fragment of *Bna.A07.WRKY70* was also cloned into pGreen-35S-eGFP to produce a fusion of *GFP-Bna.A07.WRKY70* under the control of the 35S promoter. To obtain the construct of *pBna.A07.WRKY70:GUS*, the 2600 bp 5 regulatory region upstream of the ATG start codon, as the *Bna.A07.WRKY70* promoter region was amplified and cloned into pHY107 [43]. The CDS of *Bna.A07.WRKY70* was cloned into the pGBKT7 vector containing the GAL4 DNA binding domain to form a construct of *pBD-Bna.A07.WRKY70*. Eight single colonies of each plasmid were selected randomly and sequenced by Sangon Biotechnology (Shanghai, China). Primers used for plasmid construction are listed in Supplementary Table S1.

#### *4.4. Subcellular Localization of Bna.A07.WRKY70-GFP Protein*

The *35S:GFP-Bna.A07.WRKY70* construct was transformed into *Agrobacterium tumefaciens* strain GV3101 and transiently expressed in the leaves of transgenic tobacco (*Nicotiana benthamiana*) carrying a nuclear localization signal as previously described [44]. Images of fluorescent signals were detected through a confocal laser scanning microscope (Leica TCS SP8 SR, Wetzlar, Germany) 72 h after agroinfiltration of the tobacco plants.

#### *4.5. Transcriptional Activation Assays in Yeast*

The construct of *pBD-Bna.A07.WRKY70* and the negative control pGBKT7 vector were transformed separately into the yeast strain Y2HGold, including the *HIS3* and *ADE2* reporter genes. The transformed strains were cultured on synthetic dropout nutrient medium without tryptophan (SD/-Trp) plates and then were spotted on SDO (SD/-Trp) and TDO (SD/-Trp/-His/-Ade) plates by diluting to different concentrations. The transcription activation activity of each construct was observed according to the growth conditions of the corresponding yeast cells after incubating for 2–3 days in a 30 ◦C incubator.

#### *4.6. Generation of A. thaliana Transgenic Plants*

The construct of *pBna.A07.WRKY70:GUS* and *35S:Bna.A07.WRKY70-GFP* was transformed into *Agrobacterium tumefaciens* strain GV3101, which was subsequently used to transform the *A. thaliana* wild type and *wrky70* mutant plants using the floral dip method [45]. The transgenic lines of *pBna.A07.WRKY70:GUS* in wild type was selected on soil using Basta® and the transgenic lines of *35S:Bna.A07.WRKY70-GFP* in *wrky70* mutants were screened by hygromycin. All the transgenic plants were genotyped according to DNA and RNA analyses and selfed until T3 generation homozygous plants, which were generated and used for subsequent experiments.

#### *4.7. RNA Extraction and RT-qPCR Analysis*

The total RNA from various tissues of *B. napus* and leaves of *A. thaliana* were extracted using the SteadyPure Plant RNA Extraction Kit (Accurate Biology, Changsha, China). The quality of RNA was assessed using 1% agarose gel electrophoresis, and the concentration was determined by spectrometry (Nano Drop; Thermo Scientific, Wilmington, MA, USA) (Supplementary Table S3). RNA was reverse transcribed by EasyScript One-Step gDNA Removal and cDNA Synthesis SuperMix (TransGen, Beijing, China) according to the manufacturer's instructions, and conditions were 37 ◦C for 15 min; 85 ◦C for 5 s, followed by maintaining at 4 ◦C. Quantitative real-time PCR (RT-qPCR) was utilized to evaluate gene expression with SYBR Green Master Mix (Cofitt, Hongkong, China) using the QuantStudioTM 7 Flex Real-Time PCR System (Thermo Scientific), which were performed by three independent biological replicates with two technical replicates for each biological replicate. Reactions were performed in a total volume of 20 μL containing 100 nM of each primer and 2 μL of diluted cDNA (50 ng/μL) templates and amplified using the following cycling conditions: 95 ◦C for 2 min, 40 cycles of 95 ◦C for 15 s, 60 ◦C for 30 s, and 72 ◦C for 30 s. *AtACTIN7* (amplified product with 161 bp) and *BnACTIN7* (amplified product with 400 bp) were used as the internal control in *Arabidopsis* and rapeseed, respectively. For each reaction run, the specificity of the amplification was validated, and the threshold cycle (Ct) above the background was calculated using Bio-Rad iCycler (Bio-Rad, Hercules, CA, USA). The relative expression levels of the target genes were calculated using a modified double delta method [46]. Primers used for RT-qPCR analyses are listed in Supplementary Table S2.

#### *4.8. Phenotypic Observation of A. thaliana Leaves*

The seeds of the Col-0, the *wrky70* mutant, and three independent lines—*wrky70 35S:Bna.A07.WRKY70-GFP #4*, *#6*, *and #12*—were germinated on 1/2 MS agar medium for one week. Subsequently, the seedlings were transplanted into 8 × 8 cm pots. When the *A. thaliana* plants grew 35 days after germination (DAG), the phenotype of leaf senescence was observed and photographed by a camera (D7500, Nikon, Tokyo, Japan).

#### *4.9. Measurement of the Chlorophyll Content*

The fifth and sixth rosette leaf samples of Col-0, *wrky70* mutants, and *Bna.A07.WRKY70* transgenic plants from the fourth, fifth, sixth, and seventh weeks were separately collected and weighed and then placed in a 1.5 mL centrifuge tube with 1 mL extraction solution (80% acetone), soaked the leaves in the dark for 24 h until they faded [47]. To calculate the chlorophyll content of leaves, the 0.2 mL supernatant was absorbed into Costar 96 Flat Transparent plate, and the absorbance values at 663 nm and 645 nm were measured using a microplate reader (Infinite M200pro, Tecan, Mannedorf, Switzerland). Each experiment was represented by three biological replicates.

#### **5. Conclusions**

As an indispensable portion, leaf senescence spans the latter half of leaf development, which is essential to guarantee better production and survival of the next generation. This study suggested that *Bna.A07.WRKY70* may act as a negative regulator to share a conserved function with *AtWRKY70* in controlling leaf senescence when it is expressed in *A. thaliana*. Thus, *Bna.A07.WRKY70* can be utilized as a potential target to genetically manipulate leaf senescence and to create new stay-green materials to improve the rapeseed yield.

**Supplementary Materials:** The following supporting information can be downloaded at: https: //www.mdpi.com/article/10.3390/plants12020347/s1, Figure S1: Percent identity of full-length protein sequences of the WRKY70 protein from *A. thaliana* and *B. napus*; Table S1: Primers used for gene cloning and various constructs in the present study; Table S2: Primers used for RT-qPCR analysis in the present study; Table S3: The amount of RNA per sample in the *Arabidopsis thaliana* and *Brassica napus*.

**Author Contributions:** Y.G. and M.C. conceived and designed the experiments. T.L. conducted the experiments and analyzed the data. Y.L., C.W., D.Z., J.L. and M.H. conducted parts of the experiments. T.L. wrote the draft of the manuscript, and M.C. and Y.G. revised the manuscript. All authors have read and agreed to the published version of the manuscript.

**Funding:** This research was financially supported by a grant from the Yang Ling Seed Industry Innovation Center (Grant no. K3031122024) and the National Natural Science Foundation of China (Grant no. 31801393).

**Institutional Review Board Statement:** Not applicable.

**Informed Consent Statement:** Not applicable.

**Data Availability Statement:** All data included in this study are available upon reasonable request by contact with the corresponding author.

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

#### **References**


**Disclaimer/Publisher's Note:** The statements, opinions and data contained in all publications are solely those of the individual author(s) and contributor(s) and not of MDPI and/or the editor(s). MDPI and/or the editor(s) disclaim responsibility for any injury to people or property resulting from any ideas, methods, instructions or products referred to in the content.

## *Article* **Genome-Wide Characterization of the Sulfate Transporter Gene Family in Oilseed Crops:** *Camelina sativa* **and** *Brassica napus*

**Parviz Heidari 1,\*, Soosan Hasanzadeh 1, Sahar Faraji 2, Sezai Ercisli <sup>3</sup> and Freddy Mora-Poblete 4,\***


**Abstract:** Sulfate transporters (SULTRs) are responsible for the uptake of sulfate (SO4 <sup>2</sup>−) ions in the rhizosphere by roots and their distribution to plant organs. In this study, SULTR family members in the genomes of two oilseed crops (*Camelina sativa* and *Brassica napus*) were identified and characterized based on their sequence structures, duplication events, phylogenetic relationships, phosphorylation sites, and expression levels. In total, 36 and 45 putative *SULTR* genes were recognized in the genomes of *C. sativa* and *B. napus*, respectively. SULTR proteins were predicted to be basophilic proteins with low hydrophilicity in both studied species. According to the observed phylogenetic relationships, we divided the SULTRs into five groups, out of which the SULTR 3 group showed the highest variation. Additionally, several duplication events were observed between the *SULTRs*. The first duplication event occurred approximately five million years ago between three *SULTR 3.1* genes in *C. sativa.* Furthermore, two subunits were identified in the 3D structures of the SULTRs, which demonstrated that the active binding sites differed between *C. sativa* and *B. napus*. According to the available RNA-seq data, the *SULTRs* showed diverse expression levels in tissues and diverse responses to stimuli. *SULTR 3* was expressed in all tissues. *SULTR 3.1* was more upregulated in response to abiotic stresses in *C. sativa*, while *SULTR 3.3* and *SULTR 2.1* were upregulated in *B. napus*. Furthermore, *SULTR 3* and *SULTR 4.1* were upregulated in response to biotic stresses in *B. napus*. Additionally, the qPCR data showed that the *SULTRs* in *C. sativa* were involved in the plant's response to salinity. Based on the distribution of cis-regulatory elements in the promoter region, we speculated that *SULTRs* might be controlled by phytohormones, such as ABA and MeJA. Therefore, it seems likely that *SULTR* genes in *C. sativa* have been more heavily influenced by evolutionary processes and have acquired further diversity. The results reveal new insights of the structures and functions of SULTRs in oilseed crops. However, further analyses, related to functional studies, are needed to uncover the role of SULTRs in the plants' development and growth processes, as well as in their response to stimuli.

**Keywords:** bioinformatics; biotic stresses; regulatory mechanisms; protein structure; gene expression; evolutionary analysis

#### **1. Introduction**

Sulfur (S) is a macronutrient that is required for the biosynthesis of amino acids (such as cysteine (Cys) and methionine (Met)), vitamins, cofactors, and glutathione (GSH), as well as secondary metabolites; therefore, (S) is a vital element for plant growth, development, and stress response [1–3]. Root cells take up sulfate (SO4 <sup>2</sup>−) in the form of S through a proton codependent process. The uptake and assimilation of sulfate resources that are available in the environment produce essential sulfur (S) metabolites that are crucial for development and stress responses, which is critical for plants and microbes [4]. The soil sulfate content can be modified by various factors, such as the dissimilation of soil microbes,

**Citation:** Heidari, P.; Hasanzadeh, S.; Faraji, S.; Ercisli, S.; Mora-Poblete, F. Genome-Wide Characterization of the Sulfate Transporter Gene Family in Oilseed Crops: *Camelina sativa* and *Brassica napus*. *Plants* **2023**, *12*, 628. https://doi.org/10.3390/ plants12030628

Academic Editors: Mingxun Chen, Lixi Jiang and Yuan Guo

Received: 3 January 2023 Revised: 23 January 2023 Accepted: 28 January 2023 Published: 31 January 2023

**Copyright:** © 2023 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (https:// creativecommons.org/licenses/by/ 4.0/).

the weathering of S-containing minerals, human activities that modify the deposition of S into the ecosystem, and climate change [1]. Therefore, the available sulfate content in soil can also be altered because of the ability of plant root systems to absorb nutrient compounds according to their requirements and material accessibility. It has been reported that in comparison to other micronutrients, sulfate only has a gentle and limited effect on root structures [5]. To meet demands of S required for the S-containing metabolite synthesis, plant membrane transport systems and their related metabolic enzymes optimize sulfate uptake, acquisition, storage, and use [6]. The uptake and distribution of sulfate in plants are facilitated by networks of sulfate transporters (SULTRs), which are encoded by a multigene family [7]. The H+/SO4 <sup>2</sup><sup>−</sup> co-transporter SULTRs have been reported to contain 12 transmembrane domains, along with a carboxyl-terminal region, i.e., the so-called STAS (sulfate transporter/anti-sigma factor), which is suggested to play an important role in transporters' activity and their interactions with other proteins [1,8].

The involvement of SULTRs in the transportation of S within plants was first reported by Smith et al. [9]. SULTRs are characterized by 12 transmembrane domains (TMDs) and an anti-sigma factor antagonist (STAS) domain at the C-terminus, which is critical for sulfate transporter activity [10]. The genomes of higher plants, such as *Arabidopsis thaliana*, rice, wheat, sorghum, and apple, have been reported to have 12, 12, 11, 10, and 9 *SULTR* genes, respectively [11–14]. The SULTR family has been well characterized in *Arabidopsis*, and sulfate transporters can be divided into four main groups based on their sequence resemblance, function, and location. The first group includes AtSULTR 1.1, AtSULTR 1.2, and AtSULTR 1.3, which are all high-affinity S transporters [15]. AtSULTRs 1.1 and 1.2 are co-localized in root hairs and epidermal and cortical cells in roots, and they are both responsible for the uptake of sulfate from soil [16,17]. Nevertheless, despite their common function, AtSULTR 1.1 predominantly operates under the conditions of S deficiency, while AtSULTR 1.2 operates efficiently under the conditions of either sulfur abundance or sulfur deficiency [18]. AtSULTR 1.3 is localized in the phloem, and cooperates in the source-sink distribution of sulfate [19]. The second group consists of two low-affinity transporters, AtSULTR 2.1 and AtSULTR 2.2, which are responsible for the transportation of sulfate from root to shoot [20]. The third group comprises five members (AtSULTR 3.1-5) and is the largest group. However, the precise functions of these members have not been fully established. It has been reported that SULTR 3.1, which transports sulfate to chloroplasts, could have a role in helping plants to withstand abiotic stresses [21]. Additionally, SULTR 3.5 has been reported to co-express with SULTR 2.1 to enhance the uptake of sulfate and facilitate its transportation from root to shoot under conditions of S deficiency [22,23]. The fourth group of transporters, SULTR 4.1 and SULTR 4.2, have been demonstrated to be tonoplast localized transporters that release sulfate from vacuoles into the cytosol [24,25]. As well as the study on *A. thaliana*, many other studies have been conducted to functionally characterize SULTRs in crops. For instance, 14 putative *SULTR* genes have been identified in rapeseed (*Brassica napus*), among which only those from group 1 and group 4 were induced in response to S deficiency [26]. In another study, 28 putative *SULTR* genes were identified in the soybean (*Glycine max*) genome and *GmSULTR 1.2b* was confirmed to have important roles in sulfate uptake and improving plants' tolerance to sulfur deficiency [27]. In the potato (*Solanum tuberosum*) genome, 12 *SULTR* genes have been identified and the members of group 3 (StSULTR3s) were potentially proven to be involved in biotic/abiotic stress responses through MYB TFs, which play crucial roles in the modulation of StSULTR3s under these circumstances [28]. The maize (*Zea mays* L.) genome has been shown to include eight putative *SULTR* genes, all of which were induced by drought and heat stresses, except for *ZmSULTR 3.3* [29]. In addition, various studies have confirmed that *SULTRs* can be responsive to heavy metal exposure [30,31]. Despite the progress that has been made in the functional characterization of plant *SULTRs*, there are still more important crops that need to be investigated. *Camelina sativa* is an oilseed crop from the Brassicaceae family that has many qualities, including low inputs, great adaptation and resistance abilities, short life cycles, and easy genetic transformation, which have turned *C. sativa* into an ideal

model plant [32,33]. Moreover, *C. sativa* is becoming more important as a biofuel [34,35]. Although oilseed plants typically have very high S demands [36], a study on the response of *C. sativa* to various fertilizers showed that the seed yields and oil contents of *camelina* seeds were not affected by sulfur fertilization [37]. In order to develop S-efficient crops and crop varieties that are tolerant to S deficiency, it is necessary to identify and characterize SULTRs, especially in low-input crops, such as *C. sativa*. To the best of our knowledge, there are no available reports on the genome-wide analysis of *SULTR* genes in *C. sativa*, except for one study that reported the upregulation of *SULTR 3.4* in *C. sativa* under salinity stress [38]. In this study, resources were employed to distinguish the regulation roles of *SULTR* genes in various cellular processes, especially in response to stimuli. *B. napus* is another well-known oilseed plant containing appreciable amounts of erucic acid. In the present study, we focused on SULTR sequences in the *C. sativa* and *B. napus* genomes, and compared and discussed their adjustments and their possible engagement in protection mechanisms against unfavorable environmental stimuli. We also highlighted the potential properties of these genes that could help to facilitate sulfate uptake.

#### **2. Results**

#### *2.1. SULTR Properties in Camelina sativa and Brassica napus*

In the current study, 36 and 45 putative SULTR genes were recognized in the genomes of *C. sativa* and *B. napus*, respectively (Table S1). The SULTRs of the two oilseed crops were characterized and compared according to their coding DNA sequences (CDS) and protein lengths, exon numbers, isoelectric points (pIs), molecular weights (MWs), grand average of hydropathy (GRAVY) values, and instability indices (Table S1 and Table 1). Our results showed that the physicochemical properties of the SULTR proteins in the two studied plants were almost identical to each other. For instance, the MWs ranged from 29.07 to 91.99 kDa in *C. sativa*, and from 28.94 to 83.86 kDa in *B. napus*. Additionally, the pI values ranged from 7.41 to 9.93 in C. sativa, and from 7.11 to 10.71 in B. napus. Moreover, the GRAVY values varied from 0.271 to 0.624 in C. sativa, and from 0.108 to 0.621 in B. napus. Based on the instability indices, 83% and 73% of SULTR proteins were predicted to be stable proteins in *C. sativa* and *B. napus*, respectively. In addition, the exon numbers varied from 4 to 20 in *C. sativa* and from 4 to 19 in B. napus (Figure 1 and Table 1). Overall, the SULTR proteins were predicted to be basophilic proteins with low hydrophilicity.

**Table 1.** Summary of SULTRs properties in *Camelina sativa* and *Brassica napus*. Full details of SULTRs properties are provided in Table S1.


#### *2.2. Phylogenetic Analysis and Classification of the SULTR Gene Family*

In the present study, a phylogenetic tree of the SULTR proteins was created, including 45 SULTRs from *B. napus*, 36 SULTRs from *C. sativa*, 28 SULTRs from *Glycine max*, 12 SULTRs from *Oryza sativa*, and 12 SULTRs from *Arabidopsis thaliana* (Figure 1). The studied SULTRs were classified into five main groups: 16 SULTRs from SULTRs 4.1 and 4.2 were categorized into group 1; SULTRs 2.1 and 2.2 were clustered into group 2; 30 SULTRs from SULTRs 1.1, 1.2, and 1.3 were assigned to group 3; 28 proteins from SULTRs 3.3 and 3.4 were included in group 4; 34 SULTRs from SULTRs 3.1, 3.2, and 3.5 were located in group 5 (Figure 1). The SULTRs from monocot model plant (rice) were very different from the dicot samples. Moreover, the SULTRs from *C. sativa* and *B. napus* were evaluated and compared according

to the conserved motifs. Overall, 10 conserved motifs were recognized in the protein sequences of the SULTRs, among which motif 6 was not observed in the SULTRs in group 1 (Figure 2). Additionally, 10 conserved motifs were identified in SULTRs 2.1 and 2.2, except the SULTR 2.1 from *C. sativa* only showed eight conserved motifs. Furthermore, SULTRs 1.1, 1.2, and 1.3 and 3.1, 3.2, and 3.5 were very diverse, according to the patterns of their motif distributions (Figure 2). Motifs 7 and 2 were frequently observed in the SULTR proteins and showed potential as screening markers for members of this family.

**Figure 1.** The phylogenetic tree of the SULTRs from *Camelina sativa, Brassica napus, Arabidopsis thaliana, Glycine max*, and *Oryza sativa*. The exon numbers for the SULTR coding genes are shown in the blue bar (more details related to the gene structures are provided in Table S1).

#### *2.3. Evolutionary Processes in the MGT Genes of Citrullus lanatus and Cucumis sativus*

In this study, to investigate the duplication events that have occurred in the SULTR gene family in *C. sativa* and *B. napus*, the synonymous (Ks), non-synonymous (Ka), and Ka/Ks values of each duplicated gene pair were calculated (Figure 3 and Table S2). The Ks values of the *SULTRs* in *C. sativa* were frequently between 0.6 and 1.0 (Figure 3a), while the Ka/Ks values were frequently between 0.7 and 0.9 (Figure 3b). In contrast, the Ks and Ka/Ks values of the *SULTRs* in the *B. napus* genome differed from those in *C. sativa*, with the Ks values frequently being between 1.2 and 1.6 (Figure 3c) and the Ka/Ks values frequently ranging from 0.3 to 0.5 (Figure 3d). In *C. sativa*, the first duplication event was predicted to have occurred around five million years ago (MYA) between three *SULTR*

*3.1 genes*, including *Csa06g026100-Csa04g037720* and *Csa09g058940-Csa04g037720*, while the first duplication event in *B. napus* occurred approximately three MYA between two *SULTR 3.1 genes*, *BnaA03g41530* and *BnaA09g35200* (Table S2). Several synteny blocks were observed between the *SULTRs* from *C. sativa* and *B. napus* (Figure S2). Additionally, three *SULTR* 1.3 genes (*Csa17g029070, Csa14g027370*, and *Csa03g026040*), four *SULTR* 3 genes (*Csa13g054450*, *Csa08g050710, Csa02g005990*, and *Csa08g012360*), and a *SULTR 1.1* gene (*Csa08g034630*) from *C. sativa* showed fewer synteny relationships with *SULTRs* from *B. napus* (Figure 4).

**Figure 3.** The frequency of Ks and Ka/Ks values in the *SULTRs*: (**a**) the frequency of Ks values in the *SULTRs* of *C. sativa* (Cs); (**b**) the frequency of the Ka/Ks values in the *SULTRs* of *C. sativa* (Cs); (**c**) the frequency of Ks values in the *SULTRs* of *Brassica napus* (Bn); (**d**) the frequency of the Ka/Ks values in the *SULTRs* of *Brassica napus* (Bn). The full details of the duplicated *SULTRs* are provided in Table S2.

#### *2.4. Transmembrane Structures of SULTRs*

The SULTR proteins from different groups were compared based on their transmembrane structures *in C. sativa* and *B. napus* (Figure 5). In group 1, 12 transmembrane helices and 11 pores were identified in all SULTRs. However, the SULTRs in *B. napus* showed similar structures based on the positions of the transmembrane helices while the structures in *C. sativa* were diverse. Additionally, the number of transmembrane helices in the group 2 SULTRs ranged from 10 to 12 in *B. napus* and from 8 to 10 in *C. sativa*. Most of the SULTRs in *B. napus* showed 10 transmembrane helices with nine pores (except for BnaC07g18000D with seven transmembrane helices), while the number of transmembrane helices in *C. sativa* varied between 8 and 11. In group 4, the number of transmembrane helices in the SULTRs of *B. napus* ranged from 6 to 11, while the number of transmembrane helices in *C. sativa* ranged from 9 and 13. The SULTRs in group 5 were very diverse in terms of their transmembrane structures, in which between 4 and 14 transmembrane helices were observed.

**Figure 4.** The synteny relationships between the *SULTRs* from *Camelina sativa* and *Brassica napus*.

#### *2.5. 3D Structure Analysis of SULTRs*

Our analysis of the 3D structures revealed that the SULTRs in *C. sativa* and *B. napus* had two domains and that the active binding sites could be located in small or large subunits (Figure 6). These results showed that the SULTRs in *C. sativa* were different from those in *B. napus* (Figure 6). In the group 1 SULTRs, the valine (VAL), proline (PRO), phenylalanine (PHE), asparagine (ASN), lysine (LYS), glycine (GLY), and serine (SER) amino acids were frequently observed in the binding sites of SULTRs from *C. sativa*, while PHE, GLY, and leucine (LEU) were frequently observed in the binding sites of SULTRs from *B. napus* (Figure 6). In the group 2 SULTRs, PHE, GLY, and alanine (ALA) were more frequently observed in the binding sites of *C. sativa*, while PHE, SER, and isoleucine (ILE) were frequently observed in the binding sites of *B. napus*. Additionally, six amino acids, including SER, aspartic acid (ASP), LYS, ILE, ALA, and tyrosine (TYR), were more frequently observed in the binding sites of group 3 SULTRs in *C. sativa*, while PHE and threonine (THR) were frequently observed in the binding sites of *B. napus*. In the group 4 SULTRs, SER, GLY, histidine (HIS), and TYR were more commonly identified in the binding sites in *C. sativa*, while LEU, ILE, glutamate (GLU), and arginine (ARG) were frequently observed in the binding sites of *B. napus.* In the group 5 SULTRs, SER, PHE, ILE, ALA, VAL, LEU, and TYR were more frequently observed in the binding sites in *C. sativa*, while ALA, ILE, methionine (MET), VAL, and THR were frequently observed in the binding sites of *B. napus*.

**Figure 5.** The transmembrane structures of the SULTRs in *C. sativa* and *B. napus*. The grouping was based on the phylogenetic tree.

#### *2.6. SULTR Expression Analysis*

In this study, the expression patterns of *SULTRs* in *C. sativa* and *B. napus* were evaluated in different tissues and in response to stress (Figures 7 and 8). We found that two *SULTR* 3.5 genes (*Csa20g030350* and *Csa13g022560*) and two *SULTR* 1.2 genes (*Csa09g084780* and *Csa07g050670*) were expressed more in the roots of *C. sativa*, while three *SULTR* 3.1 genes (*Csa06g026100*, *Csa09g58940*, and *Csa04g0377720*) and three *SULTR 2.1* genes (*Csa13g011940*, *Csa08g054410*, and *Csa20g015450*) were highly expressed in stem tissues (Figure 7a). In the leaf tissues of *C. sativa*, three *SULTR* 3.3 genes (*Csa17g030170*, *Csa14g030330*, and *Csa03g026970*), two *SULTR* 2.2 genes (*Csa16g042230* and *Csa09g084770*), and a *SULTR 4.1* gene (*Csa20g018910*) were highly expressed (Figure 7a). In response to abiotic stresses, *SULTR 3.1* was induced in *C. sativa* (Figure 7b). For example, *Csa06g026100* and *Csa04g037720* were expressed more in response to cold and salt stresses, while *Csa09g058940* was expressed more in response to drought, cold, and cadmium stresses. In addition, *Csa20g018910* (which is a chloroplast *SULTR 4.1*) was expressed more under cold stress (Figure 7b). Additionally, the *SULTRs* of *B. napus* showed diverse expression levels in tissues and in response to abiotic and biotic stresses (Figure 8). We found that two *SULTR 2.1* genes (*BnaA02g00410D* and *BnaC02g00440D*), a *SULTR 3.4* gene (*BnaC01g35550D*), and a *SULTR 3.5* gene (*BnaC02g08870D*) were highly expressed in the root tissues of *B. napus*, while two *SULTR 3.2* genes (*BnaC09g00110D* and *BnaA09g01000D*), two *SULTR 3.1* genes (*BnaA03g41530D* and *BnaC07g32580D*), a *SULTR 3.3* gene (*BnaC05g18450D*), and a *SULTR 2.2* gene (*BnaC06g38470D*) were expressed in seeds (Figure 8a). In the stem tissues of *B. napus*, two *SULTR* 3 genes (*BnaA03g41530D* and *BnaC04g28500D*) were highly expressed, while three *SULTR* 3 genes (*BnaA09g32410D*, *BnaA07g10140D*, and *BnaC07g13290D*), a *SULTR 2.1* gene (*BnaC09g46440D*), and a *SULTR 4.1* gene (*BnaA03g04410D*) were expressed in leaf tissues (Figure 8a). Furthermore, two *SULTR 3.3* genes (*BnaC05g18450D* and *BnaA09g30120*) and two *SULTR 2.1* genes (*BnaA10g22050D* and *BnaC09g46440D*) were more upregulated in response to PEG, NaCl, and ABA treatment (Figure 8b). Interestingly, two *SULTR 2.1* genes (*BnaC06g38470D* and *BnaA07g33850D*) were differentially expressed in response to cold stress in *B. napus*. However, *BnaA07g10140D* (which is a *SULTR 3.3*) and *BnaC09g46440D* (which is a *SULTR 2.1*) were also upregulated under cold stress. In response to biotic stresses, two *SULTR 4.1* genes (*BnaC03g05940D* and *BnaA03g04410D*) were upregulated in response to the fungal pathogen *Leptosphaeria maculans*. In addition, a *SULTR 3.4* gene (*BnaC01g3550D*) and a *SULTR 3.3* gene (*BnaA07g10140D*) were expressed more in response to *Sclerotinia sclerotiorum* and *Bacillus thuringiensis* strain 4f5, respectively (Figure 8b).

**Figure 6.** The three-dimensional docking analysis of the SULTRs in *C. sativa* and *B. napus*. The ligand binding sites are highlighted in red and lists of the binding sites are provided next to the protein structures.

**Figure 7.** The expression levels of the *SULTRs* in *C. sativa*, based on the available RNA-seq data: (**a**) in different tissues; (**b**) in response to abiotic stresses.

**Figure 8.** The expression levels of the *SULTRs* in *B. napus*, based on the available RNA-seq data: (**a**) in different tissues; (**b**) in response to abiotic and biotic stresses.

#### *2.7. SULTR Phosphorylation Prediction*

The potential phosphorylation sites of the SULTRs in *C. sativa* and *B. napus* were predicted based on serine, threonine, and tyrosine amino acids (Figure 9). The potential phosphorylation sites in the SULTRs ranged from 3 (in Csa13g054450, which is a SULTR 3.2) to 21 (in Csa08g005450, which is a SULTR 4.1 from group 1), with an average of 10.28 sites per protein in *C. sativa* (Figure 9a). Interestingly, SULTR 4.1 showed a high potential for phosphorylation events in *C. sativa*. Additionally, the potential phosphorylation sites in the SULTRs in *B. napus* ranged from a site in BnaC07g18000D (which is a SULTR 1.1) to 23 sites in BnaA10g19810D (which is a SULTR 4.1), with an average of 9.71 sites per protein (Figure 9b). In addition, more phosphorylation sites were predicted in SULTR 4.1 in *B. napus*.

#### *2.8. Distribution of Cis-Regulatory Elements in Promoter Sites*

In this study, the distribution of cis-regulatory elements in the promoter sites of the *SULTRs* in *C. sativa* and *B. napus* was investigated (Figure 10, Figure S3, and S4). The *SULTRs* in *C. sativa* and *B. napus* were compared based on the cis-regulatory elements that were related to their responses to stress and hormones (Figure 10). The cis-regulatory elements associated with auxin, ABA, MeJA, GA, and SA responses were observed in the promoter regions of the *SULTRs*. The results revealed that the cis-regulatory elements of the ABA response were frequently distributed in the *SULTRs* from *C. sativa*, while the MeJA response elements were more commonly observed in *B. napus* (Figure 10). Additionally, the cis-regulatory elements related to biotic and cold stresses were more frequently observed in the *SULTRs* from *B. napus*, while those related to drought stress were more commonly observed in the promoter sites of the *SULTRs* from *C. sativa*.

**Figure 10.** A comparison between the *SULTRs* from *C. sativa* and *B. napus* based on the number of cis-regulatory elements related to hormone and stress responses in promoter sites. More details are provided in Figures S3 and S4.

#### *2.9. Expression Patterns of SULTRs in Camelina in Response to Salinity Stress*

To understand the potential roles of the *SULTR* genes in camelina plants, the expression levels of five selected genes were analyzed in response to salt stress (i.e., 200 mM of NaCl). The camelina *SULTR* genes illustrated different expression patterns under salinity (Figure 11). For instance, *Csa01g013600* (which is a SULTR 4.2) was downregulated after 6 h of salinity stress, while its expression was upregulated after 24 h. Moreover, *Csa16g042230* (which is a SULTR 2.2) and *Csa06g026100* (which is a SULTR 3.1) had similar expression patterns. Both genes were upregulated in response to salt stress and the maximum expression was observed after 72 h. In contrast, *Csa07g050670* (which is a SULTR 1.2) was not induced by salinity stress. The expression levels of one SULTR 3.4 gene (*Csa15g020720*) were significantly reduced after 24 h and 72 h of salt stress. Overall, these data showed that some SULTR family members were involved in the response to salt stress.

**Figure 11.** The expression levels of the *SULTRs* in *C. sativa* in response to salinity stress (i.e., 200 mM of NaCl) at three timepoints (6, 24, and 72 h after salt stress) and under control conditions (C, i.e., irrigation without NaCl), based on the qPCR data. Note: \* and \*\* indicate significant differences between the expression levels following the salt treatment and those under normal conditions (based on a Student's t-test) at *p* < 0.05 and *p* < 0.01, respectively.

#### **3. Discussion**

The uptake and distribution of sulfate in plants are facilitated by networks of sulfate transporters, which are encoded by a multigene family (SULTRs) [7]. Due to the important role of sulfate in plants, the SULTRs in several plant species have been studied. For instance, the genomes of higher plants, such as *Arabidopsis thaliana*, rice (12 *SULTRs*), wheat (11 *SULTRs*), sorghum (10 *SULTRs*), and apple (9 *SULTRs*), have been identified [11–14]. In this study, we identified and characterized 36 and 45 putative *SULTR* genes in the genomes of *C. sativa* and *B. napus*, respectively (Table S1). More members of this gene family could be associated with changes in ploidy levels and genome sizes in *C. sativa* and *B. napus*, as well as duplication events in evolutionary processes [35,39]. Our investigations revealed that the SULTR proteins in the two studied plants had the same ranges for their physicochemical properties, i.e., MWs, pIs, GRAVY values, and instability indices. In addition, the exon numbers ranged from 4 to 20 in *C. sativa* and from 4 to 19 in *B. napus*. The similarities in their gene structures could indicate that significant evolutionary events have occurred in the plant genomes [40,41]. Our findings also suggested that the exon/intron patterns could provide new insights into the evolutionary relationships among the members of the gene family and that they could have originated from a common ancestor. Moreover, it has been reported that the exon number can affect expression levels, and that genes with lower exon numbers can be expressed quickly in response to environmental stresses [42,43]. SULTRs have been divided into four main classes based on their locations and functions [4]. In this study, the different SULTR classes were further separated based on our phylogenetic analysis. The SULTR 4 genes were very distinct from the other classes, while the SULTR 3 members varied significantly (Figure 1). Differences have also been observed between the SULTRs in the model monocot plant, rice, and dicot plants, indicating that the diversity in the SULTR gene family has occurred after the divergence of monocots and dicots [44,45]. According to our results for the conserved motifs in the SULTRs, some conserved sites were common between SULTR groups, which could be used to distinguish between specific groups.

According to our phylogenetic results, the camellia SULTRs were similar to the SULTRs of *B. napus*, although their evolutionary trends were different. Based on the Ka/Ks indices, the first duplication events in the SULTR genes in *C. sativa* occurred about five million years ago, while those in *B. napus* occurred three million years ago. Furthermore, it seemed that other members of the SULTR gene family originated from SULTR 3. Additionally, the Ka/Ks values revealed that the duplicated SULTRs in *B. napus* occurred under purifying (negative) selection, while both adaptive (positive) selection and purifying selection were observed in the *SULTRs* of *C. sativa* [46]. This suggested that the duplicated genes with conserved functions, pseudogenization, or both were possibly produced via purifying selection [47]. Interestingly, the results of our comparative synteny analysis revealed that several *SULTRs* from *C. sativa*, including three *SULTR 1.3* genes (*Csa17g029070*, *Csa14g027370*, and *Csa03g026040*), four *SULTR 3* genes (*Csa13g054450*, *Csa08g050710*, *Csa02g005990*, and *Csa08g012360*), and a *SULTR 1.1* gene (*Csa08g034630*), had fewer synteny relationships with the *SULTRs* from *B. napus* (Figure 4). It seemed that these genes could have been specifically developed during the evolution of the camellia, although more research is needed to determine their functions.

SULTRs can be classified into four groups based on their sequence structures, locations, and functions [48]. For instance, the genes in group 1 and group 2 are expressed more in root cells and vacuolar tissues, respectively [48,49]. In this study, the *SULTRs* in *C. sativa* and *B. napus* showed diverse expression levels in different tissues and in response to stresses. In the roots of *C. sativa*, two *SULTR 1.2* genes and two *SULTR 3.5* genes were expressed more, while two *SULTR 2.1* genes (*SULTR 3.4*, and *SULTR 3.5)* were highly expressed in the root tissues of *B. napus*. In the shoot tissues, *SULTRs 2, 3*, and 4 were expressed more. Interestingly, *SULTR* 3 showed a diverse range of functions and was expressed in all tissues, indicating that the members of this class were not specific to a tissue or organ. In addition, the members of *SULTR 3* varied greatly in terms of their transmembrane structure. Moreover, different expression patterns were observed between the members of the SULTR

gene family in *B. napus* and camellia in response to stimuli. The *SULTR 3.1* genes were expressed more in response to abiotic stresses in *C. sativa*, while the *SULTR 3.3* and *SULTR 2.1* genes were more upregulated in *B. napus*. Several members of SULTR 3 play multiple roles and interact with abscisic acid (ABA) metabolism [21–23]. In the present study, SULTR 3 and SULTR 4.1 were upregulated in response to biotic stresses in *B. napus*, including bacterial and fungal pathogens. Additionally, the cis-regulatory elements related to ABA and MeJA responses were frequently observed in the promoter sites of the SULTRs. We concluded that the SULTRs could be controlled by phytohormones, especially the hormones related to stress, such as ABA and MeJA. These interactions could effectively induce the expression of the members of this gene family in response to stress. It can also be stated that the expression levels of different *SULTRs* could be correlated with hormone and stress response elements observed in the promoter regions. Additionally, the real-time PCR data revealed that the *SULTRs* in *C. sativa* had diverse expression patterns and were involved in the response to salt stress. This indicates that SULTRs could possibly interact with some transcription factors, such as MYB, and be indirectly involved in responses to abiotic stresses [28]. The prediction of the 3D structures revealed two subunits in the SULTRs and that the active binding sites differed between the subgroups (Figure 6). PHE, ALA, ILE, and VAL were identified as key amino acids in the binding sites, playing critical roles in the function and regulation of the SULTRs. Post-translational phosphorylation modifications can affect the function and possible interaction of proteins [50,51]. The prediction of the phosphorylation sites in the SULTRs revealed that the SULTR 4.1 genes had a high potential for influencing post-translation modifications, such as phosphorylation. The SULTR 4.1 and SULTR 4.2 genes have been reported to be tonoplast transporters, which allow sulfate to leave vacuoles to reach cytosol [24,25]. It seems that phosphorylation modifications play key roles in the activity of these transporters.

#### **4. Materials and Methods**

#### *4.1. Identification of SULTR Genes in C. sativa and B. napus*

To identify all sequences related to the SULTR family, the amino acid sequences of two conserved domains, including Sulfate\_transp (PF00916) and STAS (PF01740), were used as queries in a BLASTP search of Ensembl Plants (https://plants.ensembl.org/index.html, accessed: 20 September 2022) in the protein databases of *C. sativa* and *B. napus*. Furthermore, orthologue genes were identified by following the same procedure for *Arabidopsis thaliana*, *Oryza sativa*, and *Glycine max*. All collected sequences were checked using the NCBI Conserved Domain Database (CDD) [52] and the Pfam database [53] to confirm the presence of domains related to the SULTRs [54]. The physiochemical properties, including molecular weight (MW), instability index, isoelectric point (pI), and GRAVY value, of the SULTRs were predicted using the ProtParam tool [55]. The TMHMM version 2.0 server was used to predict the transmembrane structures of the SULTRs in *C. sativa* and *B. napus* [56].

#### *4.2. Phylogenetic and Conserved Motif Analyses*

The amino acid sequences of all the identified *SULTRs* from five plant species, i.e., *C. sativa*, *B. napus*, *A. thaliana*, *O. sativa*, and *G. max*, were aligned using the online tool Clustal-Omega [57]. The entire phylogenetic relationships were constructed using the maximum likelihood (ML) method with 1000 bootstrap replicates using the IQ-TREE server [58]. Finally, a phylogenetic tree was created using the interactive tree of life tool (iTOL version 5) [59]. The conserved protein motifs in the SULTRs in *C. sativa* and *B. napus* were identified using the Multiple Expectation Maximization for Motif Elicitation program (MEME version 5.0.5) [60].

#### *4.3. Promoter Analysis*

In this study, 1500 bp upstream of the start codon in the *SULTRs* was selected as the promoter site, and these regions in *C. sativa* and *B. napus* were downloaded from Ensembl Plants. The sequence of each promoter site was screened to identify the conserved cisregulatory elements using the PlantCARE tool [61]. Then, the cis-regulatory elements were classified based on their functions.

#### *4.4. Ka/Ks Ratio and Duplication Analysis*

In the present study, pairs of *SULTR* genes from each species (*C. sativa* and *B. napus*) with similarities of more than 85% were considered to be duplicated genes [62]. Additionally, the synonymous (Ks) and non-synonymous (Ka) indices were calculated for all gene pairs using the MEGAX software [63]. The time of divergence of the duplicated *SULTR* genes was estimated using the following equation: T = (Ks/2λ) × <sup>10</sup>−6. (<sup>λ</sup> = 6.5 × <sup>10</sup><sup>−</sup>9) [64]. In addition, the synteny relationships between the SULTRs in each species, and between the orthologous genes of *C. sativa* and *B. napus*, were drawn using the Circos tool [65].

#### *4.5. Gene Expression Analysis*

In this study, the available RNA-seq data for *C. sativa* and *B. napus* were screened to extract the expression levels of the *SULTR* genes. In total, four RNA-seq datasets for *C. sativa*, including SRR935368 (root tissue), SRR935362 (leaf tissue), SRR935365 (stem tissue), and SRR935369 (flower tissue) were retrieved from the NCBI gene bank and analyzed. To extract the expression patterns of the *SULTRs* in response to stresses, the RNA-seq datasets related to salt stress (SRR935382), drought stress (SRR935380), cadmium stress (SRR935383), cold stress (SRR935372), and control conditions (SRR935385) were used. For the raw data analysis, we used FastQC software (version 0.11.6) (http://www.bioinformatics.babraham. ac.uk/projects/fastqc/, accessed: 20 September 2022) to check the quality of the data and HISAT [66] to map the sequences. The FPKM (fragments per kilobase of exon model per million mapped reads) metric was used to evaluate the transcription levels of each *SULTR* gene in *C. sativa*. To illustrate the expression levels of the *SULTRs* in *B. napus*, we utilized RNA-seq data for the rapeseed cultivar ZhongShuang11 (ZS11), which were related to 18 tissues and responses to biotic and abiotic stresses, from the Brassica Expression Database [67]. The expression patterns of the target genes were extracted based on their FPKM values. Finally, heatmaps were constructed using the log2 transformed method in TBtools software (version 0.665) [68].

#### *4.6. Prediction of 3D Structures, Modeling, Binding Sites, and Phosphorylation*

In this study, five proteins from each species (*C. sativa* and *B. napus*) were selected, based on the phylogenetic tree. Additionally, the three-dimensional structures of 10 SULTRs were predicted using the Phyre2 server [69]. In the next step, the predicted structures were checked using a Ramachandran plot analysis [70]. The binding sites of each model were highlighted on the predicted structures. The NetPhos 3.1 server [71], with a potential value of more than 0.90, was used to predict the phosphorylation sites of the SULTRs in *C. sativa* and *B. napus*.

#### *4.7. Expression Patterns of SULTR Genes in C. sativa under Salinity Stress*

Sterilized camelina seeds were planted at a depth of 2 cm in pots containing peat moss and were kept under the conditions of 16 h of light and a temperature of 25 ◦C with irrigation every three days. Then, the five-week-old seedlings were treated with salt (200 mM of NaCl) via irrigation, which was repeated after 24 h. After the salt treatment, leaves were collected at different time points (after 6, 24, and 72 h). The total RNA samples were extracted using an RNX kit (Sinaclon, Iran) and the cDNA was synthesized using a reverse transcriptase kit (Roche, Germany), according to manufacturer protocols. In the present study, five members of the SULTR family were selected for real-time PCR analysis. The genes were selected based on the phylogeny analysis. In addition, *actin-2* (Csa15g026420) was used as a reference gene to normalize the expression data. Specific primers were designed using the Primer3 online software (version 4.1.0) [72], based on the coding sequences of the selected SULTR genes (Table S3). The expression patterns of

the SULTR genes were analyzed using a Maxima SYBR Green/ROX qPCR Master Mix kit (Thermo Fisher, France) and the ABI Step One, according to manufacturer protocols. The expression levels of each SULTR gene were calculated using the delta Ct method [73], using three biological replicates.

#### **5. Conclusions**

In this study, we identified and characterized 36 and 45 putative *SULTR* genes in two important oilseed crops, *Camelina sativa* and *Brassica napus*, respectively. We found that the first duplication event occurred in the *SULTR* genes of *C. sativa* and that members of this family showed diverse structures and functions. Additionally, several *SULTR* genes in *C. sativa* were uniquely developed under evolutionary processes. SULTR 3 was identified as the class of sulfate transporter family genes with the highest diversity. Overall, our results revealed new insights into the structures and functions of SULTRs in oilseed crops. However, further functional studies are needed to evaluate the roles of SULTRs in development and growth processes, as well as in responses to stimuli. Also, investigation of upstream key proteins/enzymes that affect the activity of SULTRs, can reveal the pathways linked to SULTR.

**Supplementary Materials:** The following are available online at https://www.mdpi.com/article/ 10.3390/plants12030628/s1, Table S1: A list of the identified SULTRs and their characteristics in Camelina sativa and Brassica napus, Table S2: The predicted Ka/Ks values in the duplicated gene pairs from the sulfate transporter family in the Camelina sativa and Brassica napus genomes, Table S3: A list of the primers for the camelina SULTR genes that were used in our real-time PCR, Figure S1: The logos of 10 conserved motifs in the sulfate transporter family proteins in Camelina sativa and Brassica napus, Figure S2: A synteny analysis of the SULTR genes in the (a) Camelina sativa and (b) Brassica napus genomes, Figure S3: The distribution of cis-regulatory elements in the SULTR promoter site of Camelina sativa, Figure S4: The distribution of cis-regulatory elements in the SULTR promoter site of Brassica napus.

**Author Contributions:** Conceptualization, P.H., S.H. and S.F.; methodology, P.H. and S.F.; software, P.H., S.E. and F.M.-P.; writing—original draft preparation, P.H. and F.M.-P.; writing—review and editing, P.H. and F.M.-P. All authors have read and agreed to the published version of the manuscript.

**Funding:** This research received no external funding.

**Institutional Review Board Statement:** Not applicable.

**Informed Consent Statement:** Not applicable.

**Data Availability Statement:** Not applicable.

**Acknowledgments:** F.M.-P. acknowledges the support from ANID FONDECYT grant No. 1201973.

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

#### **References**


**Disclaimer/Publisher's Note:** The statements, opinions and data contained in all publications are solely those of the individual author(s) and contributor(s) and not of MDPI and/or the editor(s). MDPI and/or the editor(s) disclaim responsibility for any injury to people or property resulting from any ideas, methods, instructions or products referred to in the content.
