Transcriptome profiling of Brassica napus stem sections in relation to differences in lignin content

Plant materials and growth conditions

Brown seeded B. napus DH12075 (designated DH), yellow seeded B. napus YN01–429 (designated YN), alfalfa (Medicago sativa L.) and Arabidopsis thaliana ecotypes and mutant lines were grown in soil in a controlled greenhouse environment (16 h light/8 h dark, 20 °C/17 °C). For B. napus stem samples, plants with well-developed siliques were harvested and 40 cm of the stem was selected (excluding the lower most 10 cm). Collected stem samples were divided into 4 segments, 10 cm each, and apical (section 1) and basal sections (section 4) were used for lignin, microarray and qPCR analyses. For Arabidopsis lignin analysis, Arabidopsis plants (42-day-old) were cut at the base of the inflorescence stem and leaves and thin branches were removed, keeping the stem and main/robust branches for analysis.

Lignin analysis

Lignin content of stem material was determined by thioglycolic acid (TGA) method [24]. Fifty mg of stem tissue ground in liquid N2 from 3 different plants per cultivar (i.e., 3 biological replicates) was used for lignin content determination.

DNA microarray

Total RNA from stem sections of both DH and YN were extracted as described in Carpenter and Simon [25] followed by clean-up using the commercial RNeasy mini kit (Qiagen, Valencia, CA, USA). Two biological replicates were collected for each sample. RNA amplification, labeling with cy3- or cy5-dCTP dyes (GE Healthcare, Buckinghamshire, UK), and probe fragmentation was carried out using Ambion AminoAllyl MessageAmp II RNA amplification kit according to the manufacturer’s instructions (Ambion, Austin, TX, USA). A spotted 15 K 50-mer B. napus oligo array, previously developed at the Saskatoon Research and Development Centre, Agriculture and Agri-Food Canada [26], was hybridized with the cy5- and cy3-labelled probe pairs. Cy3 dye was used for labeling the YN RNA and Cy5 for DH RNA. A dye swap (cy3/cy5) labelling experiment was performed for each biological replicate. Labeling, hybridization and post- hybridization washing were conducted according to the protocol for Corning epoxide slides (Corning Inc., Lowell, MA, USA). Hybridization was carried out at 37 °C for 17 h in a MAUI hybridization station (BioMicro Systems, Salt Lake City, UT, USA). Following the post-hybridization washes, slides were scanned with a VersArray ChipReader laser scanner (Bio-Rad Laboratories, Hercules, CA, USA).

Image files were transferred to ArrayPro Analyzer software (Media Cybernetics Inc., Bethesda, MD, USA) for image analysis, spot verification, normalization, filtering and feature extractions. A standard statistical program GeneSpring GX (Agilent Technologies, Santa Clara, CA, USA) was used to determine that the data set would fit a normal distribution (log 10) and to check the validity of expression data by assessing the variation between control spots. Data files for the duplicates on individual slides and the dye swap data files for each experiment were merged and average spot intensities used to reduce experimental bias. Normalization was performed on merged data using the Lowess (sub-grid) method, and the local background was subtracted from the values of each spot. The intensity of each spot at λ549 nm (Cy5) and λ647 nm (Cy3) was analyzed using a BASE plug-in and finally transformed into a DH/YN ratio value that denoted the most upper (section 1) and lowest stem sections (section 4), including DH1:YN1, DH4:YN4, DH1:DH4 and YN1:YN4. Additional data processing was performed using tools available in BASE ( http://base.thep.lu.se ). Gene ontology (GO) analysis was conducted using the TAIR database [27]. Microarray ratios with transcript changes > 2-fold upregulated or downregulated between two different B. napus stem sections were considered for all ratios where p ≤ 0.05 (Additional files 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11 and 12: Table S1A and B to S6A and B). Data in these 12 Additional files were screened using either Excel or Microsoft Access to identify individual or groups of B. napus genes expressed in multiple DH and YN ratios in the main categories of cell wall, other carbohydrate genes and transcription factors, and in additional subcategories. Ten major transcription factor families, Arabidopsis homologue IDs and corresponding mutant lines for selected B. napus genes were retrieved from http://www.arabidopsis.org/browse/genefamily/index.jsp [28].

Quantitative real-time RT-PCR

For quantitative real time reverse transcription PCR (qRT-PCR) experiments, total RNA was extracted from stem sections using a TRIzol reagent ( http://www.invitrogen.com ) following the manufacturer’s instructions. Three biological replicates (independent RNA preparations) were collected and used with three technical reps for qRT-PCR analysis. RNA concentration was determined using a nanodrop spectrophotometer (Fisher Scientific, Canada). Twenty micrograms of total RNA (per sample in each replicate) was treated with Turbo DNase ( www.ambion.com ) as per the manufacturer’s instructions to eliminate trace amounts of genomic DNA. Reverse transcription reactions were performed with Superscript™ III Reverse Transcriptase (http:// www.invitrogen.com ) according to the manufacturer’s instructions using 1.5 μg RNA per reaction. Two RT reactions per sample were done and reactions were pooled after 10-fold dilution with nuclease-free water. Polymerase chain reactions were carried out using 96-well plates in a LightCycler® 480 II ( http://www.roche-applied-science.com/lightcycler-online ) using SYBR® Green. Reactions contained 10 μl of 2× SYBR® Green Master Mix (Roche), 5 μl of diluted cDNA and 200 nM of gene/EST-specific primer (Additional file 13: Table S7) in a final volume of 20 μl. For data normalization, five reference genes (actin, adenine phosphorybosyl transferase, β-tubulin, cyclophilin and elongation factor 1-α) were included in the experiment and two genes with stable expression, namely β-tubulin and actin, were selected using geNorm software [29]. PCR efficiency, in the range of 85 to 100%, was determined from amplification plots using the program LinRegPCR [30].

Plant transformation

RNAi and over-expression constructs were developed using either the closest Arabidopsis homologue according to Hossain et al. [31, 32] or the closest alfalfa stem cDNA according to Li et al. [33]. Alfalfa primers were identified by Laberge (personal communication). Floral bud transformation was conducted in Arabidopsis according to Hossain et al. [31, 32]. Alfalfa transformation was conducted using Agrobacterium tumefaciens according to Aung et al. [34].

DH and YN stems differ in lignin and exhibit different gene expression profiles

To validate our introductory hypothesis, we first determined the stem lignin content of DH and YN plants. Total lignin content showed higher levels in both apical (section 1) (DH1) and basal sections (section 4) (DH4) of DH compared to their counterparts in YN cultivar sections (YN1 and YN4) (Fig. 1a). Within stem sections of each cultivar, we found the highest lignin content to be in the basal sections and the lowest in the apical sections. Increased lignin accumulation is the result of increased enzyme activity and/or increased gene expression. Hence, we compared global gene expression levels between the two cultivars and also between these sections using a 15 K B. napus microarray that was used earlier to investigate gene expression in developing B. napus seed [35] and diseased stem sections [26]. In detail, we related lignin levels with transcriptome changes between DH and YN (both apical and basal sections for each cultivar (Fig. 1b; Additional files 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11 and 12: Tables S1A and B to S6A and B). The greater upregulation of genes in the apical stem sections compared with basal stem section expression pattern is reasonable because the apical region of the plant is developmentally more active compared to the basal region.

Functional classification of genes up-regulated in high lignin stems of B. napus

Gene annotation of sequences differentially expressed using the 15 K B. napus array was conducted using The Arabidopsis Information Resource ( http://www.arabidopsis.org ). Annotated genes were grouped into three major functional categories: cellular component, molecular function and biological process, and then divided into subcategories (Fig. 2; Table 1). For ease of presentation, data was expressed as the average of genes in each subcategory in DH vs. the average in YN. Within the category of biological processes in DH4 vs. YN4, genes involved in defense responses were 7.6% of the total, whereas stress responsive genes were 7.4%; genes involved in developmental processes were 4.9 and 4.7% were involved in transport; 2.35% of genes were involved in transcription and 0.6% in DNA or RNA metabolism. In the cellular component category, 6.4% genes were related to plasma membrane, and nucleus genes accounted for 6.02%; cell wall genes were 2.12 and 1.6% genes were related to cytoplasm functions. In the molecular function category, 7.4% genes were involved in protein binding, and 5.7% genes accounted for DNA or RNA binding activity; 4.8% genes were involved in transporter activity, and 3.4% genes showed transcription factor activity (Fig. 2; Table 1).

	DH4:YN4		DH1:YN1		DH1:DH4		YN1:YN4
	Up	Down	Up	Down	Up	Down	Up	Down
Cellular Component
Other intracellular components	754	693	668	580	475	396	377	205
Chloroplast	602	567	581	462	413	320	308	159
Other cytoplasmic components	636	469	531	403	361	321	280	161
Other membranes	531	482	466	483	374	290	284	156
Unknown cellular components	607	607	513	373	463	352	347	187
Plastid	291	343	260	303	174	151	246	122
Plasma membrane	438	381	421	326	262	240	133	79
Nucleus	392	227	331	200	247	203	181	107
Other cellular components	284	285	184	153	225	162	163	103
Mitochondria	198	191	210	188	125	94	97	55
Cell wall	127	123	115	88	75	66	64	40
Cytosol	121	102	106	80	66	66	56	31
Extracellular	110	82	82	62	66	62	44	33
ER	90	76	81	60	60	48	44	33
Golgi apparatus	49	43	40	37	31	26	18	12
Ribosome	23	23	23	18	16	13	10	8
Molecular Function
Other enzyme activity	515	695	439	383	492	295	383	216
Unknown molecular functions	650	471	578	515	341	376	259	156
Other binding	525	469	455	393	327	274	248	150
Transferase activity	295	308	243	229	211	186	188	100
Hydrolase activity	341	328	275	299	240	199	151	98
Protein binding	351	260	306	230	223	199	176	115
Nucleotide binding	279	132	231	217	131	96	139	53
Kinase activity	137	245	114	110	175	143	85	88
DNA or RNA binding	257	182	160	183	83	125	121	73
Transporter activity	164	201	197	131	153	93	70	56
Other molecular functions	185	157	143	110	122	93	103	45
Transcription factor activity	165	149	130	115	117	93	79	47
Nucleic acid binding	119	111	98	102	76	64	61	28
Structural molecule activity	35	22	28	25	26	16	20	8
Receptor binding or activity	15	18	17	12	16	10	9	5
Biological Process
Other cellular processes	1254	1145	1066	943	803	710	641	398
Other metabolic processes	1139	1016	974	860	754	664	597	364
Unknown biological processes	763	760	637	575	583	418	444	233
Response to abiotic or biotic stimulus	387	373	350	279	221	197	200	115
Response to stress	381	342	337	263	206	188	181	115
Other biological processes	353	334	325	269	215	197	183	105
Protein metabolism	375	359	333	317	247	229	130	127
Developmental processes	251	261	222	181	197	158	189	83
Transport	266	232	232	193	165	124	139	85
Cell organization and biogenesis	171	154	150	137	127	110	97	56
Transcription	159	146	130	119	104	88	73	42
Electron transport or energy pathways	81	74	61	54	51	31	45	17
Signal transduction	75	52	66	46	39	38	37	22
DNA or RNA metabolism	35	27	24	21	25	20	18	14

NB data was ranked in descending order using the DH4-YN4 pair wise comparison. Up., upregulated. Down, down regulated

Annotation for genes differentially expressed in DH1 vs. YN1, DH1 vs. DH4 and YN1 vs. YN4 pairwise comparisons was also carried out. In all three functional categories, genes belonging to the above three groupings showed a similar pattern of distribution to those of DH4 vs. YN4.

Validation of microarray data by quantitative RT-PCR

To validate the results of microarray experiments, transcript levels of seven genes were examined by qRT-PCR using total RNA prepared from stem sections. Emphasis was placed on functional diversity in selecting the following genes for qRT-PCR experiments; ANAC062 (Accession # ES901271), PHYTOCHOME INTERACTING FACTOR 4 (PIF4; EV181356), GLYCINE-RICH RNA BINDING PROTEIN (RBP/ CAA78513.1), DIMINUTO 1, LIPASE 1 (AY866419), FIBRILLIN (AF084554.1) and GDSL-MOTIF CONTAINING LIPASE/HYDROLASE (CN730052) (Fig. 3a and b; Additional file 13: Table S7). From the microarray data, it was determined that these selected genes were upregulated in the DH cultivar with high lignin content. The expression profiles of most of the genes tested by qRT-PCR were consistent with our microarray data, although the difference in expression was less in the latter data set. Expression of DIMINUTO 1 was significantly higher in the basal section of DH (DH4). Expression of ANAC062 was also higher in the basal section of DH4, whereas PIF4 showed higher expression in both apical (DH1) and basal sections (DH4). Expression of LIPASE1 was also higher in both sections of DH compared to those of YN. Transcripts of RBP exhibited slight numerical increase in both sections of DH, although the difference was not statistically significant. Because of their extremely low expression levels relative to other tested genes, the two genes encoding FIBRILLIN and GDSL-MOTIF PROTEIN were assessed separately (Fig. 3c).

Both of the latter genes showed significantly higher expression in the two sections of DH relative to their counterparts in YN.

Upregulation of monolignol biosynthetic genes

Cinnamoyl-CoA reductase (CCR) and cinnamyl alcohol dehydrogenase (CAD) also play important roles in monolignol biosynthesis. The cinnamoyl-CoA esters, precursors of monolignol biosynthesis, are generated by the general phenylpropanoid pathway and then converted into monolignols by CCR and CAD [7]. CCR1 was among the most upregulated genes in the top and basal sections of DH relative to those of YN, and the same gene was upregulated in the basal section relative to the top section of YN in the array (Table 2). CCR2 was also upregulated in the top section of DH relative to the top section of YN. Another major gene change in the phenylpropanoid pathway was coded by p-coumarate 3-hydroxylase (C3H), also known as REDUCED EPIDERMAL FLUORESCENCE (REF8) was upregulated in the top and basal sections of DH relative to those of YN (Table 2).

The DE-ETIOLATED 3 (DET3), At1g12840, was upregulated in the top and basal section of DH relative to the top and basal section of YN (Table 2). This gene encodes the C-subunit of the V-type ATPase [37] implicated in the maintenance of pH homeostasis in plants [38]. Moreover, the FAH1 gene was upregulated in the basal section of DH4 relative to the basal section of YN4, and the same gene is upregulated in the basal section relative to the top section of both cultivars (Table 2). The basal section is more mature and richer in lignin content relative to the top section (Fig. 1a), suggesting the involvement of this gene with lignin accumulation in Brassica. The ferulate-5-hydroxylase (F5H) enzyme is a key regulatory point in the determination of lignin monomer composition [39] by catalyzing the conversion of ferulic acid to sinapic acid and syringyl lignin.

Differential expression of glucosyltransferase, peroxidase and laccase genes

It is well established that monolignols are synthesized in the cytoplasm and then exported into the cell wall for polymerization [7]. Cellular homeostasis of these monomers is regulated through glucosylation [40]. This modification step increases the solubility and stability of monomers and provides access to the membrane transport systems for transportation and storage. Glucosyltransferases play an important role in these processes [7]. After transport of the monolignols to the cell wall, lignin is formed through dehydrogenative polymerization of the monolignols [7]. The dehydrogenation to monolignol radicals has been attributed to different classes of oxidative proteins, such as peroxidases, laccases, polyphenol oxidases and coniferyl alcohol oxidase. However, it is still unknown which particular peroxidase and laccase genes are directly involved in this process during secondary wall formation.

Genes involved in carbohydrate biosynthesis

Other genes involved in carbohydrate synthesis were also upregulated in DH stems (Table 5). One of these, RADIAL SWELLING 10 (RSW10) (BN17967), was upregulated in the top section of DH relative to the top section of YN (Table 5). RSW 10 encodes a ribose 5-phosphate isomerase involved in the formation of uridine used for the synthesis of UDP-sugars.

Plant primary carbohydrate metabolism is complex, yet flexible, and is regulated at many levels. To gain molecular insights into the carbohydrate metabolic process in B. napus stem, we retrieved genes involved in cellular carbohydrate metabolic process from www.arabidopsis.org . A total of 419 loci were retrieved and compared to identify upregulated genes in the arrays. In total, 85 genes were found to be upregulated of which 52, 33, 35 and 29 belonged to DH4:YN4, DH1:YN1, DH4:DH1 and YN4:YN1, respectively (Additional file 14: Table S8).

Targets for transcriptional regulation of cell wall formation

The largest group of transcription factors with differential expression (25 gene loci; 13 upregulated and 12 downregulated) consisted of Cys2His2-like (C2H2) zinc finger genes. C2H2 zinc fingers are well known and display a wide range of functions, from DNA and RNA binding to the involvement in protein-protein interactions. In our array data, the numbers of differentially expressed members (> 2-fold) for this family in DH4:YN4, DH1:YN1, DH1:DH4 and YN1:YN4 were 38, 25, 11 and 11, respectively (Table 6; Additional files 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11 and 12: Tables S1A and B to S6A and B). BN21811 related to At1g14580, a member of this family, was over-expressed in all 4 groups whereas the other 3 family members were upregulated in 3 groups. Another gene, BN15654 related to At5g04240, was upregulated in the basal section of DH relative to the basal section of YN. This latter gene is a member of the jumonji group of C2H2 transcription factors (Table 6; Additional files 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11 and 12: Tables S1A and B to S6A and B).

C3HC4 zinc finger genes formed the next largest gene category (Table 6). One of these genes BN22477 (related to At5g05830) was strongly upregulated in three of the four tissue comparisons and another in two comparisons. In contrast, BN17974 (related to At3g09770) was strongly down regulated in the DH4:YN4 and DH4:DH1 comparisons. B helix loop helix transcription factors were the third most differentially expressed genes in DH1:YN1, DH1:DH4 and YN1:YN4, but to a lesser extent than in DH4:YN4 (Table 6). In more detail, BN155636 related to At2g18300 (AtbHLH64) and BN12702 related to At1g05805 (AtbHLH128) were over-expressed in three groups of our array (Additional files 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11 and 12: Tables S1A and B to S6A and B). AtbHLH64 is known to play a role in Arabidopsis cytokinin signaling [43]. Though the role of bHLH TFs in lignin regulation is correlative rather than firmly established, we cannot rule out their possible involvement considering the common precursors of the lignin and flavonoid pathways and the known influence of bHLH on the flavonoid pathway genes, as well as the known interaction between bHLHs and MYB TFs in other physiological and developmental systems such as Brassica trichome development [44].

Plant basic region/leucine zipper motif (bZIP) transcription factors regulate processes including pathogen defense, light and stress signaling, seed maturation and flower development in Arabidopsis [45]. The numbers of bZIP family members differentially expressed more than 2-fold in DH4-YN4, DH1-YN1, DH1-DH4 and YN1-YN4 were 24, 11, 11 and 11, respectively in the array data (Table 6). BN15088 related to At2g18160 (ATBZIP2) were upregulated in 3 groups.

AP2 (APETALA2) and EREBPs (ethylene responsive element binding proteins) are a family of transcription factors unique to plants and share a highly conserved region of about 60 to 70 amino acids (the so-called AP2 domain) with no apparent similarity outside this domain. The numbers of differentially expressed members (> 2-fold) for this family in DH4-YN4, DH1-YN1, DH1-DH4 and YN1-YN4 were 16, 10, 7 and 2, respectively, in the array data (Table 6). In detail, BN16694 related to At1g68550 (AtERF118), was upregulated in 3 groups in our microarray experiments. Another of these genes, BN16341 related to At5g67190 (DREB2) was upregulated in the basal section relative to the top section of DH relative to YN. DREB2 contains the conserved DNA binding domain found in ERF (ethylene response factor) and AP2 (APETALA2).

NAC proteins constitute one of the largest families of plant-specific transcription factors, and are expressed in various developmental stages and tissues [46]. The NAC domain was originally characterized from consensus sequences from petunia NAM and from Arabidopsis ATAF1, ATAF2 and CUC2. In our array data, the number of differentially expressed members (> 2-fold) for this family in DH4-YN4, DH1-YN1, DH1-DH4 and YN1-YN4 were 15, 12, 8 and 5, respectively (Table 6). One member of this family, BN13387 related to At1g61110, was over-expressed in all 4 groups of our array and was previously reported to be expressed in stamens [47]. Other NAC transcription factors were only upregulated in 3 groups (Table 6). BN24136 was highly upregulated and is related to At3g49530 (ANAC062), a gene that plays an integrative role in plant responses to stress [48]. BN17366 related to At1g01720 is another Arabidopsis NAC-domain containing protein 1 (ATAF1) upregulated in the basal section relative to the top section of DH stems.

In addition to the six classes of over-represented transcription factors, MYB factors were more than 2-fold upregulated in our microarray comparisons of lignified and less lignified tissues. MYB factors represent a family of proteins with a conserved R2R3 (stimulating) or R3 (inhibiting) MYB DNA-binding domain. They often bind to bHLHs. In our array data, the numbers of overexpressed members of this family in DH4-YN4, DH1-YN1, DH1-DH4 and YN1-YN4 were 6, 3, 8 and 7, respectively (Table 6). BN17384 related to At5g60890 (AtMYB34) was expressed in all 4 groups of the array data while 2 other members of this family were upregulated in 3 groups (Table 6). BN26939 related to At2g47460 (AtMYB12) was also strongly upregulated in the array data.

In order to highlight the promising genes for future functional characterization from these elaborated sets, 25 TFs with representation from families mentioned above, plus six biochemical or physiological genes, were selected and their Arabidopsis knockdown mutants retrieved from TAIR (Additional files 14, 15, 16 and 17: Tables S8 to S11). In our analysis on a subsample of these lines, 4 TF lines and 4 non-TF lines showed lower lignin content than the Col-0 control line, while 2 TF lines and 3 non-TF lines showed enhanced lignin content (Additional files 16 and 17: Tables S10 and 11).

Biochemical genes

Biochemical genes specifying monolignol formation, such as COMT [10], CCoAOMT 1 [49], SAMS [50], CAD [51], REF8 (C3H) [52] and F5H (FAH5) [53], were more highly expressed in the more lignified stem tissues we examined. Upregulation in more lignified tissue also occurred for genes specifying monolignol glucosylation, and monolignol polymerization, including peroxidases and laccases, when pairs of tissue sections were compared. Cellulose biosynthesis, the entire cell wall and plant development appeared to be affected by some of these differentially expressed Brassica genes. For example, we found that Brassica CES genes [54, 55] and peroxidases [56] were also more upregulated in more lignified tissues. Passardi et al. [56] showed that two highly homologous Arabidopsis peroxidases, AtPrx33 and Atprx34, were involved in cell elongation, and the knockout mutant, cesa2, also had severe defects in cell wall formation and microtubule orientation [57]. However, the numbers of cellulose biosynthetic genes upregulated in DH were fewer than lignin biosynthetic genes. Hence, our comparative transcriptome analyses of lignified and less lignified material could serve as a source of differentially expressed biochemical genes that impact on cell wall recalcitrance and could be useful to enable plant molecular breeders to fine tune crop residue composition for industrial applications.

Targets for transcription regulation

Microarray analysis of differentially lignified tissues resulted in six transcription factor families being more prominent in expression changes than other major TFs. Transcription factors have recently received attention from plant molecular breeders because of their ability to modify entire biochemical pathways or developmental families of genes. For example, overexpression of the Arabidopsis GLABRA3 (AtGL3) bHLH transcription factor in B. napus enhanced trichome coverage and insect resistance [44, 58].

Here, we screened our microarray data for the top 10 TF families present in lignified B. napus stems with the intent of beginning the process of characterizing the function of a subset of regulatory genes using Arabidopsis knockdown mutants. The most modified TF family in our stem microarrays was the zinc finger family. Several zinc finger proteins (ZFPs) in plants, e.g. Arabidopsis and petunia, have already been found to be involved in a variety of processes such as the regulation of floral organogenesis, leaf initiation, lateral shoot initiation, gametogenesis and stress response [59]. The C2H2 group was the most highly upregulated set of zinc finger families in our analyses. This latter group is known to be involved in the regulation of Arabidopsis flowering time [60]. Prigge and Wagner [61] reported that SERRATE (At2g27100), another member of the C2H2 family, plays a role in embryogenesis and is transcribed in shoot meristem and in emerging organ primordia throughout development. Dong and co-workers [62] also showed that SERRATE is required for the accurate processing of microRNA (miRNA) precursors in the plant cell nucleus. Moreover, some zinc finger genes highly expressed in trichomes have been reported to play a role in cell wall biosynthesis [63]. Thus, selections from among these differentially expressed zinc finger genes could target B. napus regulatory genes that can modify entire developmental programs.

BHLH TFs are one of the largest TF families in Arabidopsis. Phylogenetic analysis, divergence in binding site specificity and their varied ability to engage in homodimerization and heterodimerization events strongly support their roles in a multiplicity of transcriptional programs [64, 65]. BHLH TFs are known to play important roles in the flavonoid branch of the phenylpropanoid pathway; hence it makes sense that we would find them so plentiful in the closely related lignified tissue. In Arabidopsis seeds, the MYB transcription factor TRANSPARENT TESTA 2 (TT2) forms a complex together with the bHLH transcription factor TT8 and a WD40 scaffold protein to control the expression of BANYULS, a proanthocyanidin biosynthetic gene [66]. In a leaf study, TT8 and its closest homolog GL3 and ENHANCED GLABRA (EGL3) regulate flavonoid biosynthesis through interaction with two homologous MYB proteins PRODUCTION OF ANTHOCYANIN PIGMENTS 1 (PAP1, MYB75) and PAP2 (MYB90) [67]. In our stem microarrays, several bHLH genes were also differentially expressed along with a number of MYB genes. These genes should be paired and tested to determine whether any are binding partners and whether their modification affects lignin in transgenic plants or mutants. Already, the alfalfa bHLH TT8 homologue has been shown in RNAi studies to downregulate carbohydrate and dry matter accumulation, but not lignin, in alfalfa forage [68].

AP2/EREBP genes represent the fifth most differentially expressed genes in the B. napus stem microarrays. These genes play a variety of roles throughout the plant life cycle; from being key regulators of several developmental processes, like floral organ identity determination or control of leaf epidermal cell identity, to forming part of the mechanisms used by plants to respond to various types of biotic and environmental stresses [69]. NAC genes ranked sixth among our most differentially expressed TFs in the stem microarrays. Mitsuda et al. [70] showed that NAC SECONDARY WALL THICKENING PROMOTING FACTOR1 (NST1) and NST3 are key regulators of the formation of secondary walls in woody tissues of Arabidopsis. Lu et al. [71] showed that ATAF1 mRNA expression is strongly induced by dehydration and abscisic acid (ABA) treatment and they suggested a general role of this protein as a repressor in drought stress responses. Wu et al. [72] also implicated ATAF1 in diverse biotic and abiotic stress responses including drought, high-salinity, ABA, methyl jasmonate, mechanical wounding and Botrytis cinerea infection. Considering the expression pattern of NAC transcription factors in the array data and their role in lignin biosynthesis as well as in diverse biotic and abiotic stress responses, this family is a promising target for lignin modulation of the plant cell wall.

Although not among the most differentially expressed B. napus stem TFs, MYB-related transcription factors are an important group since they regulate different branches of secondary metabolism, in addition to the identity and fate of plant cells [73]. Studies on MYB transcription factors in Arabidopsis, pine and eucalyptus confirmed regulatory roles that MYB TFs play in lignin biosynthesis [14, 15, 74, 75]. MYB58 and MYB63 are transcriptional regulators specifically activating lignin biosynthetic genes during secondary wall formation in Arabidopsis [75]. Dominant repression of their functions led to a reduction in secondary wall thickening and lignin content, while overexpression of MYB58 and MYB63 resulted in specific activation of lignin biosynthetic genes and concomitant ectopic deposition of lignin in cells that are normally unlignified [75]. Newman and co-workers [13] suggested a role for another Arabidopsis MYB gene, AtMYB61, in cell wall deposition in normal vascular development [76]. However, transcription of these three MYB genes was unaltered in our array experiments. Instead, MYB12 and MYB34 were differentially regulated in the arrays. MYB34 plays a key role in the regulation of indole glucosinolate homeostasis in Arabidopsis [77], and MYB12 acts as a flavonol-specific activator of flavonoid biosynthesis [78]. Luo et al. [79] showed that MYB12 can also activate the caffeoyl quinic acid (CQA) biosynthetic pathway when expressed in a tissue-specific manner in tomato. The role of these two MYB TFs in secondary metabolism pathways and the correlation of their expression patterns with stem lignification establish this family as an important target in the modulation of lignin biosynthesis.

Analysis of Arabidopsis mutants for lignin, eight mutant lines showed lower lignin content than the Col-0 control line, while five showed enhanced lignin content. These verified lignin genes can now be used to modify Brassica lignin using our microarray oligonucleotides. Moreover, the Brassica and Arabidopsis sequences can be used to find homologues in other species. In fact, DIM1 and eEF-1Bβ1 have already been shown to modify lignin and carbohydrate in Arabidopsis [31, 32] and TT8 and HB12 have been shown to modify carbohydrate structure in alfalfa [34, 68]. In addition, several other genes (HB5, DIM1, ZFP2, eEF-1Bβ1, ANAC62) have also been used to modify alfalfa (Amyot, personal communication).