Analysis of the Small Auxin-Up RNA (SAUR) Genes Regulating Root Growth Angle (RGA) in Apple

Small auxin upregulated RNAs (SAURs) are primary auxin response genes; the function of regulating root growth angle (RGA) is unclear in the apple rootstock. We firstly identified 96 MdSAUR genes families from new apple genome GDDH13 using the resequence database of ‘Baleng Crab (BC)’ and ‘M9’. A total of 25 MdSAUR genes, regulating the formation of RGA, were screened for the expression profiles in stems and roots and the allelic variants of quantitative trait loci (QTL). Finally, through the joint analysis of network and protein–protein interaction, MdSAUR2, MdSAUR29, MdSAUR60, MdSAUR62, MdSAUR69, MdSAUR71, and MdSAUR84 were screened as the main candidate genes for regulating RGA. This study provides a new insight for further revealing the regulatory mechanism of RGA in apple dwarf rootstocks.


Introduction
Apple (Malus domestica Borkh.) is one of the most popular perennial tree fruits in the world [1]. Currently, the dwarf form and close planting is the main cultivation mode, which is used to develop the apple industry. However, dwarfing rootstock is the limiting factor that restricts the popularization of dwarf and close planting [2]. An ideal dwarfing rootstock with a relatively deep root system is pivotal for water uptake, nutrient uptake, and adaptability [3]. Root growth angle (RGA), defined as the angle between the direction of root growth and the horizontal axis, determines the depth of root architecture and the development of high-density planting of apple [4]. The bending of RGA results from differential elongation of the cells in the distal elongation zone, which results from the differential accumulation of auxins on opposite sides of the root [5,6]. However, the underlying mechanisms of auxin in RGA is unclear.
Auxin mainly includes three different types of early auxin-responsive gene families, namely, Aux/indole-3-acetic acid (IAA), Gretchen Hagen 3 (GH3), and small auxin-up RNA (SAUR) [7]. Auxin signal transduction by numerous transcription factors that act downstream of IAA and auxin response factor (ARF) proteins leads to a variety of auxininduced morphological changes in plants, including gravitropic bending [8,9]. Precise auxin levels are vitally important to plants, which have many effective mechanisms to maintain auxin homeostasis [10]. The GH3 gene family, a supervisor of the fluctuation of auxin, could maintain auxin homeostasis by a mechanism in which an amino acid of some GH3 proteins conjugates to excessive IAA and inactivates IAA [11]. The SAUR gene family comprises the most rapid auxin-responsive genes related to the auxin signaling pathway [7,12]. The expression of SAURs can be induced within minutes by the active auxin, implying that auxin plays a critical role in their gene transcription [13]. However, many SAURs are regulated post-transcriptionally as a conserved downstream destabilizing element in the 3 -untranslated region that confers high mRNA instability [14]. Therefore, regulation of SAURs may occur at the transcriptional, posttranscriptional, and protein levels [15,16].
After the first SAUR gene was identified in soybean [17], members of SAUR gene family have been identified in a wide range of plants such as rice [18], moso bamboo [19], and poplar [20]. In apple, a total of 80 MdSAUR genes have been successfully identified in the draft of the apple genome [21,22]. Now, a new, high-quality, de novo assembly of the apple genome, named GDDH13, has been updated [23]. More accurate gene sequences and different designations were identified in the new apple genome GDDH13 [24]. Therefore, it is necessary to further identify the MdSAUR gene family in this new apple genome. In our previous study, the BSA-seq and RNA-seq of hybrids from 'M9' and 'BC' were analyzed, which contributes to further screening of MdSAUR genes regulating RGA from expression and allelic variation [25].
In this study, a genome-wide analysis for the homologues SAUR genes in the new apple genome, GDDH13, was carried out. Subsequently, the expression profiles of Md-SAURs in the stems and roots of 'M9' and 'BC' after cutting were measured by quantitative real-time PCR (qRT-PCR). Then, the allelic variants of MdSAURs located in QTLs were predicted by the resequence database of 'M9' and 'BC'. Finally, these network analyses of MdSAURs were analyzed to further screen MdSAURs genes regulating RGA.

Identification of SAUR Gene Family in Apple Genome GDDH13
According the resequence database of 'BC' and 'M9', a total of 96 SAUR genes in the new apple genome GDDH13 were identified and designated as MdSAUR 1-96, ordered by their location in the chromosome from the top to bottom of the apple chromosomes 0-17 ( Table 1). The number of MdSAUR genes per chromosome ranged from 1 (chromosome 6,11,17) to 25 (chromosome 10). The protein lengths of coding sequences range from 51 (MdSAUR 5) to 219 amino acids (MdSAUR 71), with an average sequence length of 124 amino acids. The lengths of MdSAUR genes in the genome vary from 156 to 660 bp. According to the results of the predicted protein localization, we found that 93 SAUR proteins possess signal sequences targeting in the cytoplasmic. MdSAUR 33 and MdSAUR 71 were only located in periplasmic. MdSAUR 56 was in the outer membrane and periplasmic.

Phylogenetic Relationships of SAUR Genes in Arabidopsis and Apple
A phylogenetic tree was constructed on the basis of 231 putative nonredundant SAUR protein sequences from Arabidopsis, rice, and apple. All SAUR proteins were clustered into five groups ( Figure 1, Table S1). The number of SAUR proteins per group ranged from 26 (group V) to 64 (group I). Groups II, III, and V contained more SAUR proteins in apple than in Arabidopsis and rice, whereas the opposite was found in group I. Group IV contained 14 MdSAUR proteins, which was more than the number of OsSAURs (9), and less than the number of AtSAURs (20).

Expression Analysis of MdSAUR Genes in Stems and Roots
A total of 36 MdSAURs genes had no expression level in stems of hybrids from 'M9′ and 'BC' at 7 d after cutting by RNA-seq analysis (Table S2). Of 60 expressed MdSAURs, 18 had a higher transcripton level in hybrids of large RGA than those of small RGA (Table  S2; Figure 2a). At 14 d after cutting, roots of hybrids from 'M9′ and 'BC' owned 53 expressed MdSAURs, but didn't contain the different expression genes (Table S2; Figure 2b). These candidate MdSAUR genes were further screened by using qRT-PCR analysis of 'M9′ and 'BC' (Figure 3). The expression level of 18 MdSAUR genes in stems of 'BC' (large RGA)

Expression Analysis of MdSAUR Genes in Stems and Roots
A total of 36 MdSAURs genes had no expression level in stems of hybrids from 'M9' and 'BC' at 7 d after cutting by RNA-seq analysis (Table S2). Of 60 expressed MdSAURs, 18 had a higher transcripton level in hybrids of large RGA than those of small RGA (Table S2;

Co-Localization of MdSAURs with QTL for RGA
To better understand the potential function of SAUR genes related to RGA, we colocalized the SAUR genes with reported RGA QTL. As a result, a total of 24 MdSAUR genes were mapped in QTL M10.1 and B13.2 by BAS-seq analysis of hybrids with largesmall RGA from 'BC' and 'M9' (Figure 4). To understand if the co-localized SAURs are genetically associated with RGA, allelic variations were analyzed to gain an insight into the MdSAUR genes in QTLs. In the QTL M10.

Co-Localization of MdSAURs with QTL for RGA
To better understand the potential function of SAUR genes related to RGA, we colocalized the SAUR genes with reported RGA QTL. As a result, a total of 24 MdSAUR genes were mapped in QTL M10.1 and B13.2 by BAS-seq analysis of hybrids with largesmall RGA from 'BC' and 'M9′ ( Figure 4). To understand if the co-localized SAURs are genetically associated with RGA, allelic variations were analyzed to gain an insight into the MdSAUR genes in QTLs. In the QTL M10.

Network Analysis of Candidate MdSAURs
We predicted potential proteins which may regulate 25 candidate MdSAUR genes by network analysis in the AppleMDO online database. Based on the Pearson correlation coefficient (PCC) ≥ 0.9 and the numbers of networks ≥ 10, six MdSAURs were predicted as primary candidate genes (Figure 7, Table S3). The numbers of potential proteins with

Network Analysis of Candidate MdSAURs
We predicted potential proteins which may regulate 25 candidate MdSAUR genes by network analysis in the AppleMDO online database. Based on the Pearson correlation coefficient (PCC) ≥ 0.9 and the numbers of networks ≥ 10, six MdSAURs were predicted as primary candidate genes (Figure 7, Table S3). The numbers of potential proteins with MdSAURs ranged from 21 (MdSAUR2) to 67 (MdSAUR71). MdSAUR2, MdSAUR69, and MdSAUR84 formed three independent networks with different proteins of 21, 31, and 10, respectively. In addition, MdSAUR60, MdSAUR62, and MdSAUR71 constituted three interleaved networks containing 20 mutual proteins and 144 individual proteins.

Protein-Protein Interaction Analysis of Candidate MdSAURs
The corresponding functional protein-protein interaction networks were reconstructed using the 25 candidate MdSAUR proteins to explore the functional protein-protein interactions and the gene regulatory relationships among these MdSAUR proteins (Figure 8). MdSAUR2, MdSAUR29, and MdSAUR71 proteins were regarded as crucial node proteins with big nodes and tight coordination for further research. In three specific proteins (or nodes), MdSAUR29 was the core hub-protein with big nodes and tight coordination in the mapped network. However, another 22 proteins had no interaction or connection to other proteins in the network.

Function Analysis of Primary Candidate MdSAURs
Through the network analysis and protein-protein interaction analysis, a total of seven SAURs were selected as primary candidate genes for regulating RGA formation for future study (Figure 6). MdSAUR2, MdSAUR29, MdSAUR60, and MdSAUR62 had a higher expression level in stems of 'BC' than in those of 'M9′ ( Figure 3). MdSAUR60, MdSAUR62, MdSAUR69 and MdSAUR71 were mapped in QTL M10.1, and MdSAUR84 was located in QTL B13.2 ( Figure 4). The previous kompetitive allele-specific PCR (KASP) assay showed [25], the RGA of hybrids with the A:G genotype of marker Z3658 in QTL M10.1 were large, but had no significant difference from those with the A:A genotype ( Figure S1). However, significantly larger RGA values were detected in hybrids with T:G than in hybrids with T:T genotype of the marker b13 in QTL B13.2 ( Figure S1).

MdSAUR Gene Family in Apple Genome GDDH13
In this present study, we successfully identified 96 SAUR genes in the new apple genome GDDH13, which has more SAUR genes and a more accurate prediction sequence than the old apple genome (80 members) [21,22]. Compared with the SAUR gene family in Arabidopsis (72), rice (58), and sorghum (71), many more members were identified in apple, suggesting that the SAUR family in apple experienced severe expansion during the long evolutionary history [26][27][28]. Some SAUR proteins were localized in the nucleus, cytoplasm, or plasma membrane [16,29,30]. Our results indicated that more than 97% of the SAUR proteins localized in the cytoplasm (Table 1). This feature has also been found in

Function Analysis of Primary Candidate MdSAURs
Through the network analysis and protein-protein interaction analysis, a total of seven SAURs were selected as primary candidate genes for regulating RGA formation for future study ( Figure 6). MdSAUR2, MdSAUR29, MdSAUR60, and MdSAUR62 had a higher expression level in stems of 'BC' than in those of 'M9' (Figure 3). MdSAUR60, MdSAUR62, MdSAUR69 and MdSAUR71 were mapped in QTL M10.1, and MdSAUR84 was located in QTL B13.2 ( Figure 4). The previous kompetitive allele-specific PCR (KASP) assay showed [25], the RGA of hybrids with the A:G genotype of marker Z3658 in QTL M10.1 were large, but had no significant difference from those with the A:A genotype ( Figure S1). However, significantly larger RGA values were detected in hybrids with T:G than in hybrids with T:T genotype of the marker b13 in QTL B13.2 ( Figure S1).

MdSAUR Gene Family in Apple Genome GDDH13
In this present study, we successfully identified 96 SAUR genes in the new apple genome GDDH13, which has more SAUR genes and a more accurate prediction sequence than the old apple genome (80 members) [21,22]. Compared with the SAUR gene family in Arabidopsis (72), rice (58), and sorghum (71), many more members were identified in apple, suggesting that the SAUR family in apple experienced severe expansion during the long evolutionary history [26][27][28]. Some SAUR proteins were localized in the nucleus, cytoplasm, or plasma membrane [16,29,30]. Our results indicated that more than 97% of the SAUR proteins localized in the cytoplasm (Table 1). This feature has also been found in other species, such as cotton, watermelon, and maize [31][32][33]. SAUR clusters were reported in rice, tomato, and maize, in which most genes tended to be grouped together by phylogenetic analysis [27,31,34]. The SAUR proteins in the old apple genome were compared with those in Arabidopsis and rice, which constructed a phylogenetic tree containing six groups [21,22]. In our study, SAUR genes of the new apple genome GDDH13 were also grouped together with Arabidopsis and rice by phylogenetic analysis, which contained five groups and indicated that species specificity commonly existed in the SAUR family.

Differential Expression of MdSAUR Genes
The high expression divergence in different tissues reflected the complexity of gene family functions [35]. Given that the expression profiles of SAUR genes could provide important clues for gene function, the expression of SAUR genes in different tissues (leaf, stems, root, and flowers) were analyzed in maize, watermelon, and cotton [31][32][33]. Previous reports found that SAUR genes participating in the auxin signal pathway were upregulated or repressed to some extent following an auxin treatment [18,36]. In our study, expression profiles of MdSAUR genes were investigated in stems and roots after cutting by RNA-seq and qRT-PCR analysis. In fact, 61 MdSAUR genes were found to be expressed in the stems or roots, suggesting that the MdSAUR gene family might play a major role in the RGA development in apple. A total of 18 MdSAUR genes were found to have differential expression in stems, which acted as candidate genes for regulating RGA formation. Furthermore, the expression levels of MdSAUR70 and MdSAUR71 were high in stems, suggesting that they played key roles in RGA development. Similar to other primary auxin-responsive gene families (AUX/IAA, GH3), some members of SAUR gene family also exhibited tissue-specific expression patterns [34].

Allelic Variations and Gene Structures of MdSAUR Genes in QTLs
Auxin is essential for root growth and development, regulating stem cell specification and division, meristem size, cell elongation, and differentiation [37]. With SAURs being the most rapid auxin-responsive genes, most studies have focused on expression patterns and functional studies in shoots; several SAUR overexpression lines also exhibit root phenotypes [18]. Many QTLs of RGA have been reported in diverse species, such as sorghum (Sorghum bicolor L.), rice (Oryza sativa L.), and wheat (Triticum aestivum L.), which contributed to screening RGA by allelic variations [38][39][40]. In our study, 13 MdSAURs of QTL were screened by the allelic variations analysis based on the resequence database of 'BC' and 'M9', suggesting a likely genetic involvement of these SAUR genes in RGA formation.
The majority of SAURs lacked introns. For example, most of the SAUR genes in Solanaceae species possessed no introns and only 9 out of 99 SAURs in tomato and 3 out of 134 SAURs in potato had introns in their coding regions [34]. This phenomenon also exists in the SAUR genes of other species [28,33,41]. A similar phenomenon was found in the old apple genome, such that approximately 55% of the SAUR genes had no intron [22]. In our study, of 13 candidate SAURs in QTLs, only MdSAUR69 had an intron in the new apple genome. As the occurrence of alternative splice in intronless genes is usually low, the function of certain SAUR family genes is likely stable.

Primary Candidate MdSAURs for RGA Formation
Network analysis of AppleMDO webtools could predict some proteins that positively or negatively regulate the target protein, which contributes to finding the primary target proteins with more regulatory proteins [25,42]. In our study, 25 MdSAURs, screened from expression and allelic variations, were further analyzed based on the network analysis of AppleMDO webtools. MdSAUR2, MdSAUR60, MdSAUR62, MdSAUR69, MdSAUR71, and MdSAUR84 were screened as primary candidate genes because these MdSAURs regulated more proteins with PPC ≥ 9.0.
The STRING database aims to integrate all known and predicted associations between proteins, including both physical interactions and functional associations [43]. In the old apple genome, a total of 24 specific MdSAUR proteins were the core hub-proteins that were regarded as crucial nodes or core hub proteins with big nodes and tight coordination for further research [22]. In our study, these 25 key SAUR proteins were used to further identify the regulatory relationships with these putative functional protein-protein interactions using the STRING online database. Only three SAUR proteins (MdSAUR2, MdSAUR29 and MdSAUR71), which had big nodes and tight coordination in the mapped network, were selected as primary candidate genes.
Natural variation affects RGA formation, which has been verified in rice, corn, and wheat [43][44][45]. In our research, the allelic variations of MdSAURs (locatingin QTL) probably affected the gene function. The allelic variations of CDS in MdSAUR62, MdSAUR69, MdSAUR71, and MdSAUR84 affected the amino acid sequences, and the allelic variants of MdSAUR62 upstream could result in the differential expression between 'BC' and 'M9' in stems. Furthermore, significantly larger RGA values were detected in hybrids with T:G than in hybrids with T:T genotype of the marker b13 in QTL B13.2.

Sequence Retrieval and Characterization Analysis
MdSAUR family numbers were achieved using the resequence database of 'BC' and 'M9' [25,46]. The locations and sequences of all MdSAURs were found via searching against the gene database of AppleMDO (http://bioinformatics.cau.edu.cn/AppleMDO/index. php (accessed on 22 October 2020) using the identity document (ID) of the gene [23]. Subcellular localization prediction of each of the MdSAUR family genes was carried out using the CELLO v2.5 server (http://cello.life.nctu.edu.tw/ (accessed on 22 October 2020) [47].

Phylogenetic Analysis in Arabidopsis, Rice, and Apple
Multiple sequence alignments for all available SAUR full-length protein sequences of Arabidopsis and apple were performed using the gene database of AppleMDO and TAIR (https://www.arabidopsis.org/index.jsp (accessed on 2 November 2020)). The rice SAUR genes were searched according to the previous report [27]. A neighbor-joining (NJ) phylogenetic tree of full-length sequences of SAURs in apple, rice, and Arabidopsis was constructed by the software MEGA 5.0 [48].

Gene Expression and qRT-PCR Analysis
In RNA-seq of hybrids from 'BC' and 'M9', the expression patterns of MdSAURs were analyzed at 7 d (stems) and 14 d (roots) after cutting [25]. To further characterize the expression of selected MdSAUR genes, tissue samples were also collected from 'BC' and 'M9' at 7 d (stems) and 14 d (roots) after cutting. Total RNA of stems and roots was extracted by using the CTAB method [49]. The first strand of cDNA was synthesized using the PrimerScript™ RT Reagent Kit with the gDNA Eraser (Perfect Real Time) kit from Takara (Tokyo, Japan). Specific primers were designed for the selected MdSAUR genes using Primer 5.0 software, and the actin gene of apple was used as a standardized internal reference (Table S4). Expression levels of MdSAURs were determined by real-time quantitative PCR (RT-qPCR) [50]. Statistical analysis was performed using SPSS 17 software (Armonk, NY, USA).

Gene Structure and Allelic Variations Analysis
Mapchart 2.2 software was used to generate the chromosomal location image for these MdSAUR genes in QTL [51]. By comparing the resequencing database of 'M9' and 'BC', MdSAUR genes with SNPs or SVs, affecting cis-elements of the promoter and amino acids, were screened as candidate genes. The gene structures of these selected MdSAURs were found through the apple genome (GDDH13): JBrowse (https://www.rosaceae.org/tools/ jbrowse (accessed on 1 December 2021) [23].

Network Analysis
Candidate MdSAUR genes, regulating RGA, were screened by the joint analysis of expression and allelic variation. To further screen the primary MdSAUR genes, network analysis of these candidate MdSAURs was analyzed using AppleMDO webtools [24,42]. Moreover, amino acid sequences of these MdSAURs were also used to analyze the protein-protein interactions using the STRING (https://cn.string-db.org/ (accessed on 12 May 2022)) [22].

Statistical Analysis
Expression differences of MdSAUR genes in stems and roots of 'BC' and 'M9' were analyzed by one-way analysis of variance (ANOVA) through Dunnett's multiple comparison at a significance level of α = 0.05.

Conclusions
This study provides the phylogenetic analysis of the SAUR gene family in the new apple genome for the first time. In total, 25 MdSAURs regulating RGA were selected by the analysis of RNA-seq and BSA-seq. Then, combining the network analysis and protein-protein interaction analysis, seven SAURs (MdSAUR2, MdSAUR29, MdSAUR60, MdSAUR62, MdSAUR69, MdSAUR71, and MdSAUR84) were selected as primary genes for regulating RGA formation, to further study.
Supplementary Materials: The following supporting information can be downloaded at: https: //www.mdpi.com/article/10.3390/genes13112121/s1: Table S1 Predicted SAUR protein sequences in apple, Arabidopsis, and rice; Table S2 Tissue-specific expressions of apple SAUR genes in RNA-seq database of hybrids from 'BC' and 'M9' after cutting; Table S3 Co-expressed network analysis between MdSAURs and regulatory proteins; Table S4 Primers of apple on 18 SAUR genes used for qRT-PCR; Figure S1 Box plots showing differences in RGA between hybrids with genotypes of marker in QTL. Numbers of the hybrids are presented in parentheses. Asterisks represent p < 0.05 by Dunnett's multiple comparison.  Data Availability Statement: All RNA-seq reads and BSA-seq data have been deposited in the NCBI SRA under the accession number PRJNA655783. The DNA resequencing reads of 'Baleng Crab' and 'M9' are freely available in the NCBI SRA under the accession numbers SRR12234087 and SRR12233711, respectively. The data presented in this study are available in the graphs and tables provided in the manuscript.