Genome-Wide Analysis of Genes Involved in the GA Signal Transduction Pathway in ‘duli’ Pear (Pyrus betulifolia Bunge)

Gibberellic acid (GA) is an important phytohormone that regulates every aspect of plant growth and development. While elements involved in GA signaling have been identified and, hence, their functions have been well studied in model plants, such as Arabidopsis and rice, very little is known in pear. We, therefore, analyzed the genes related to GA signaling from the recently sequenced genome of the wildtype ‘duli’ pear (Pyrus betulifolia Bunge), a widely used rootstock for grafting in pear cultivation in China due to its vigorous growth and resistance to abiotic and biotic stress. In total, 15 genes were identified, including five GA receptors PbGID1s (GA-INSENSTIVE DWARF 1), six GA negative regulators, PbDELLAs, and four GA positive regulators, PbSLYs. Exogenous application of GA could promote the expression of PbGID1s but inhibit that of PbDELLAs and PbSLYs in tissue culture ‘duli’ pear seedlings. The expression profiles of these genes in field-grown trees under normal growth conditions, as well as in tissue-cultured seedlings treated with auxin (IAA), GA, paclobutrazol (PAC), abscisic acid (ABA), and sodium chloride (NaCl), were also studied, providing further evidence of the involvement of these genes in GA signaling in ‘duli’ pear plants. The preliminary results obtained in this report lay a good foundation for future research into GA signaling pathways in pear. Importantly, the identification and preliminary functional verification of these genes could guide molecular breeding in order to obtain the highly desired dwarf pear rootstocks for high-density plantation to aid easy orchard management and high yielding of pear fruits.


Introduction
Gibberellins or gibberellic acids (GAs), tetracyclic diterpenoid phytohormones, are widely distributed in higher plants [1,2]. They play an important role in seed germination, hypocotyl elongation, flowering time transition, seed and fruit development, and abiotic and biotic stress response [3][4][5][6]. Since the discoveries of the GA-insensitive dwarf mutants of wheat and rice, one of the most decisive players in 'Green Revolution' in the late 1960s which led to a massive increase in grain production, saving millions of lives in developing countries, much effort has been put in place to unravel the molecular mechanism of the GA signal transduction pathway in plants. As a result, the key genes involved in the GA signaling pathway have been isolated, including those encoding for GID1s, DELLAs, and SLYs from various plants [4,7,8].
GID1, a cytosolic gibberellin receptor, was first discovered in rice (OsGID1) for its ability to sense and bind active GAs, transmit GA signals, and induce downstream reactions [9]. Subsequently, many other GID1s were also successfully isolated and studied in other plant species [3,10], the first of which were the three OsGID1 homologous genes from Arabidopsis, AtGID1a, AtGID1b, and AtGID1c. All three AtGID1s could strongly bind active GA and physically interact with DELLA proteins in the presence of GA [11]. They are also functionally redundant because the single mutant has no obvious GA-insensitive phenotype, whilst GA-GID1-DELLA complex in turn allows the C-terminus of DELLA to interact with the SCF SLY1/GID2/SNE complex. This leads to the ubiquitination and subsequent degradation of DELLA by the proteasome complex, resulting in the release of GA from the GA-GID1-DELLA complex. GA can then enter the nucleus to activate the expression of downstream genes, thereby promoting GA-related growth and development [33].
While GA signaling in the regulation of architecture, fertility, and stress is well characterized in Arabidopsis, rice, and other crop plant species, this knowledge in woody plants, such as pear tree, is very limited. Given the important role of DELLAs and other GA signaling components in the regulation of plant height and the lack of dwarf pear rootstocks, we aimed to first identify GID1, DELLA, and SLY homologous genes from the wildtype 'duli' pear genome. Subsequently, the evolution and expression profiles of these genes in various tissues of 5 year old 'duli' pear trees under normal growth conditions, as well as in tissue-cultured seedlings treated with IAA, GA, PAC, ABA, and NaCl, were assessed in order to dissect their biological functions. The preliminary data obtained in this study could provide valuable information for further research of these genes, which could guide molecular breeding to create the much-desired dwarf 'duli' pear rootstocks in future.
the GA-GID1-DELLA complex in turn allows the C-terminus of DELLA to interact with the SCF SLY1/GID2/SNE complex. This leads to the ubiquitination and subsequent degradation of DELLA by the proteasome complex, resulting in the release of GA from the GA-GID1-DELLA complex. GA can then enter the nucleus to activate the expression of downstream genes, thereby promoting GA-related growth and development [33].
While GA signaling in the regulation of architecture, fertility, and stress is well characterized in Arabidopsis, rice, and other crop plant species, this knowledge in woody plants, such as pear tree, is very limited. Given the important role of DELLAs and other GA signaling components in the regulation of plant height and the lack of dwarf pear rootstocks, we aimed to first identify GID1, DELLA, and SLY homologous genes from the wildtype 'duli' pear genome. Subsequently, the evolution and expression profiles of these genes in various tissues of 5 year old 'duli' pear trees under normal growth conditions, as well as in tissue-cultured seedlings treated with IAA, GA, PAC, ABA, and NaCl, were assessed in order to dissect their biological functions. The preliminary data obtained in this study could provide valuable information for further research of these genes, which could guide molecular breeding to create the much-desired dwarf 'duli' pear rootstocks in future.

PbGID1s, PbDELLAs, and PbSLYs Contain the Typical Conserved Domains and Motifs
As shown in Figure 3a, all PbGID1s contained two conserved domains, the co-esterase (pfam00135) and α/β hydrolase-3 (pfam07859) belonging to the α/β hydrolase family. These PbGID1s also contained 10 conserved motifs (Supplementary Materials Figure S1a  Analysis of the six PbDELLAs showed that they all contained the characteristic conserved DELLA (pfam12041) domain in their N-termini and the GRAS (pfam03514) domain in the C-termini (Figure 3b). Note that the DELLA (motif 5) showed some variations, where the amino acids EL in PbRGLa and PbRGLb were replaced by GC and GY ((DGCLA and DGYLA instead of DELLA), respectively. It is also noteworthy that the sequence of PbGAI1b was different from the other five PbDELLAs in that it contained a putative transmembrane domain between amino acids 2 and 24, implying that it may have a different function (Figure 3b). With the exception of PbGAI2b that did not have motif 9 SAW close to the C-terminus, a typical repressor motif of the GRAS family proteins, all other PbDELLAs had the 10 conserved motifs (Supplementary Materials Figure  S1b).
Conserved domains and motifs of PbSLYs were also analyzed, and the results are shown in Figure 3c. Only one conserved domain, the F-box (pfam00646) domain, was present. Conserved motifs 1, 2, 3, and 6 were found in the sequences of all four PbSLYs, whilst motif 4 and motif 5 were only found in PbSLY2s and PbSLY1s, respectively, suggesting that these two types of PbSLYs may have different roles in GA signaling. Analysis of the six PbDELLAs showed that they all contained the characteristic conserved DELLA (pfam12041) domain in their N-termini and the GRAS (pfam03514) domain in the C-termini ( Figure 3b). Note that the DELLA (motif 5) showed some variations, where the amino acids EL in PbRGLa and PbRGLb were replaced by GC and GY ((DGCLA and DGYLA instead of DELLA), respectively. It is also noteworthy that the sequence of PbGAI1b was different from the other five PbDELLAs in that it contained a putative transmembrane domain between amino acids 2 and 24, implying that it may have a different function (Figure 3b). With the exception of PbGAI2b that did not have motif 9 SAW close to the C-terminus, a typical repressor motif of the GRAS family proteins, all other PbDELLAs had the 10 conserved motifs (Supplementary Materials Figure S1b).

Synteny and Duplication Analysis of PbGID1s, PbDELLAs, and PbSLYs of 'duli' Pear
Conserved domains and motifs of PbSLYs were also analyzed, and the results are shown in Figure 3c. Only one conserved domain, the F-box (pfam00646) domain, was present. Conserved motifs 1, 2, 3, and 6 were found in the sequences of all four PbSLYs, whilst motif 4 and motif 5 were only found in PbSLY2s and PbSLY1s, respectively, suggesting that these two types of PbSLYs may have different roles in GA signaling.

Synteny and Duplication Analysis of PbGID1s, PbDELLAs, and PbSLYs of 'duli' Pear
Understanding the gene duplication events occurring in the genome and the synteny between different genomes can help understand gene evolution and function. Gene duplication, including tandem duplication, segmental duplication, and whole-genome duplication (WGD), is the major driving force in plant evolution. A multigene family is the result of region-specific duplication of the genome or WGD [39]. To evaluate the gene duplication events for PbGID1s, PbDELLAs, and PbSLYs, synteny and selective pressure analyses were carried out, and the results are presented in Figure 4a. PbGID1c-1-1 had a syntenic relationship with PbGID1c-1-2 and PbGID1c-2, as did PbGID1b-1 with PbGID1b-2. For PbDELLAs, PbGAI1a and PbGAI1b had a syntenic relationship, as did PbRGLa and PbRGLb. For PbSLYs, a syntenic relationship was found between PbSLY1-1 and PbSLY1-2, as well as between PbSLY2-1 and PbSLY2-2. All these gene pairs were the result of WGD or segmental duplications.
Understanding the gene duplication events occurring in the genome and the synteny between different genomes can help understand gene evolution and function. Gene duplication, including tandem duplication, segmental duplication, and whole-genome duplication (WGD), is the major driving force in plant evolution. A multigene family is the result of region-specific duplication of the genome or WGD [39]. To evaluate the gene duplication events for PbGID1s, PbDELLAs, and PbSLYs, synteny and selective pressure analyses were carried out, and the results are presented in Figure 4a. PbGID1c-1-1 had a syntenic relationship with PbGID1c-1-2 and PbGID1c-2, as did PbGID1b-1 with PbGID1b-2. For PbDELLAs, PbGAI1a and PbGAI1b had a syntenic relationship, as did PbRGLa and PbRGLb. For PbSLYs, a syntenic relationship was found between PbSLY1-1 and PbSLY1-2, as well as between PbSLY2-1 and PbSLY2-2. All these gene pairs were the result of WGD or segmental duplications. Next, the syntenic relationships between the genes of each family of GID1s, DELLAs, and SLYs were compared between 'duli' pear and Arabidopsis. As shown in Figure 4b, PbGID1c-1-1, PbGID1c-1-2, and PbGID1c-2 had synteny with Arabidopsis AtGID1c, of which PbGID1c-1-2 and PbGID1c-2 also had synteny with AtGID1a. Both PbGID1b-1 and PbGID1b-2 were syntenic with AtGID1b. For the PbDELLAs, PbGAI1a and PbGAI1b had a syntenic relationship with AtGAI and AtRGA, respectively, as did PbGAI2a with AtRGL1. PbSLY2-2 had a syntenic relationship with AtSLY2/AtSNE. Next, the syntenic relationships between the genes of each family of GID1s, DELLAs, and SLYs were compared between 'duli' pear and Arabidopsis. As shown in Figure 4b, PbGID1c-1-1, PbGID1c-1-2, and PbGID1c-2 had synteny with Arabidopsis AtGID1c, of which PbGID1c-1-2 and PbGID1c-2 also had synteny with AtGID1a. Both PbGID1b-1 and PbGID1b-2 were syntenic with AtGID1b. For the PbDELLAs, PbGAI1a and PbGAI1b had a syntenic relationship with AtGAI and AtRGA, respectively, as did PbGAI2a with AtRGL1. PbSLY2-2 had a syntenic relationship with AtSLY2/AtSNE. Therefore, on the basis of the syntenic relationship between these genes in 'duli' pear and Arabidopsis, they are divided into two groups. The first group, including PbGAI2b, PbRGLa, and PbRGLb of the PbDELLA family and PbSLY1-1, 1-2, and 2-1 of the PbSLY family, was due to 'dispersed' duplication, indicating a separation by other sequences where the genes may have arisen from transposition, such as 'replicative transposition', 'non-replicative transposition', or 'conservative transposition'. The second group originated from WGD or segmental duplication, and it included all PbGID1s, PbGAI1a, PbGAI1b, and PbGAI2a of the PbDELLAs, as well as PbSLY2-2.
Ka/Ks represents the ratio between the nonsynonymous substitution rate (Ka) and synonymous substitution rate (Ks) of protein-coding genes. While synonymous codons produce the same amino acids (synonymous changes), nonsynonymous changes result in different amino acids to be translated. Therefore, the Ka/Ks ratio is used to determine whether there is a selective pressure acting on the protein-coding gene; we found that the values of Ka/Ks for the GA signal transduction-related genes in both 'duli' pear and Arabidopsis were all less than 1 (Supplementary Materials Table S1). This clearly demonstrates that these genes of 'duli' pear have been subjected to purification selection during the evolution process; hence, the individual gene pairs identified above may have similar functions between Arabidopsis and 'duli' pear.

Three Types of cis-Acting Elements Were Present in the Promotor Regions of PbGID1s, PbDELLAs, and PbSLYs
Spatial and temporal expression of a gene is very important for its proper biological function. This is determined by the cis-acting elements in its promoter region. As shown in Figure 5 and Supplementary Materials Figure S2, three types of cis-acing elements were present in the promoter regions of PbGID1s: growth and development regulatoryresponsive, hormone-responsive, and stress responsive elements, named here as class I, II, and III, respectively. While the types and numbers of cis-elements in each type were very similar, the light-responsive element GATA-motif, the 6K protein-binding site Unnamed-1, and the palisade tissue mesophyll cell differentiation element HD-zip1 in class I, as well as the TC-rich repeats in class III, were unique to PbGID1c-1-1 and PbGID1c-1-2. PbGID1c-2 was also very different in that the zeatin metabolism regulatory element O2-site and the auxin response element TGA-element in class II were not present. Interestingly, the promoter of PbGID1c-2 contained unique cis-acting elements, the dehydration, low temperature, and salt stress response elements DRE1 and WRE3, and the wound and pathogen response element W box. For the promoter sequences of PbGID1b-1 and PbGID1b-2, there were also some unique elements found, such as (1) the regulatory element A-box of class I present in both, (2) the F-box and the rhythm response element circadian present in class I, the GA response element GARE-motif present in class II, and the drought-induced response element MBS in class III unique to PbGID1b-2, and (3) the light-responsive element AT1-motif, AE-box, and ACE in class I and the salicylic acid-induced cis-acting element TCA-element in class II unique to PbGID1b-1.
The common cis-acting elements found in the promoter regions of all six PbDELLAs were the G-Box, ABRE, ARE, and W box although the numbers were slightly different ( Figure 5b). However, there were some unique elements in each of the proPbDELLAs. For example, the AT1-motif of class I and TGA-element of class II were found in the promoter region of PbGAI1b, while the p-Box was found in PbGAI2a. The GCN4-motif and TATC-box were unique to proPbGAI2b, the 3-AF1 binding site and the cis-acting element DRE were unique to proPbRGLa, the TCCC-motif and O2-site were unique to PbRGLb, and chs-CMA2a and WRE3 were unique to proPbRGLa and proPbRGLb. The common cis-acting elements found in the promoter regions of all six PbDELLAs were the G-Box, ABRE, ARE, and W box although the numbers were slightly different ( Figure 5b). However, there were some unique elements in each of the proPbDELLAs. For example, the AT1-motif of class I and TGA-element of class II were found in the promoter region of PbGAI1b, while the p-Box was found in PbGAI2a. The GCN4-motif and TATC-box were unique to proPbGAI2b, the 3-AF1 binding site and the cis-acting element DRE were unique to proPbRGLa, the TCCC-motif and O2-site were unique to PbRGLb, and chs-CMA2a and WRE3 were unique to proPbRGLa and proPbRGLb.
The cis-acting elements found in the promoter regions of all PbSLYs were Box 4, GT1-motif, and as-1 of class I and ABRE and TGACG-motif of class II, although their numbers were slightly different (Figure 5c). The TCT-motif and CAT-box in class I, GARE-motif in class II, and MBS in class III were shared by PbSLY1-1 and PbSLY1-2, while chs-CMA2a and RY-element in class I were present on both PbSLY2-1 and PbSLY2-2. Unique elements, such as Gap-box, I-box, AT1-motif, F-box and transcription factor-binding site AP-1 in class I, as well as LTR and WRE3 in class III, were only found in the promoter region of PbSLY1-1, while chs-CMA1a and GCN4-motif of class I and TGA-element of class II were only in that of PbSLY2-1. Furthermore, the TATC-box, p-box, and TCA-element in class II and the WUN-motif in class III were unique to PbSLY1-2, while the HD-zip1 of class I and TC-rich repeats of class III were unique to PbSLY2-2. The cis-acting elements found in the promoter regions of all PbSLYs were Box 4, GT1motif, and as-1 of class I and ABRE and TGACG-motif of class II, although their numbers were slightly different (Figure 5c). The TCT-motif and CAT-box in class I, GARE-motif in class II, and MBS in class III were shared by PbSLY1-1 and PbSLY1-2, while chs-CMA2a and RY-element in class I were present on both PbSLY2-1 and PbSLY2-2. Unique elements, such as Gap-box, I-box, AT1-motif, F-box and transcription factor-binding site AP-1 in class I, as well as LTR and WRE3 in class III, were only found in the promoter region of PbSLY1-1, while chs-CMA1a and GCN4-motif of class I and TGA-element of class II were only in that of PbSLY2-1. Furthermore, the TATC-box, p-box, and TCA-element in class II and the WUN-motif in class III were unique to PbSLY1-2, while the HD-zip1 of class I and TC-rich repeats of class III were unique to PbSLY2-2.
Taken together, it was found that the common cis-acting elements were present in the promoter regions of the genes in the same family. However, there were also unique ones found in some genes within the same family. This suggests that, while genes in the same family play similar roles, each gene also has its unique function in growth, development, hormone response, and stress response.

PbGID1s, PbDELLAs, and PbSLYs Were Ubiquitously Expressed with Tissue-Specific High-Level Expression for Some Genes in Each Family
In order to understand the expression profiles of PbGID1s, PbDELLAs, and PbSLYs, real-time RT-PCR was carried out using total RNA isolated from roots, shoots, leaves, flowers, and young fruits of 5 year old trees. The results showed that the transcripts of all candidate genes were present in all the tissues tested; however, their expression levels were different ( Figure 6). For PbGID1s, with the exception of PbGID1c-1-2 which was expressed in all tissues at similar levels, other PbGID1s were expressed at very low levels in the roots and flowers compared to in other tissues. The highest expression levels of PbGID1c-1-1/1c-2 (Type A) were found in leaves and fruits while those of PbGID1b-1/1b-2 (Type B) were found in new shoots. However, PbDELLAs showed opposite trends to PbGID1s, whereby their expression levels were much higher in roots and flowers, but lower in new shoots. PbGAI1a, PbGAI1b, and PbGAI2b were also expressed at high levels in leaves, while PbGAI2b and RGLb were expressed at low levels in young fruits. The main expression site for PbSLYs was in the reproductive tissues. The expression levels of PbSLYs in flowers were similar to those of PbDELLAs, but opposite to PbGID1s.
Taken together, it was found that the common cis-acting elements were present in the promoter regions of the genes in the same family. However, there were also unique ones found in some genes within the same family. This suggests that, while genes in the same family play similar roles, each gene also has its unique function in growth, development, hormone response, and stress response.

PbGID1s, PbDELLAs, and PbSLYs Were Ubiquitously Expressed with Tissue-Specific High-Level Expression for Some Genes in Each Family
In order to understand the expression profiles of PbGID1s, PbDELLAs, and PbSLYs, real-time RT-PCR was carried out using total RNA isolated from roots, shoots, leaves, flowers, and young fruits of 5 year old trees. The results showed that the transcripts of all candidate genes were present in all the tissues tested; however, their expression levels were different ( Figure 6). For PbGID1s, with the exception of PbGID1c-1-2 which was expressed in all tissues at similar levels, other PbGID1s were expressed at very low levels in the roots and flowers compared to in other tissues. The highest expression levels of PbGID1c-1-1/1c-2 (Type A) were found in leaves and fruits while those of PbGID1b-1/1b-2 (Type B) were found in new shoots. However, PbDELLAs showed opposite trends to PbGID1s, whereby their expression levels were much higher in roots and flowers, but lower in new shoots. PbGAI1a, PbGAI1b, and PbGAI2b were also expressed at high levels in leaves, while PbGAI2b and RGLb were expressed at low levels in young fruits. The main expression site for PbSLYs was in the reproductive tissues. The expression levels of PbSLYs in flowers were similar to those of PbDELLAs, but opposite to PbGID1s. . Quantitative RT-PCR was carried out using total RNA isolated from different tissues as indicated. PbACTIN was used as the reference gene. The relative expression level of each gene was calculated relative to the transcript level of the roots (set as 1). Each bar represents the mean ± SE (n = 3). Different letters Quantitative RT-PCR was carried out using total RNA isolated from different tissues as indicated. PbACTIN was used as the reference gene. The relative expression level of each gene was calculated relative to the transcript level of the roots (set as 1). Each bar represents the mean ± SE (n = 3). Different letters above the bars indicate significant differences at p < 0.05 (n = 3) analyzed by Duncan's multiple range test.

Changes in Expression Levels of PbGID1s, PbDELLAs, and PbSLYs in Tissue-Cultured Seedlings Treated with GA 3 , PAC, IAA, ABA, and NaCl
To see if/how the expression profiles of PbGID1s, PbDELLAs, and PbSLYs changed in response to phytohormones and salt stress, tissue-cultured seedlings were treated with GA 3 , PAC, IAA, ABA, and NaCl for different lengths of time. The expression levels of PbGID1s, PbDELLAs, and PbSLYs were monitored by qRT-PCR. As shown in Figure 7, in GA 3 -treated seedlings, the highest expression levels of PbGID1c-1-1, 1c-1-2, and 1c-2 were reached at 12, 120, and 12 h, with 1.8-, 4.2-, and 3.6-fold increases compared to 0 h, respectively. The expression levels of PbDELLAs and PbSLYs were increased first and then decreased, with PbGAI2a/b decreasing to the lowest levels at 3 h after treatment. Similarly, the expression levels of all PbSLYs were also the lowest at 3 h, with PbSLY2-1 and 2-2 being the most affected by GA 3 .
When the seedlings were treated with PAC, an inhibitor of GA, all PbGID1s, with the exception of PbGID1c-1-2, were downregulated first and then upregulated, with PbGID1c-1-1 and PbGID1c-2 being the most sensitive and reaching the lowest expression levels with 5.3-fold and 6.8-fold decreases at 6h and 3h compared that at 0 h, respectively. On the contrary, PbDELLAs and PbSLYs showed opposite expression patterns to PbGID1s, i.e., their expression levels increased first and then decreased in PAC-treated seedlings. The highest expression levels of PbGAI2a and PbRGLb were reached at 3 h while those of all four PbSLYs were increased at 72 h, with PbSLY2-1 and 2-2 showing the most profound change with increases by 6.3-and 33.6-fold, respectively.
In IAA-treated seedlings the expression levels of all PbGID1s were increased except for PbGID1b-1. PbGID1c-1-1 and PbGID1b-2 were quickly affected, with the highest expression levels with 1.3-and 2.3-fold increases achieved after 3 h, whilst PbGID1c-1-2 was the most affected with a 3.3-fold increase after 120 h. The expression levels of all six PbDELLAs were decreased to start with, before increasing, with PbGAI1a, 1b, and RGLb being the most sensitive to IAA treatment with the lowest expression levels with 10.4-, 1.8-, and 5.5-fold decreases shown after 6 h. All PbSLYs were upregulated when treated with IAA, with PbSLY1-1/1-2 responding the fastest and reaching the highest expression levels at 6-12 h, while PbSLY2-1 and 2-2 were the most affected, with their expression levels increased by 7.4-and 28.1-fold, respectively.
To see if the expression of PbGID1s, PbDELLAs, and PbSLYs was affected by stress, ABA and NaCl were supplemented in the media. The effects of ABA and NaCl on the expression levels of these genes were very similar. In the case of PbGID1s, the initial downregulation of PbGID1b-1 and PbGID1b-2 was the most significant with 4.2-and 6.0-fold (ABA) and 7.1and 3.8-fold (Nacl) decreases. However, PbDELLAs showed opposite trends, whereby their expression levels first increased and then decreased. This was particularly true for PbGAI2b and PbRGLa/b. All PbSLYs, except for PbSLY1-2 whose transcript level was increased throughout the treatment period, showed an increase then decrease.
Therefore, these combined data indicate that GA signaling plays an important role in the phytohormonal and abiotic stress response during 'duli' pear seedling development. . Quantitative RT-PCR was carried out using total RNA isolated from tissue-cultured 'duli' seedlings treated with GA3 (2mg/L), PAC (2 mg/L), IAA (0.2 mg/L), ABA (4mg/L), and NaCl (0.6%) for 0, 3, 6, 12, 72, and 120 h. PbACTIN was used as the reference gene. The relative expression level of each gene was calculated relative to the transcript level of that gene at 0 h (set as 1). Each bar represents the mean ± SE (n = 3). Different letters above the bars indicate significant differences at p < 0.05 (n = 3) analyzed by Duncan's multiple range test.
When the seedlings were treated with PAC, an inhibitor of GA, all PbGID1s, with the exception of PbGID1c-1-2, were downregulated first and then upregulated, with PbGID1c-1-1 and PbGID1c-2 being the most sensitive and reaching the lowest expression levels with 5.3-fold and 6.8-fold decreases at 6h and 3h compared that at 0 h, respective- . Quantitative RT-PCR was carried out using total RNA isolated from tissue-cultured 'duli' seedlings treated with GA 3 (2mg/L), PAC (2 mg/L), IAA (0.2 mg/L), ABA (4 mg/L), and NaCl (0.6%) for 0, 3, 6, 12, 72, and 120 h. PbACTIN was used as the reference gene. The relative expression level of each gene was calculated relative to the transcript level of that gene at 0 h (set as 1). Each bar represents the mean ± SE (n = 3). Different letters above the bars indicate significant differences at p < 0.05 (n = 3) analyzed by Duncan's multiple range test.

Discussion
The 'duli' pear is a wild pear species native to China, which is widely distributed in the northern region of the country. While the fruits are small and not edible, its rather well-developed root system, vigorous growth, high tolerance to biotic and abiotic stress, and good compatibility with European and Asian pear trees in grafting make it an ideal rootstock. In particular, a dwarfed rootstock is very desirable in order to achieve high density, uniformed planting, and ease of pesticide spraying and harvest.
Previous studies in wheat, rice, and other crop plants clearly showed that the plant hormone GA plays important roles in regulating plant height and stature [18,[21][22][23]. As such, many GA insensitive dwarf mutants were isolated and used in agriculture to achieve high yield. Therefore, studying the GA signal transduction-related genes, expression patterns, and functions in 'duli' pear could provide valuable information for the molecular breeding of dwarf rootstocks.

The Genome of 'duli' Pear Encodes More Members of Each Family of Genes Related to GA Signal Transduction Pathway Than Other Plant Species
A total of 15 GA signal transduction-related genes were identified from the 'duli' pear genome ( Figure 1, Table 1). Compared to Arabidopsis, rice, and many other plants, the number of genes in each family was increased. This seems to be correlated to a genomewide doubling event [40]. AtGID1s in Arabidopsis have 13 domains, TWVLIS, LDR, FFHGGSF, HS, IYD, YRR, DGW, GDSSGGNI, GNI, MF, LDGKYF, WYW, and GFY, that are responsible for binding to GA or DELLA [41]. Similarly, our analysis of PbGID1s also identified TWVLIS, FFHGGSF, YRR, GDSSGNI, LDGKYF, and GFY motifs in their sequences (motifs 3, 7, 2, 5, 1, and 4, respectively, Supplementary Materials Figure S1a), indicating that PbGID1s are most likely typical GID1s, and that these important functional motifs are highly conserved among plant GID1s. It is noteworthy, however, that PbGID1c-1-1 and PbGID1c-1-2 shared the highest similarity of~94%; the only differences between them were as follows: (1) PbGID1c-1-1 had 20 more amino acids than PbGID1c-1-2 at the Nterminus; (2) only three amino acids were different within the remainder of the sequences. Both genes were also located on the same (12th) chromosome ( Figure 1). Therefore, they are most likely the same gene that has been misannotated. Further cloning and sequencing will confirm this conclusion.
Compared to five in Arabidopsis, only one in rice, and one in tomato, the 'duli' pear genome contained six PbDELLAs. Further analysis showed that they belonged to three groups with two in each group, which is very similar to the DELLAs in the closely related apple genome [38]. It is worth noting that the predicted protein sequence of PbGAI1a had two DELLA sequences, which were also longer than others ( Figure 3). We were puzzled by this, as this has not been previously reported in the literature. We, therefore, carried out RT-PCR, cloning, and sequencing of PbGAI1a. The results showed that the deduced protein sequence contained only one DELLA, not two DELLAs (data not shown). Therefore, the two DELLA domains of PbGAI1a in the database were most likely caused by error when the 'duli' pear genome sequence was assembled.
The genome of the 'duli' pear encodes four SLYs compared to two in Arabidopsis and one in rice, respectively [28,31]. All four PbSLYs contained the F-box conserved domains and the specific motifs 1, 2, and 3 (F-box, LSL, and GGF, respectively), which were shown to be essential for the function of SLY1/GID2 in Arabidopsis [30] and rice [28]. Therefore, the presence of the same motifs in PbSLYs suggests that they are also important for their function in pear. Motif 5 was only found in PbSLY1-1/2 (SLY-I) while motif 4 was found in PbSLY2-1/2 (SLY-II). The study of OsGID2/SLYI in rice showed that it contains a unique VR1 motif in its N-terminus, which is not shared with AtSLY1. However, this motif is not required for its function because both WT-OsGID2 and OsGID2-∆VR1 can rescue the phenotype of gid2 [28]. Therefore, whether the two groups of PbSLYs function differently in 'duli' pear is worth investigating in the future.
The combined results clearly demonstrated that the 'duli' pear genome encodes all three families of GA-related signaling proteins. Despite their small variations in some specific motifs/domains were found compared to their homologous proteins in other plant species, they all contained the conserved domains and motifs important for their proper function in GA signaling.

Putative Functions of GA Signal-Related Genes in 'duli' Pear
In general, GAs are actively synthesized in young developing tissues, such as new shoots and small fruits, while well-developed mature tissues contain low levels of GAs. In line with this, we found that PbGID1s were expressed at low levels in roots and flowers but high levels in new shoots and young fruits, while PbDELLAs and PbSLYs (except for roots) showed an opposite trend ( Figure 6). This indicates that GA can promote the expression of PbGID1s but inhibit that of PbDELLAs and PbSLYs. Similarly, the highest expression level of rice OsGID2/SLY1 was found in unopened flowers, where the highest level of active GAs was detected [28]. Different expression levels of PbGID1s, PbDELLAs, and PbSLYs were found when 'duli' pear seedlings were treated with IAA, GA 3 , and PAC, as well as ABA and NaCl (Figure 7). The Type A PbGID1s were upregulated by IAA and GA3 but downregulated by PAC, indicating that they function in the growth and development of 'duli' pear seedlings. Similar results were found in tomato plants, where the three GID1s were expressed in all tissues but only GID1a responded to GA treatment, confirming the important regulatory role of GID1a in seed germination, stem elongation, and leaf development [16]. Studies in Arabidopsis and peach showed that the Arabidopsis gid1a/gid1c double mutant and peach single mutant ppgid1c were all severely dwarfed, indicating that GID1a and/or GID1c play important roles in regulating stem elongation and, hence, plant height [13,17]. These results, together with the results from our homology analysis of these genes (Figure 2), suggest that the Type A PbGID1s are most likely involved in the regulation of vegetative growth and plant stature in 'duli' pear. On the other hand, the Type B PbGID1s were downregulated in seedlings treated with ABA and NaCl; thus, they may play critical roles in abiotic stress. In supporting this, Illouz-Eliaz et al. found that tomato GID1b, a Type B GID1 maintained its stability when tomato grew in uncontrollable and unfavorable field conditions [16].
Different DELLA proteins have different roles. For example, in Arabidopsis, AtGAI and AtRGA play a major role in vegetative growth, AtRGA, AtRGL1, and AtRGL2 play a major role in seed germination and flower development, and AtRGL3 plays a major role in stress response [26,27]. The three PbDELLAs, PbGAI1a, PbGAI1b, and PbGAI2b, were highly expressed in roots and leaves but lowly expressed in new shoots, while PbGAI2a, PbGAI2b, PbRGLa, and PbRGLb were highly expressed in flowers but lowly expressed in young fruits except PbRGLa ( Figure 6). Treatment with GA 3 resulted in the downregulation of PbGAI2a and PbGAI2b (3 h), while PAC caused the upregulation of PbGAI2a and PbRGLb (3 h) (Figure 7). PbDELLAs may also play important roles in stress response. This is because, when seedlings were treated with stress hormone ABA and salt, PbRGLa and PbRGLb showed the fastest and most increased expression levels, followed by PbGAI2a and PbGAI2b, while PbGAI1a and PbGAI1b remained unchanged (Figure 7). Therefore, these expression data, combined with the results from homology and evolutionary analysis (Figure 2), indicate that PbGAI2a, PbGAI2b, PbRGLa, and PbRGLb function in the regulation of reproduction, as well as in stress, while PbGAI1a, 1b and 2a play important roles in stem elongation, as well as leaf development, in 'duli' pear. It is worth noting that PbGAI2a and PbGAI2b were the only pair of genes without a syntenic relationship (Figures 2 and 4a), indicating that there may be some functional differentiation of these two genes.
The expression profiles of PbSLYs in IAA-, GA3-, PAC-, ABA-, and NaCl-treated tissuecultured seedlings showed similar trends to those of PbDELLAs but opposite trends to those of PbGID1s. Most interestingly, the type II PbSLYs (PbSLY2-1 and PbSLY2-2) showed the most significant changes, while type I PbSLYs (PbSLY1-1 and PbSLY1-2) showed the fastest changes following different treatments (Figure 7). AtSLYs in Arabidopsis also showed dif-ferent functions because the sly1sne double mutant a showed more severe dwarf phenotype and significantly reduced fertility compared to the single mutant sly1 [32]. Furthermore, overexpression of SNE (SLY2) can partially restore the GA-insensitive phenotype of sly1 mutants [42,43]. Therefore, SNE, a type II SLY, can replace SLY1, i.e., SLY2 functions redundantly with SLY1. However, overexpression of SLY1 in the sly1 mutant can lead to the degradation of AtRGL1 and AtRGL2, while excessive SNE/SLY2 results in no change in RGL2 and only a reduction in RGL1 levels. This indicates that the functions of SNE/SLY2 and SLY1 are not exactly the same [31]. Whether the two types of PbSLYs function in a similar fashion to those in Arabidopsis remains to be explored.
In summary, we identified 15 GA-related signaling genes from the wildtype 'duli' pear genome. Further bioinformatic analysis and expression studies of these genes showed that they are typical and comparable to those characterized from Arabidopsis and other plant species. The type A PbGID1s play a major role in the regulation of growth and development, while type B PbGID1s play a major role in stress. The members of type I and II PbDELLAs (PbGAI1a, PbGAI1b, PbGAI2a, and PbGAI2b) may be involved in stem elongation and leaf (vegetative) development, while type II and III PbDELLAs (PbGAI2a, PbGAI2b, PbRGLa, and PbRGLb), as well as all four PbSLYs, may be involved in reproduction and stress response (Supplementary Materials Figure S3). As such, these preliminary data could guide further research into GA signaling pathways and the identification of desirable genes for 'duli' pear rootstock development.

Plant Materials
Roots, shoots, leaves, flowers, and young fruits of 'duli' pear (Pyrus betulifolia Bunge) were sampled from three trees (5 years old) maintained in the 'duli' Resource Orchard of Hebei Agricultural University, Hebei, China. The stem segments (about 1-2 cm) of the same trees were removed and used as explants to generate tissue-cultured seedlings for hormonal and abiotic stress treatments, as described below.

Identification and Chromosome Location of the Candidate GA Signaling-Related Genes
To identify the genes in the GA signaling pathway in 'duli' pear, the known amino-acid sequences of Arabidopsis and/or rice GID1, DELLA, and SLY were used as references to search the whole genome sequence of 'duli' pear (GWHA,m>AYT00000000, the Genome Warehouse of the BIG Data Center, Beijing Institute of Genomics (BIG), Chinese Academy of Sciences) [44,45]. This was followed by sequence analysis to identify and verify the Pfam and conserved domains via the CDD (https://www.ncbi.nlm.nih.gov/cdd/, accessed on 6 December 2021) [46] and SMART (http://smart.embl.de/smart/batch.pl, accessed on 7 December 2021) databases, as well as conserved motifs via MEME (http://memesuite.org/tools/meme, accessed on 7 December 2021) [47]. The molecular weight (MW) and isoelectric point (pI) of each of the encoded proteins were predicted using the online software ExPASy (https://www.expasy.org/, accessed on 6 February 2022). The location of these genes on the specific chromosomes was mapped and marked using the software TBtools [45].

Evolutionary Analysis of the Putative GA Signaling-Related Genes
The amino-acid sequences of the above identified proteins in 'duli' pear were compared with those from other plants using Clustal (https://www.ebi.ac.uk/Tools/msa/clustalo/, accessed on 7 February 2022). A phylogenetic tree was constructed using the maximum likelihood method by applying the JTT matrix model and the MEGA6.0 software with default parameters and1000 bootstraps. This was imported to the online tool iTOLs (https: //itol.embl.de/itol.cgi, accessed on 7 February 2022) to produce the figure.

Synteny and Gene Duplication Analysis
The syntenic relationship of the genes with a role in the GA signaling between Arabidopsis and pear was analyzed. The orthologous and paralogous relationships between them were drawn using TBtools [45]. The Ka, Ks, and Ka/Ks values for the individual gene pairs were calculated by the adjoint method using the software KaKs_Calculator 2.0 (( https://sourceforge.net/projects/kakscalculator2/, accessed on (8 March 2022) [48].

Prediction and Analysis of cis-Elements in the Promoter Regions of the GA Signaling Genes of 'duli' Pear
In order to identify relevant cis-acting elements, the DNA sequences about 2000 bp upstream of the transcriptional start codon ATG of the PbGID1, PbDELLA, and PbSLY genes of 'duli' pear were analyzed using the PlantCARE online database (http://bioinformatics. psb.ugent.be/webtools/plantcare/html/, accessed on 12 December 2021) [49].

Total RNA Extraction and First-Strand cDNA Synthesis
The total RNA was isolated from roots, young shoots, leaves, flowers, fruits, and tissue-cultured 'duli' seedlings of pear plant using the RNAprep Pure Plant Kit according to the manufacture's protocol (TIANGEN, Beijing, China). The first-strand cDNA was synthesized from 1 µg of total RNA using the First-Strand cDNA Synthesis Kit (TIANGEN, Beijing, China).

qRT-PCR Assay
The transcript levels of the identified genes were analyzed by qRT-PCR using the above-synthesized cDNAs as templates, 2×TransStart ® Top Green qPCR SuperMix (Trans-Gen Biotech, Beijing, China), and gene-specific primers (Supplementary Materials Table S2) in a 20 µL reaction mix. The reaction was run on a Step-Two Plus real-time PCR system (Applied Biosystems, San Francisco, CA, USA). PbACTIN was used as the reference gene. At least three replicates were included in each experiment, and experiments were repeated three times.

Conclusions
Bioinformatic analysis showed that there were five GA receptors (PbGID1s), six GA negative regulators (PbDELLAs), and four GA positive regulators (PbSLYs) in the 'duli' pear genome. The expression profiles of these genes showed that, although they were ubiquitously expressed, their expression levels and trends were different. Exogenous GA and IAA application could increase the transcript levels of PbGID1s but decrease those of PbDELLAs and PbSLYs, while PAC, ABA, and NaCl showed the opposite effect on the expression of these genes. Future research should focus on the detailed functional analysis of these genes so that their involvement in plant height, reproduction, and stress response can be revealed.