Integrated transcriptome and proteome analysis reveals brassinosteroid-mediated regulation of cambium initiation and patterning in woody stem.

Abstract Wood formation involves sequential developmental events requiring the coordination of multiple hormones. Brassinosteroids (BRs) play a key role in wood development, but little is known about the cellular and molecular processes that underlie wood formation in tree species. Here, we generated transgenic poplar lines with edited PdBRI1 genes, which are orthologs of Arabidopsis vascular-enriched BR receptors, and showed how inhibition of BR signaling influences wood development at the mRNA and/or proteome level. Six Populus PdBRI1 genes formed three gene pairs, each of which was highly expressed in basal stems. Simultaneous mutation of PdBRI1–1, −2, −3 and − 6, which are orthologs of the Arabidopsis vascular-enriched BR receptors BRI1, BRL1 and BRL3, resulted in severe growth defects. In particular, the stems of these mutant lines displayed a discontinuous cambial ring and patterning defects in derived secondary vascular tissues. Abnormal cambial formation within the cortical parenchyma was also observed in the stems of pdbri1–1;2;3;6. Transgenic poplar plants expressing edited versions of PdBRI1–1 or PdBRI1–1;2;6 exhibited phenotypic alterations in stem development at 4.5 months of growth, indicating that there is functional redundancy among these PdBRI1 genes. Integrated analysis of the transcriptome and proteome of pdbri1–1;2;3;6 stems revealed differential expression of a number of genes/proteins associated with wood development and hormones. Concordant (16%) and discordant (84%) regulation of mRNA and protein expression, including wood-associated mRNA/protein expression, was found in pdbri1–1;2;3;6 stems. This study found a dual role of BRs in procambial cell division and xylem differentiation and provides insights into the multiple layers of gene regulation that contribute to wood formation in Populus.


Introduction
Wood consists essentially of cellulose, hemicellulose and lignin and is an important natural and renewable resource for the timber, fiber and bioenergy industries. Wood (secondary xylem) formation involves a complex developmental process starting from the initiation of a vascular cambium, which is secondary meristematic tissue in stems [1,2]. The proliferation and differentiation of the cambial derivative cells generate two separate multicellular tissues: xylem and phloem. Understanding wood developmental programs, particularly (pro)cambium formation, is needed to improve wood biomass production.
Brassinosteroids (BRs), which are plant steroid hormones, play an important role in wood formation. Exogenous application of BR reduces lignification and alters cell wall carbohydrate biosynthesis in the secondary xylem of Liriodendron tulipifera [3]. Transgenic poplar plants overexpressing the BR biosynthesis-related genes DWF4 and BR6OX or the basic helix-loop-helix transcription factor BEE3 (BRASSINOSTEROID ENHANCED EXPRESSION 3) result in increased secondary xylem formation [4][5][6], while inhibition of BR synthesis by propiconazole treatment results in decreased growth and secondary vascular differentiation in poplar [6]. However, key questions, such as which aspects of wood formation are dependent on BRs in tree species, remain unanswered.
Primary mechanisms for the biosynthesis, signal transduction and response of BRs have been well characterized in Arabidopsis [7][8][9]. BR biosynthesis is mediated by several cytochrome p450 monooxygenases, including CPD, DWF4 and BR6OX2 types. BRs are perceived by the cell surface receptor kinase BRASSINOSTEROID-INSENSITIVE 1 (BRI1) and its coreceptor BRI1-ASSOCIATED RECEPTOR KINASE 1 (BAK1), which activate an interconnected signal transduction cascade, resulting in the transcriptional regulation of BRresponsive genes. Loss-of-function BR mutants such as cpds and bri1s display pleiotropic phenotypes, including severe dwarfism, dark-green leaves, and abnormal patterning of stem vascular bundles [10][11][12]. BR signaling coupled to auxin maxima controlled by polar transport is required to establish periodic patterns of stem vascular bundles [12]. BRs facilitate root xylem differentiation by inducing the action of the BIN2-BES1 regulatory module, which is a crucial downstream component of the TDIF signaling pathway [13]. Currently, the mechanisms of BR action during wood development in tree species remain unclear.
Transcriptome analysis is a widely used approach for elucidating the mechanisms underlying the regulation of genes involved in wood formation. For certain processes, the expression of mRNA and protein is concordant. However, mRNA expression does not always reflect the final outcome of protein expression; this is defined as discordant expression hereafter. Thus, when only transcriptome analysis is used, some important genetic information can be missed. In contrast, the integration of multiomic data is an effective method for clarifying developmental processes. For instance, recent advances in integrated multiomic analysis of lignin biosynthesis have led to an in-depth understanding of wood properties in Populus [14]. To date, little is known about the multidimensional regulation of wood development by endogenous BRs in tree species.
Arabidopsis BRI1 has three homologs. BRI1-like 1, 2 and 3 (BRLs); only BRL1 and BRL3 behave as BR receptors, and both of these proteins function redundantly in vascular differentiation with BRI1 [8,15]. In this study, we report the results of experiments showing the dual roles of BR signaling in early xylem development in Populus. Genomic editing of PdBRI1, 2, 3, and 6, the orthologs of three Arabidopsis BR receptors, impairs BR signaling and inhibits the initiation of the vascular cambium, leading to patterning defects in derived secondary vascular tissues in Populus stems. Integrated transcriptome and proteomic analysis of pri1-1;2;3;6 stems revealed concordant and discordant regulation of wood-associated gene and protein expression by BRs at the mRNA and protein levels. Our results demonstrate the similar but distinct roles of BRs in xylem development between Arabidopsis and poplar and could lead to a greater understanding of the molecular mechanisms underlying wood formation in trees.

Identification of six BRI1/BRL orthologs in Populus
Six PtrBRI1 proteins were identified in the Populus genome through sequence similarity searches using Arabidopsis BRI1, BRL1, BRL2 and BRL3 as bait samples. These six PtrBRI1s were named PtrBRI1-1 to PtrBRI1-6 based on their position on Populus chromosomes. Phylogenetic analysis of the ten Arabidopsis and Populus BRI1/BRL proteins yielded a well-defined tree with three clades, which allowed us to identify three pairs of paralogous genes at the terminal nodes (Fig. 1a). PtrBRI1-2 and PtrBRI1-3 are the closest homologs of Arabidopsis BRI1. Synteny analysis showed that two paralogous members of PtrBRI1-1/PtrBRI1-6 or PtrBRI1-4/PtrBRI1-5 may have originated from a chromosome duplication event (Fig. 1b). The tissue expression patterns of the six PdBRI1 genes were subsequently investigated in the poplar genotype Nanlin 895 using qRT-PCR. As shown in Fig. 1c, these PdBRI1s were highly expressed in the shoots and basal regions of poplar stems during secondary growth.  (Fig. 1a), to determine whether PdBRI1s affect vascular development. The expression pattern of PdBRI1-1, −2, −3 and − 6 in the stems was first analyzed using in situ hybridization. As shown in Supplemental Fig. S1, these four genes were mainly expressed in wood-forming cells, including cambial zone, primary xylem parenchyma, and developing secondary xylem/phloem cells. No signals were observed in wood cells of the control sections hybridized with sense probes.

Simultaneous mutation of PdBRI1-1, PdBRI1-2, PdBRI1-3 and PdBRI1-6 results in abnormal cambium formation and alterations in the patterning of secondary vascular tissues
To investigate how mutation of PdBRI1(s) affects wood development, we analyzed the vascular phenotype of stem sections from 45-day-old control and PdBRI1-edited lines. As shown in Fig. 3a, the stems of the control plants had a continuous cambial layer that produced secondary xylem toward the center of the stem and secondary phloem toward the outside of the stem. Lignified phloem fibers differentiated into groups at the periphery of the phloem. No significant difference in the overall vascular morphology was observed in the pdbri1-1 and pdbri1-1;2;6 plants compared with the controls (Fig. 3b, c). Occasionally, the continuity of the cambial zone was interrupted in the pdbri1-1;2;6 stems ( Fig. 3c). In contrast, stem sections from the pdbri1-2;3;6 and pdbri1-1;2;3;6 plants revealed significant defects in cambium initiation and patterning of derived secondary vascular tissues (Fig. 3d, e). The cambial activity of these mutant lines showed a periodically disjointed distribution in the stems (Fig. 3d, e). The xylem and phloem cells derived from the cambium were absent in the discontinuous cambial region (Fig. 3e, f). Abnormal cambial activity in the cortex gave rise to secondary xylem toward the epidermis (Fig. 3e, g). Examination of wood cell morphology revealed a reduction in the cell size and cell wall thickness of xylary fibers and vessels in the pdbri1-1;2;3;6 plants compared with the controls (Fig. 3h, i). The number of cambial cell layers in pdbri1-1;2;3;6 was identical to that in the controls (Fig. 3h, i), implying that inhibition of BR signaling in poplar may not influence cambial division during radial stem growth.
Statistical analysis showed that the radial widths of phloem and xylem were reduced by >34% and 40%, respectively, in the stems of the pdbri1-2;3;6 and pdbri1-1;2;3;6 plants compared with the controls after 45 days of growth ( Fig. 3j). At this stage, there were slight reductions in phloem width (< 14%) and xylem width (< 8%) in the pdbri1-1 and pdbri1-1;2;6 stems. Consistent with these findings, the ratio of phloem width to xylem width was markedly increased in the pdbri1-2;3;6 and pdbri1-1;2;3;6 plants compared with the controls, but it was not significantly different among the pdbri1-1, pdbri1-1;2;6 and control plants (Fig. 3k). Taken together, these results indicate that mutation of PdBRI1s alters the differentiation of xylem and phloem tissues from procambium cells in poplar stems.

PdBRI1-1;2;6
The growth status of PdBRI1-edited lines growing in soil was further investigated. The pdbri1-2;3;6 and pdbri1-1;2;3;6 plants were not able to survive due to severe physiological defects. The pdbri1-1 plants were similar in height to the controls, and they were > 20% taller than the pdbri1-1;2;6 plants after 4.5 months of growth (Supplemental Fig. S3a, b). However, the stem diameter of pdbri1-1, similar to that of pdbri1-1;2;6, was >18% thinner than that of the control plants at this stage (Supplemental Fig. S3c). The radial width of the xylem was reduced by at least 29% in pdbri1-1 and pdbri1-1;2;6 compared to the control plants (Supplemental Fig. S2d). No significant difference in phloem radial width was observed among the three genotypes. Thus, a larger ratio of phloem width to xylem width occurred in the stems of both the pdbri1-1;2;6 and pdbri1-1 plants relative to the controls after 4.5 months of growth (Supplemental Fig. S3e), suggesting that PdBRI1-1 and PdBRI1-6 may promote xylem differentiation, similar to the case for their Arabidopsis orthologs BRL1 and BRL3 [15].

Differential expression of genes and proteins in the basal stems of PdBRI1-1;2;3;6-edited plants compared with the controls.
To explore BR-mediated regulatory networks, we profiled the transcriptome and proteome of the control and pdbri1-1;2;3;6 stems using RNA sequencing (RNA-seq) and liquid chromatography-mass spectrometry (LC-MS), respectively. Hierarchical cluster analysis of the DEGs and DEPs in the two genotypes indicated good repeatability among three biological samples for each genotype (Supplemental Fig. S4). Interestingly, discordant regulation of mRNA and protein expression was observed in the stems of pdbri1-1;2;3;6 plants. In total, 35 793 genes and 6654 proteins were detected in Populus stems after stringent quality checks and data analysis. Of these genes and proteins, 5556 (16%) differentially expressed genes (DEGs; filter criteria: |log 2 [FC]| > 1, padj < 0.05) and 1229 (18%) differentially expressed proteins (DEPs, |log 2 [FC] | > 0.26, p value < 0.05) were identified in the pdbri1-1;2;3;6 lines compared with the controls (Supplemental Fig. S5a, b; Table S1). The ratio of upregulated genes to all DEGs was 73%, which was much higher than that (27%) of downregulated genes to DEGs (Supplemental Fig. S5a). This finding implies that the upregulation of specific mRNAs may be a key response to BR signaling inhibition in poplar stems. In contrast, similar ratios of upregulated proteins (48%) and downregulated proteins (52%) were observed for all DEPs (Supplemental Fig. S5b).
Gene Ontology (GO) analysis (p < 0.05) of the 5556 DEGs and 1229 DEPs revealed strong bias involving representation within the biological process category at both the transcriptome and proteome levels (Supplemental Fig. S5c-f; Table S2). Genes whose expression was induced in pdbri1-1;2;3;6 plants were enriched in protein phosphorylation, transcriptional regulation, and organ morphogenesis (Supplemental Fig. S5c). The downregulated genes mainly participated in cell wall polysaccharide biosynthesis and biogenesis (Supplemental Fig. S5d), which was highly correlated with the inhibited wood development observed in the pdbri1-1;2;3;6 stems ( Fig. 3). According to the proteomic data, the upregulated proteins were enriched in nucleic acid metabolic processes and chromatin organization (Supplemental Fig. S5e). Most of the downregulated proteins were involved in energy reserve-related metabolism, polysaccharide metabolism and nucleotide metabolism (Supplemental Fig. S5f). Our multiomic dataset revealed BR-mediated regulation of various biological processes at the mRNA and protein levels in poplar stems.

Comparative analysis of the transcriptome and proteome of PdBRI1-1;2;3;6-edited plants.
To determine the regulatory patterns of DEGs and DEPs in the pdbri1-1;2;3;6 stems, 5556 DEGs and 1229 DEPs were selected for correlation tests. The expression changes of 1502 mRNAs and proteins were moderately correlated, with R (Spearman) = 0.3721-0.7984 ( Fig. 4a; Supplemental Table S3). Of the 1502 DEGs/DEPs, 234 (16%)  Table  S3. b, c GO terms (p < 0.05) enriched in the group of proteins that were significantly changed at both the mRNA and protein levels or proteins that were significantly changed at the mRNA or protein level.
showed changes in expression at both the transcriptome and proteome levels, and 1268 (84%) showed discordant expression. A total of 964 proteins (64%) were differentially expressed at the protein level but not at the mRNA level. Similarly, a total of 290 genes (20%) showed altered expression in an mRNA-specific manner, with 197 genes upregulated and 93 genes downregulated. GO analysis (p < 0.05) of the 1502 DEGs/DEPs showed various enrichments in biological processes, which were more significant at the mRNA level than at the protein level (Fig. 4b). These findings suggest that inhibition of BR signaling in poplar may cause larger changes at the transcriptional level. We further investigated the top 3 enriched biological processes among these DEGs and DEPs. As shown in Fig. 4c, the DEGs/DEPs at both the mRNA and protein levels were enriched in processes such as oxidation-reduction, carbohydrate metabolism, and cell wall organization or biogenesis. The genes that were differentially expressed specifically at the transcriptional level were enriched in protein phosphorylation, transcriptional regulation and chemical homeostasis. In contrast, DEPs specifically identified at the protein level were mainly involved in carbohydrate catabolism, nitrogen compound regulation, and DNA packaging. Collectively, these results reflected the multifaceted control of gene expression in the stems of pdbri1-1;2;3;6.  Table S3). The Arabidopsis class II KNOX gene KNAT7 negatively regulates secondary cell wall thickening in stems and is functionally conserved in Populus [16]. A Populus ortholog of KNAT7 was identified to have mRNA and protein expression levels that were higher in pdbri1-1;2;3;6 than in the controls (Fig. 5a), which was consistent with the thinner xylary fiber and vessel cells in these mutant lines compared with the controls (Fig. 3h, i). BXL4/XYL4 possesses β-D-xylosidase activity, which is responsible for cleavage of the xylan backbone in stem tissues of Arabidopsis [17]. Higher mRNA and protein expression levels were observed for a Populus ortholog of BXL4 in pdbri1-1;2;3;6 ( Fig. 5a). PME catalyzes the demethylesterification of pectin, a key component of polysaccharides in cell walls, and some PMEs are important for developmental patterns of cell expansion in developing wood cells of Populus [18]. We found that 4 of 5 Populus PMEs showed reduced mRNA and protein expression levels in pdbri1-1;2;3;6 ( Fig. 5a). Orthologs of CCoAOMT, F5H1, C3H1, CAD4, and COMT2 in Populus are mainly expressed in maturing xylem and are required for the synthesis of lignin in the stem [14,19]. Their mRNA and protein expression levels were lower in pdbri1-1;2;3;6 than in the controls (Fig. 5a). The auxin efflux carriers PINFORMED 1 (PIN1) and PIN3 function during the establishment and activity of the interfascicular cambium in Arabidopsis [12,20]. Orthologs of PIN1 and PIN3 in Populus are mainly expressed during secondary vascular growth [21] and exhibited increased mRNA and protein expression levels in pdbri1-1;2;3;6 ( Fig. 5a). Nine of these DEGs/DEPs were further selected to validate the transcriptomic and proteomic data. The qRT-PCR results revealed that the transcripts of the KNAT7 and BXL4 orthologs increased, and those of the PME3, C3H1, CAD4, COMT2, PIN1, and PIN3 orthologs decreased in pdbri1-1;2;3;6 stems compared with the control stems (Fig. 5b). Moreover, parallel reaction monitoring (PRM) data showed that the protein enrichment of two PME3 orthologs was significantly lower in the pdbri1-1;2;3;6 stems than in the controls (Fig. 5c). These results verified the accuracy of our omic approaches. The other seven Populus orthologs could not be detected using PRM technology due to their low amount of expression in the stems.
Multiple wood-associated genes were differentially expressed at the transcriptional level but not at the protein level in pdbri1-1;2;3;6 stems relative to the controls ( Fig. 6a; Supplemental Table S3). These included 28 LACCASE genes that were mainly expressed in stemdifferentiating xylem or roots and may be involved in lignin polymerization [22]. Twelve Populus orthologs of Arabidopsis LAC4 and LAC17 showed reduced expression in pdbri1-1;2;3;6 stems. LAC4 and LAC17 are necessary for lignin deposition during stem vascular development [23,24]. A number of DEGs associated with cytokinins, BRs and ethylene were also identified in pdbri1-1;2;3;6 stems (Fig. 6a). These three hormones are known to regulate cell division and differentiation in vascular tissues. In Arabidopsis, CYP735A1 and ADENOSINE PHOSPHATE-ISOPENTENYLTRANSFERASE 3 (IPT3) are essential for the biosynthesis of trans-zeatin, an isoprenoid cytokinin [25,26], while CYTOKININ OXIDASE/DEHYDROGENASEs (CKXs) catalyze the irreversible degradation of cytokinins [27]. Orthologs of IPT3 and CKXs in Populus showed lower and higher transcript levels, respectively, in pdbri1-1;2;3;6 than in the controls (Fig. 6a). Moreover, Populus orthologs of Arabidopsis LONELY GUYs (LOGs), which are responsive to the direct activation pathway of cytokinins [28], were induced in pdbri1-1;2;3;6 stems. Increased expression levels of Populus orthologs of four BR biosynthesis-related genes (DWF4, ROT3, CYP90D1 and BR6OX2) in pdbri1-1;2;3;6 stems indicated feedback suppression of BR biosynthesis-related genes in the BR signaling pathway of Populus, similar to that in Arabidopsis [7]. In addition, several xylem-specific ETHYLENE RESPONSE FACTOR (ERF) genes [29], including orthologs of ERF1, ERF18 and ERF109, which are required for cell division in the cambium of Arabidopsis stems [30], showed higher expression levels in the pdbri1-1;2;3;6 stems than in the control stems. Sixteen of these DEGs were validated by qRT-PCR analysis (Fig. 6b), the data of which were in agreement with the RNA-seq data.

BRs mediate the regulation of cambium initiation and xylem differentiation in poplar stems
Stem radial growth of trees requires the proliferation of procambial and cambial cells and the continuous accumulation of their differentiation products (xylem and phloem) [1]. These two processes are finely regulated by several hormones, including auxin [31], cytokinin [32] and ethylene [33]. A recent study in Populus showed that exogenous application of propiconazole inhibits BR synthesis and reduces xylem differentiation and secondary growth [6]. However, it remains unknown whether other aspects of wood development in tree species require BRs. This study aims to explore the regulatory mechanism through which BRs affect wood development. We inhibited BR signaling in poplar by CRISPR-Cas9-mediated knockout of PdBRI1-1, −2, −3 and − 6. These four genes are orthologs of Arabidopsis BRI1, BRL1 and BRL3, which encode vasculature-enriched members of the BR receptor family [8]. Transgenic poplar lines expressing edited versions of PdBRI1-2;3;6 or PdBRI1-1;2;3;6 showed severe BR-deficient phenotypes, including periodic rupturing of the cambial ring in the stems (Figs. 2, 3), which is consistent with the reduced numbers of stem vascular bundles in Arabidopsis mutants with reduced BR synthesis or signaling, including cpd, bri1-116, bri1-301, bin2 [12] and bri1 brl1 brl3 mutants (Supplemental Fig. S6). The procambium arises from oriented and coordinated cell divisions and creates (inter)fascicular cambia, forming a cylinder of cambium in Arabidopsis and woody stems [2]. Therefore, it is possible that BRs promote early procambial divisions in poplar stems, similar to those observed in Arabidopsis stems [12]. In addition to discontinuous cambial segments at normal positions, abnormal initiation of the cambia in the cortical parenchyma was observed in the stems of pdbri1-2;3;6 and pdbri1-1;2;3;6 (Fig. 3). These morphological defects were not found in the stems of bri1 brl1 brl3, a null Arabidopsis mutant which mutations in three BR receptors [15] (Supplemental Fig. S6). This finding suggests that there are similar but distinct roles Extensive studies of Arabidopsis mutants have demonstrated that BR signaling coupled with auxin polar transport is required for cambium initiation in stems [12]. As part of this pathway, MP/ARF5 and its upstream positive regulators, PINs, play major roles in translating auxin accumulation into the establishment of procambium identity [34]. mp/arf5 is activated in preprocambial strands and exhibits severe defects in procambium formation [35]. Compared with WT controls, pin1 pin2 double mutants display a disorganized vascular pattern with an increased number of vascular bundles and xylemdifferentiated cells in stems [12]. Our transcriptome data revealed higher expression levels of the Populus orthologs of MP/ARF5, PIN1, and PIN3 in pdbri1-1;2;3;6 than in the controls ( Fig. 5; Supplemental Table S1). Taken together, these findings suggest that BR-induced vascular cambium initiation might involve MP-PIN cascade-mediated auxin maxima in poplar stems.
In addition to the DEGs involved in auxin, a number of cytokinin-associated DEGs were identified in the pdbri1-1;2;3;6 stems compared with the controls (Fig. 6). Two canonical cytokinin synthesis-related genes, the Populus orthologs of CYP735A1 and IPT3, showed reduced expression, whereas three cytokinin degradation-related genes, the Populus orthologs of CKX1, CKX3 and CKX5, showed increased expression in pdbri1-1;2;3;6. Since reduced cytokinin signaling leads to altered subcellular polarity of PIN1, PIN3 and PIN7 and a loss of PIN7 expression in the procambium of Arabidopsis roots [36], it is speculated that BRs promote procambium formation possibly by mediating cytokinin signaling in poplar stems. Interestingly, inhibition of BR signaling in poplar induced the expression of five Populus orthologs of LOGs that are responsive to the pathway through which cytokinins are directly activated [28]. These observations suggest that a complex regulatory network involving BRs, auxin, and cytokinins occurs during procambium formation.
Another striking aspect of pdbri1-2;3;6 and pdbri1-1;2;3;6 plants is the changes in a balanced differentiation between secondary xylem and phloem tissues in the stems (Fig. 3). Consistently, multiple DEGs involved in cell differentiation were identified in the stems of pdbri1-1;2;3;6 plants (Supplemental Fig. S5). Our results, together with those of others concerning various systems, including Zinnia [37], Arabidopsis [13,15] and Populus [6], indicate that BRs play a key role in the regulation of xylem differentiation. An obvious defect in stem xylem differentiation was observed in the pdbri1-1;2;6 and pdbri1-1 plants at 4.5 months of growth but not at 45 days of growth ( Fig. 3; Supplemental Fig. S3), suggesting that there is potential redundancy among and a stage preference of PdBRI1-1, PdBRI1-2, PdBRI1-3 and PdBRI1-6 in regulating xylem differentiation during secondary growth. Figure 6. Transcriptional regulation of wood/hormone-associated DEGs in pdbri1-1;2;3;6 stems. a A battery of wood-associated genes that were differentially expressed in pdbri1-1;2;3;6 at the mRNA level but not at the protein level. b Sixteen DEGs marked with red dots in (a) were selected for qRT-PCR analysis in the CK and pdbri1-1;2;3;6 stems.

Concordant and discordant regulation of wood-associated gene and protein expression at the mRNA and protein levels when BR signaling is impaired in poplar
To determine exactly how endogenous BRs influence wood development, we compared the transcriptome and proteome profiles of control and pdbri1-1;2;3;6 stems. Only 16% of mRNA and protein expression was correlated, the results of which are similar to those of numerous studies in which a poor correlation between transcript and protein abundance was observed [38][39][40]. It is possible that the role of BRs in biological processes involves various posttranscriptional regulatory events. Of the 5556 DEGs from pdbri1-1;2;3;6 stems, those involved in cell wall biogenesis and metabolism (e.g. lignin and xylan metabolic processes) constituted an overwhelming majority of the downregulated genes (Supplement Fig.  S5d), consistent with the inhibition of xylem development in pdbri1-1;2;3;6 (Fig. 3). The 3 most enriched GO terms of the DEGs/DEPs at both the mRNA and protein levels included "cell wall organization or biogenesis" (Fig. 4). Among them, the Populus orthologs of CCoAMT, F5H1, C3H1, CAD4 and COMT2, which correspond to five key enzymes in the Arabidopsis lignin biosynthesis pathway, showed reduced mRNA and protein expression in pdbri1-1;2;3;6 stems ( Fig. 5). Interestingly, twelve Populus orthologs of LAC4 and LAC17 that are responsive to lignin polymerization in Arabidopsis stems [22] showed reduced expression in pdbri1-1;2;3;6 in a mRNA-specific manner (Fig. 6). Thus, BRs may regulate lignin-related gene expression via multiple pathways, promoting lignin biosynthesis in poplar stems.
Receptor kinases mediate BR signal transduction mainly via a series of reversible phosphorylation and dephosphorylation processes in plants [7]. In our study, a large overrepresentation of genes involved in protein phosphorylation showed altered expression at the mRNA level but not at the protein level in pdbri1-1;2;3;6 stems compared with the controls (Supplemental Fig. S5c). This finding suggests that inhibition of BR signaling in poplar induces BR deficiency-related responses, likely involving transcriptional regulation of protein phosphorylation. Overall, this study reveals the dual roles of BRs in early xylem development and provides molecular insights into the multiple layers of gene regulation by BRs that contribute to wood development in Populus.

Plant materials and growth conditions
The poplar (Populus deltoides × Populus euramericana) genotype Nanlin 895 was used for gene expression analysis and genetic transformation. The poplar plants were grown in a greenhouse at 22-24 • C under a photoperiod comprising 16 h of light/8 h of darkness.

qRT-PCR analysis
Total RNA isolation and first-strand cDNA synthesis were conducted following previously described methods [42]. qRT-PCR was performed on an Applied Biosystems ® QuantStudio ® 1 System (Thermo Fisher Scientific, U.S.A.) in conjunction with a SYBR Green Master Mix Kit (Trans-Gen Biotech, China) according to the manufacturers' protocols. The 2 -CT method [43] was used to calculate the relative expression levels of PdBRI1s, with PdUBQ10 used as an internal control. The specificity of the primers (Supplement Table S4) was validated by melting profiles showing a single product-specific melting temperature.

In situ hybridization
In situ hybridization was conducted following a previously described method [42]. Gene-specific fragments from PdBRI1-1, −2, −3 and − 6 were used to synthesize digoxygenin-labeled antisense and sense RNA probes (the primers are listed in Supplemental Table S4) with DIG RNA labeling mixture (Cat.: 11175025910; Roche). The basal stems of 1-month-old poplar (Nanlin 895) plants were sectioned, fixed in 4% glutaraldehyde, and embedded in paraffin. Stem sections (10 μm) were cut with a rotary Leica RM 2235 microtome (Leica, Nussloch, Germany), mounted onto Superfrost Plus glass slides (Thermo Fisher, Wilmington, USA), and hybridized with antisense or sense RNA probes.

Microscopy and histochemistry
Basal stems of 45-day-old or 4.5-month-old poplar plants were fixed in 4% paraformaldehyde at 4 • C for 3 d and then embedded in Paraplast following a previously described method [42]. The Paraplast-embedded stems were sliced into 7-μm-thick sections using a Leica RM 2235 microtome (Leica) and adhered to Superfrost Plus microscope slides (Thermo Fisher). After staining with 0.1% (w/v) toluidine blue-O (TBO), the sections were observed and imaged using an Olympus X51 light microscope. The width of the xylem and phloem in the stems was measured using Olympus DP2-BSW software. At least three independent replicates were included for each genotype. Differences between the means of the control and transgenic lines were evaluated by Duncan's test at the 0.05 and 0.01 probability levels.
For observation of cambium development in CK and PdBRI1-1;2;3;6-edited lines, 0.5-cm stem segments of 45day-old seedlings were fixed in 1% osmic acid and then embedded in resin (SPI-Chem, Westchester, USA) following our previous method [42]. One-micrometer sections were cut with a glass knife via a Leica EM UC6 microtome and stained with TBO to observe the cambial layer.

RNA-seq analysis
Three representative PdBRI1-1;2;3;6-edited lines (BE4, BE7, and BE8) and controls were grown in bottles for 45 days. After removing the leaves, basal stems (2 cm) were sampled for RNA isolation. Thirty plants from each line were pooled as a single sample, and three biological replicates were included. RNA purification, library construction and sequencing were conducted by Shanghai Applied Protein Technology Co., Ltd. (Shanghai, China). cDNA libraries were sequenced using an Illumina HiSeq 4000 system, with 150-bp paired-end reads. Sequences for which Q was ≥20 were identified and extracted via custom Perl scripts. The filtered reads were mapped to the poplar genome (v.3.1) using HISAT2 software [45]. Transcript quantification was performed in terms of fragments per kilobase of transcript per million mapped reads (FPKM). Genes with an adjusted p value (padj) < 0.05 and log 2 [Fold Change] > 1 were defined as differentially expressed. The expression data were genewise normalized and hierarchically clustered based on Pearson coefficients with average linkage in Genesis (version 1.8) [46]. The raw data are available in the Gene Expression Omnibus database (accession No. GSE179471).

Protein extraction, TMT labeling and mass spectrometry
Total proteins were extracted from the samples that were used for RNA-seq analysis using the tricarboxylic acid (TCA)/acetone precipitation and SDT cracking method [47]. Protein concentrations were determined using a BCA protein assay kit (Bio-Rad, USA). The extracted proteins were digested into peptides according to the filter-aided sample preparation (FASP) procedure [47], and the resulting peptides were eluted for TMT labeling (Thermo Fisher Scientific). The TMT-labeled samples were then fractionated into 12 fractions with a high-pH reversed-phase fractionation kit (Thermo Fisher Scientific). Each fraction was then injected for nanoliquid chromatography tandem mass spectrometry (nano-LC-MS/MS) analysis. The peptide mixture was loaded onto a reverse-phase trap column (Thermo Scientific Acclaim Pep Map100; 100 μm × 2 cm, nanoViper C18) connected to a reversephase analytical column (Thermo Scientific Easy Column; 10 cm, 75 μm inner diameter, 3 μm resin; C18-A2) in buffer A (0.1% formic acid) and separated with a linear gradient of buffer B (84% acetonitrile and 0.1% formic acid) at a flow rate of 300 nL/min. The MS data were acquired with a Q Exactive mass spectrometer (Thermo Scientific) via a superior data-dependent method (top 10 ranking). The most abundant precursor ions were filtered and removed during the scan (300-1800 m/z), and a maximum injection time of 10 ms was used. The survey scans were acquired at a resolution of 70 000 at m/z 200, and the resolution of the HCD spectra was set to 35 000.

Quantitative proteomic analysis
Peptide recognition mode was selected, and the MS/MS spectra were filtered via the MASCOT engine (Matrix Science; version 2. 2) embedded in Proteome Discoverer 1.4 (Thermo Fisher) against the Ptrichocarpa v4.1 protein database. Protein identification was based on the following parameters: peptide mass tolerance, ± 20 ppm; fragment mass tolerance, 0.1 Da; enzyme, trypsin; maximum missed cleavages, 2; fixed modifications, carbamidomethyl (C), TMT 6/10 plex (N-term), and TMT 6/10 plex (K); variable modifications, oxidation (M) and TMT 6/10 plex (Y); FDR, ≤ 0.01; protein quantification, protein ratios calculated as the medians of only unique peptides of the protein; and experimental bias, all peptide ratios normalized by the median protein ratio, with the median protein ratio being 1 after normalization. Log 2 (FC) > 1.2 or < 0.83 and P value <0.05 (Student's t test) were used as filter criteria to identify differentially expressed proteins (DEPs). The mass spectrometry proteomic data have been deposited in Integrated Proteome Resources (iProX) database under dataset identifier IPX0003236001/PXD027114.

Gene ontology annotations
DEGs or DEPs were separately standardized using zscores, and a GO expression matrix was generated via the mean z-score for each GO term. The DEPs were subjected to GO analysis against the translated Populus transcriptome (FDR-corrected p value< 0.05). The GO process revealed that the DEGs and DEPs are associated with three functional categories: biological processes (BPs), cellular components (CCs), and molecular functions (MFs). The DEGs and DEPs were selected for an integrated analysis, and Pearson correlation coefficients and p values were calculated using the Spearman method [48].

Parallel reaction monitoring (PRM) analysis
To validate the expression levels of two DEPs through TMT-based proteomic analysis, additional quantification was performed through LC-PRM/MS analysis according to a protocol [49] adopted by Shanghai Applied Protein Technology Co., Ltd. The analysis of raw data was performed via Skyline (MacCoss lab, University of Washington) [50].