By referring to the expression of genes related to male and female flowering as described by Jia et al. (2019), three experimental materials with different flowering stages were selected from the Lanxi International Chinese bayberry Research Center(Latitude 29.30°N, longitude 119.60°E) in November 2019 (period during which floral buds can be morphologically differentiated). These included the early maturing variety “Zaojia” (ZJ), the medium maturing variety “Biqizhong” (BQ) and the late maturing variety “Dongkui” (DK) (Table 1). The ages of all selected trees were around 15 years, and they were in consistent cultivation conditions. Each variety was sampled in triplicates, and for transcriptome sequencing, young leaves (not unfolded) were taken from annual branches facing south and at 1 m above the ground.

A polysaccharide polyphenol RNA extraction kit was used for extracting total RNA from the samples (TIANGEN, Beijing), and after RNA detection, Biomarker Biotechnology Co., Ltd. was commissioned to carry out the transcriptome sequencing. For this purpose, magnetic beads with Oligo (dT) were used to enrich the total RNA of the samples before fragmenting the mRNA with the fragmentation buffer. The first cDNA strand was then synthesized with random hexamers using mRNA as template, and this was followed by the synthesis of the second cDNA strand by adding buffer, dNTPs, RNase H and DNA polymerase I. After purification with the QiaQuick PCR kit and elution with EB buffer, end repair was performed, poly (A) tails were added and sequencing adaptors were connected. Appropriate fragments were then selected by gel electrophoresis prior to PCR-based amplification. The resulting libraries were eventually sequenced on an Illumina HiSeq4000.

After gene splicing, protein sequences were aligned with those from eight public databases (COG, GO, KEGG, KOG, Pfam, Swissport, eggNOG and Nr) using a threshold of e≤e^-10. The BLAST algorithm was then used for sequence similarity comparison, with the resulting sequence similarities subsequently used for functional annotations. Relative gene expression was assessed based on RPKM (Reads Per Kilobase of exon model per Million mapped reads) where larger RPKM values were indicative of higher expression levels (Trapnell et al., 2010).

Differentially expressed genes were screened by the false discovery rate (FDR) (Zhao et al., 2020), with a |log₂ fold change|≥2 and an FDR<0.5 selected as thresholds for a gene to be considered as being differentially expressed.

The SPLs of Chinese bayberry were isolated and identified by tBLASTn analysis of AtSPL amino acid sequences obtained from the genomic data of Chinese bayberry (Ren et al., 2019). The Chinese bayberry SPLs and target sites of miRNA156 were then predicted and confirmed using Genscan Web (http://genes.mit.edu/GENSCAN.html) as well as the BLASTx algorithm (http://www.ncbi.nlm.nih.gov/BLAST). After obtaining the nucleotide and amino acid sequences of 16 Arabidopsis and 46 apple SPL family genes from the Plant Transcription Factor Database (PlnTFDB3.0) (http://plntfdb.bio.uni-potsdam.de/v3.0/), phylogenetic trees were also constructed using the NJ method in MEGA 7.0, along with full-length protein sequences and the test parameter (bootstrap) set to 1000. The exon and intron structures of Chinese bayberry SPL genes were obtained by Gene Structure Display Sever (http://gsds.cbi.pku.edu.cn/index.php).

Escherichia coli competent cell DH-5α (Shanghai Jinchao Technology Development Co., LTD.), Agrobacterium tumefaciens strain GV3101 and pCambia2301-KY vectors (Shanghai Kaiyi Biotechnology Co., LTD.) were the main requirements of the study.

Using the genome sequence of Chinese bayberry (Ren et al., 2019), specific primers for both sides of the open reading framework (ORF) of the target gene were designed with Primer Premier 5.0 software for gene cloning. Total RNA extraction was also performed on healthy Chinese bayberry leaves using the modified cetyl trimethyl ammonium bromide (CTAB) method, with the extracted RNA acting as template to synthesize cDNA according to the instructions of the HiScript 1st Strand cDNA Synthesis Kit (Vazyme). This was followed by PCR amplification with the Phanta Max ultra-fidelity DNA polymerase (Vazyme), using the cDNA as template. In this case, each reaction consisted of the following component: 1 μL of Phanta Max super-Fidelity DNA Polymease, 2 μL of cDNA, 2 μL each of both forward and reverse primers, 25 μL of 2×Phanta Max Buffer, 1 μL of dNTP Mix and ddH₂O (for making up the volume to 50 μL), while the PCR procedure involved an annealing temperature of 49 °C and an extension rate of 1 kb/min, carried out for 39 cycles. Other operations shall follow the product instructions of Vazyme Company. The amplified products were finally detected on 1.5% agarose gel, before being sent to the company for sequencing to verify the accuracy of cloning results.

The restriction enzyme BamHI (Takara Company) was first used to linearize the vector before extracting the pCambia2301-KY plasmid for digestion with the same enzyme. The overexpression vector was then constructed with Vazyme recombinant enzyme at 37 °C by using the following components: 2 μL of 5×CE II Buffer, 1 μL of Exnase II, 4 μL of linearized carrier, 1 μL of insert fragment and ddH₂O (to make up the volume to 10 μL). After 30 min of reaction, the vector was placed on ice for cooling. The cells were then transfected into competent E. coli DH-5α cells and cultured in LB medium containing 50 mg/L Kan. This enabled the selection and subsequent culture of resistant colonies for the positive detection of the gene by PCR. The amplified products were finally sent for sequencing. The positive transformer colony plasmid was extracted and transfected into Agrobacterium tumefaciens GV3101 and sterile glycerol was added to preserve the bacteria at -80 °C until required for the next transfection.

Tobacco Benn was selected for this set of transformation experiment. WT tobacco was infected with Agrobacterium carrying recombinant vector plasmids of target genes using the leaf disk method. After four times of continuous screening/subculture, resistant buds were eventually recovered and transferred to a rooting medium to induce roots. Once the root system was vigorous, healthy and completely regenerated plants were transplanted to the soil (nutrient soil-vermiculite ratio was 1:1 or 2:1) where they were maintained until the T₀ generation for seed collection.

Collected seeds were sterilized with 70% ethanol, 30% sodium hypochlorite or 40% of 84 disinfectant and rinsed with sterile water 5-6 times. Seeds were then added to 1/2 MS solid selective medium containing 80 mg/L of Kan, and vernalized at 4°C for 2 days to break dormancy. They were subsequently cultured in a light incubator of the laboratory of Zhejiang Academy of Agricultural Sciences (light 28°C, 16 h, Darkness 25°C, 8 h, humidity 50%-70%). After about a week, the seeds were transferred to the soil to maintain grow. The leaves of the transgenic resistant plants and the wild-type ones of the T₁ generation were randomly sampled and stored at -80°C after being frozen in liquid nitrogen.

Transgenic positive plants of T₁ generation were obtained through screening with 80 mg/L Kan, 1/2 MS solid selective medium and PCR. The leaves of grown plants were collected, and total RNA was extracted with the RNA simple Total RNA Kit (TIANGEN) after quick-freezing in liquid nitrogen. In addition, synthesis reactions were also performed in 10-μL reaction volumes with the first Strand cDNA synthesis kit. For this purpose, the following components were used as required by the FastFire qPCR PreMix (SYBR Green) Kit (TIANGEN): 5 μL of 2×FastFire qPCR PreMix, 1 μL of forward primer and reverse primers (10 μm) and 1 μL of cDNA template. The reaction was performed on a Light Cycler 96 real-time PCR instrument under the following conditions: 95°C for 60 s, followed by 45 cycles, each at 95°C for 5 s, 63°C for 10 s and 72°C 15 s. Three technical replicates were set for each sample. Quantitative primers were designed according to gene sequences, with Ntactin-F/R selected as the reference gene (Zhao et al., 2020), and WT tobacco acting as the control to determine the relative expression of target genes. The 2^-△△Ct method was used to process the data (Livak and Schmittgen, 2001), while the IBM SPSS Statistics 22 and Origin 2022/Microsoft Excel 2010 were used for statistical analysis and plotting respectively.

Primers were therefore designed based on the reference genome sequence (Table S1), with ORF sequences of the MrSPL4 gene in ZJ, BQ and DK subsequently obtained by PCR amplification. Therefore, it was inferred that the expression of this gene was different between the different test materials, probably due to the promoter element, but this remained to be experimentally verified. The full length and CDS of MrSPL4 were 1,664 bp and 555 bp respectively. The amplified product was first recovered, and the vector was digested with BamHI. The resulting enzymatically digested product was then recombined with the amplified product to construct a plant overexpression vector. The latter was transformed into E. coli competent cells DH-5α before identifying the transformed bacterial solution by PCR.

Agrobacterium-mediated transformation of Nicotiana benthamiana with the recombinant plasmid 35S::MrSPL4-pCambia2301-KY was performed. The tobacco leaves infected by Agrobacterium tumefaciens were directly transferred to a selective medium containing kanamycin (Kan) to induce differentiation and budding. After the buds had grown to 2-3cm, they were inserted into a rooting medium to induce the formation of roots. Once the root system was vigorous, the seedlings were then tempered and transplanted to soil for culture to obtain completely regenerated tobacco with Kan resistance. Leaf DNA from the resistant regenerated tobacco plants to be tested was used as a template for PCR-based validation.

The transcriptome sequencing of nine samples of young leaves (three biological replicates for each variety) in the floral bud morphological differentiation period was completed, and a total of 59.78 Gb of clean data, with an average GC content of 47.28% and a Q30 base ratio of 93.55%, were obtained. After comparison with the reference genome (Ren et al., 2019), the percentage of clean reads aligned to the reference genome was found to be 95.39% (Table S2, the transcriptome data of BQ and DK were uploaded to https://bigd.big.ac.cn/gsa/browse/CRA008253, and the datasets of ZJ generated and analyzed during the current study are available in the NCBI repository https://www.ncbi.nlm.nih.gov/sra/PRJNA733585). The number of differentially expressed genes between the three samples was then compared. In this case, 623 genes were differentially expressed between ZJ and DK, and of these, 476 were up-regulated and 147 were down-regulated. Similarly, 2,343 genes were differentially expressed between ZJ and BQ, with 1,385 and 958 genes being up-regulated and down-regulated respectively. Finally, the number of differentially expressed genes between DK and BQ was 1,572, with 734 and 838 being up-regulated and down-regulated respectively (Table 2).

Through KEGG enrichment analysis, it was found that the differentially expressed genes mainly involved functions such as cellular processes, environmental information processing, genetic information processing, metabolism and organic systems (Figures 1A–C). The pathways that were significantly enriched in all the three groups included phytohormone signal transduction, sulfur and carbon metabolism, fatty acid, phenylpropionic acid, pyruvic acid, α-linolenic acid metabolism, glycine, serine, threonine, arginine, proline, cyano amino acids, cysteine and methionine metabolism, terpenoid skeleton biology, carotenoid biosynthesis, glutathione and glycerophosphatide metabolism, glycolysis/gluconeogenesis, starch and sucrose, amino sugar and nucleotide sugar metabolism, protein processing in the endoplasmic reticulum. The results indicated that these pathways may participate in the regulation network of Chinese bayberry flowering or other important pathways.

The differentially expressed genes mentioned above were compared with 306 flowering genes reported in Arabidopsis (Bouché et al., 2016) and 25 genes were found to be homologous in Chinese bayberry (Figure 2A). In particular, one of the differentially expressed genes, MRNA_003335_1, was down-regulated in DK but up-regulated in ZJ and BQ. This gene also contained the SBP domain and belonged to the SPL gene family, named MrSPL4. The relative expression of MrSPL4 was therefore verified by qRT-PCR (Figure 2B), and the results showed that ZJ had the highest expression, followed by BQ, with DK showing the lowest expression level. These results were, in fact, consistent with the expression determined by the transcriptome.

Through the screening of all genes in the reported genome provided in Ren et al. (2019), 17 SPL family genes with SBP domains were found in the Chinese bayberry genome (Table 3). Of these, 12 genes including MrSPL4, contained the target site of miR156 in the CD region. The software MEGA7.0 was then used to analyze the evolution of 17 SPL genes in Myrica rubra (MrSPL), 16 SPL genes in Arabidopsis thaliana and 46 SPL genes in apple (Figure 3). In this case, it was observed that the SPL proteins could be divided into four different groups (I, II, III and IV), with each containing at least one MrSPL gene. More specifically, groups I and II contained one MrSPL each, group IV contained twelve MrSPLs and group III contained three Chinese bayberry SPL genes, including MrSPL4. MrSPL4 had the highest homology with the AT1G20980 (AtSPL14) gene in Arabidopsis thaliana, with previous studies showing that this gene (AtSPL14) not only promoted the normal growth and development of Arabidopsis thaliana, but also played a crucial role in the development of flowering as well as the transformation from a vegetative to reproductive growth. Thus, it was speculated that MrSPL4 could be playing an important role in the flowering process of Chinese bayberry.

The electrophoresis results (Figure 4A) showed that the bands were consistent with the expected amplification product size, hence indicating that the sequence of the coding region of MrSPL4 was successfully obtained. No differences in the gene sequence were noted between the three samples (Figure 4B). MrSPL4 also contained two exons and one intron (Figure 4C), along with a binding site of miRNA156 in the CDS1 region (Figure 4D).

The overexpression vector of MrSPL4 was constructed and transformed into E. coli competent cells DH-5α before identifying the transformed bacterial solution by PCR (Figure 4E). The results showed that bands 1, 3, 5, 6, 7 and 8 were consistent with the expected target fragments. In fact, preliminary results further showed that the MrSPL4 gene was successfully inserted into the vector to yield six positive transformant colonies. Two of these positive colonies were randomly selected for sequencing, with the results being still consistent. This experiment therefore showed that the overexpression vector of MrSPL4 gene was successfully constructed and named as 35S::MrSPL4-pCambia2301-KY.

The tobacco leaves infected by Agrobacterium tumefaciens were induced differentiation and budding (Figure 5A), and small seedlings were further grown and induced roots (Figures 5B, C). The regenerated tobacco with Kan resistance was finally transplanted to soil (Figure 5D). Additionally, MrSPL4-F and MrSPL4-R were used as primers, water was set as a blank control, the 35S::MrSPL4-pCambia2301-KY expression vector plasmid was used as a positive control and leaf DNA of WT plants was used as a negative control. The results (Figure 5E) showed that for all resistant regenerated tobacco, the target band was amplified, with its size being consistent with that of the positive control. In addition, no specific bands were observed in WT tobacco and the blank control. Hence, the results showed that the MrSPL4 gene was successfully transferred into tobacco.

Three positive tobacco lines (35S::MrSPL4) from the T₁ generation were randomly selected to detect the expression level of the MrSPL4 gene. Results showed that the relative expression level was significantly higher than that of WT tobacco, with the up-regulation multiple being between 3,862.0-5,938.4 (Figure 6A). This was a clear indication that the gene was overexpressed in transformed tobacco.

The plant heights for 35S::MrSPL4-transformed tobacco and WT ones were measured on 34 d, 41 d and 78 d after transplanting and found to be significantly higher for the transformed plants compared with the WT (Figure 6B). Furthermore, the growth rate was also significantly faster than for the WT. In terms of the budding stages, those of 35S::MrSPL4-transformed tobacco and WT plants were 34 days and 46 days after transplanting, respectively. Based on the above, it could be concluded that 35S::MrSPL4-transformed tobacco showed characteristics of rapid plant growth and early flowering (Figure 6C). Therefore, it was speculated that MrSPL4 gene affected the phenotype of transgenic tobacco to promote plant growth and flowering.