Abstract
Background
Gametogenesis is a key step in the production of ovules or pollen in
higher plants. The sex-determination aspects of gametogenesis have been
well characterized in the model plant Arabidopsis. However, little is
known about this process in androdioecious plants. Tapiscia sinensis
Oliv. is a functionally androdioecious tree, with both male and
hermaphroditic individuals. Hermaphroditic flowers (HFs) are
female-fertile flowers that can produce functional pollen and set
fruits. However, compared with male flowers (MFs), the pollen viability
and number of pollen grains per flower are markedly reduced in HFs. MFs
are female-sterile flowers that fail to set fruit and that eventually
drop.
Results
Compared with HF, a notable cause of MF female sterility in T. sinensis
is when the early gynoecium meristem is disrupted. During the early
stage of HF development (stage 6), the ring meristem begins to form as
a ridge around the center of the flower. At this stage, the internal
fourth-whorl organ is stem-like rather than carpelloid in MF.
A total of 52,945 unigenes were identified as transcribed in MF and HF.
A number of differentially expressed genes (DEGs) and metabolic
pathways were detected as involved in the development of the gynoecium,
especially the ovule, carpel and style. At the early gynoecium
development stage, DEGs were shown to function in the metabolic
pathways regulating ethylene biosynthesis and signal transduction
(upstream regulator), auxin, cytokinin transport and signalling, and
sex determination (or flower meristem identity).
Conclusions
Pathways for the female sterility model were initially proposed to shed
light on the molecular mechanisms of gynoecium development at early
stages in T. sinensis.
Keywords: Morphology and anatomy, Androdioecious tree, Tapiscia
sinensis, Sex differentiation, Transcriptome analysis
Background
Angiosperm plants exhibit a wide variety of breeding systems [[35]1,
[36]2]. Among them, hermaphroditism is thought to be the ancestral
breeding system in angiosperms [[37]3, [38]4], and dioecy appears to
have evolved from hermaphroditism multiple times [[39]1, [40]5–[41]7].
Although two main pathways for the evolution of dioecy have been
proposed, dioecy that evolves from an ancestral hermaphrodite or from a
monoecious species is considered to be the most common [[42]1,
[43]6–[44]9]. On this pathway, it is believed that hermaphroditic
individuals coexist with unisexual mutants that have lost either their
male (gynodioecy) or their female (androdioecy) function [[45]8]. Among
the variety of reproductive mechanisms in flowering plants, one
permutation – androdioecy (mixtures of males and hermaphrodites) – is
distinguished by its rarity and represents less than 0.005% of all
angiosperms [[46]10]. However, androdioecy is thought to be relatively
unstable; models of mating system evolution predict that androdioecy
should be a brief stage between hermaphroditism and dioecy (separate
males and females), or vice versa [[47]5, [48]8, [49]9, [50]11–[51]15].
Loss of function is not a big challenge. The biggest challenge to
female sterilty is due to natural selection. Typically if an individual
looses female fertility, its fitness will be halved, and natural
selection should act against it, eliminating the mutation.
Separate sexes have evolved on numerous independent occasions from
hermaphroditic ancestors in flowering plants [[52]16]. Mechanisms
controlling sex can be genetic, epigenetic (physiological and
environmental) [[53]17–[54]24], or plant hormones [[55]2, [56]25]: (1)
Genetic control of sex determination is driven by chromosomes and sex
determining genes in the dioecious plants. According to Heikrujam et
al. (2015), genetic sex determination may be due to a single locus or
multiple loci either unlinked or tightly linked on autosomes [[57]2].
For example, in persimmon, only a single locus is responsible for sex
determination [[58]26, [59]27] while the evolution of separate sexes in
asparagus invokes sterility mutations at two linked loci [[60]27,
[61]28]. (2) Hormonal regulation of unisexual flower development
[[62]29]. Unisexual flower development in Cucumis sativus (cucumber)
and Cucumis melo (melon) is regulated by the interaction of
environmental cues, plant hormones and genetic factors that
differentiate gender phenotypes [[63]24, [64]25, [65]30]. In kiwifruit,
a Y-encoded suppressor of feminization, Shy Girl, arose via
lineage-specific duplication of a cytokinin response regulator
[[66]29]. At present, in some plant species, such as Arabidopsis
[[67]31, [68]32], cotton [[69]33], rice [[70]34], tomato [[71]35],
Silene latifolia [[72]36, [73]37], Zea mays [[74]38–[75]40], Diospyros
lotus [[76]26], and cucumber [[77]25], significant progress in
elucidating the molecular mechanism of pistil development and sex
determination has been achieved [[78]24–[79]26]. Molecular studies in
these species have led to important advances in our understanding of
plant sex determination and differentiation [[80]25, [81]26, [82]41].
For example, ANT promotes ovule primordium growth; the CUP-SHAPED
COTYLEDON genes (CUCs) play a role in the establishment of the ovule
primordium boundaries; when the functions of stk, shp1, and shp2 were
lost in a triple mutant, fewer ovules were initiated, and ovule
development was severely disrupted. Even more remarkably, the ABCDE
model of flower determination and development has indicated that a
specific class of MADS-box genes are key regulators of pistil
development [[83]42], such as AGAMOUS (AG), APETALA 2 (AP2), BELL 1,
INNER NO OUTER (INO), AINTEGUMENTA (ANT), SPOROCYTELESS (SPL/NZZ), and
SUPERMAN (SUP) [[84]43, [85]44].
Recent advances in the evolution of androdioecy showed that androdioecy
arose from hermaphroditism because some male flowers have
non-functional pistils, and that these residual organs are common in
androdioecious and gynodioecious species that arose from
hermaphroditism [[86]13, [87]45–[88]49], and that in species where
androdioecy arose from dioecy, non-functional pistils are absent
[[89]50, [90]51]. In the functionally androdioecious Tapiscia sinensis
Oliv. (Tapisciaceae) [[91]52–[92]54], the male individuals have
pistillode, some of which have a shorter column and abortive ovule,
while the vast majority of pistils lack ovaries (Figs. [93]1 and
[94]2a). Thus, this plant might provide a perfect model of male
individuals that originated from an ancestral hermaphrodite.
Furthermore, the stage at which pistil abortion occurs and sex
differentiation is determined by sex chromosomes, sex-determining
genes, or hormone regulation remians unclear.
Fig. 1.
[95]Fig. 1
[96]Open in a new tab
Photographs showing the developmental stages of the male and bisexual
flowers. (a–e) Male inflorescences at stages 5, stage 6, stage 8, stage
10 and stage 13, respectively. (a–e) Bisexual inflorescences at stage
5, stage 6, stage 8, stage 10, and stage 13, respectively. Bar=1 cm
Fig. 2.
[97]Fig. 2
[98]Open in a new tab
Characteristics of the flowers. a Morphology of HF (left) and MF
(right). b Compared with HF (right), there are three types of gynoecia
in MF. c Scanning electron microscope (SEM) observation of HF. d-f SEM
observation of three types of MF. g Change in pistil length (PL) over
time. h Change in ovary transverse diameter (OTD) over time; i OTD in
HF - Type I, Type II, and Type III comparison. j PL in HF - Type I,
Type II, and Type III comparison. Each data point represents the mean
value of 30 technical replicates
To determine the molecular mechanism of female sterility in male
individuals, we compared the developmental anatomy and performed
transcript profiling using the RNA-Seq technique of the gynoecia of the
HF and MF. A number of candidate genes and related pathways were
revealed, which provided new insights into the genetic and biochemical
controls for early stages of pistil development, and our results will
be helpful in shedding light on the mechanism of sex differentiation in
the functionally androdioecious trees of Tapiscia.
Results
Floral morphology
The morphology of MFs is funnel-shaped, whereas that of HFs is water
goblet-shaped (Fig. [99]2a) and have a normal pistil with a stigma, a
style, and an ovary (Fig. [100]2c). In comparison with HF, the gynoecia
of MF were slow-growing and much shorter from the early development
stage to flowering (Fig. [101]2b, g). In natural populations, the MFs
have three types (Fig. [102]2b, d-f). In comparison with Type III, the
ovary transverse diameters (OTD) in Type I and Type II is smaller
(p < 0.05) (Fig. [103]2i), while no significant differences exist among
their pistil lengths (PL) (Fig. [104]2j). Correspondingly, the
flowering duration of MF plants was prolonged until complete
fertilization of HF because they had no seed setting and produced more
flowers. MF had shorter pistils, and smaller ovules than HF. At
flowering time, MF bloomed at a PL of 1.17 ± 0.17 mm (Fig. [105]2g) and
an OTD of 0.518 ± 0.05 mm (Fig. [106]2h). In contrast, HF bloomed at a
PL of 2.63 ± 0.54 mm (Fig. [107]2j) and an OTD of 1.15 ± 0.08 mm (Fig.
[108]2i). The OTD (T = -7.102, df = 12, P < 0.01) and the PL
(T = − 2.876, df = 12, P < 0.05) were clearly shorter in MF than in HF.
Remarkably, in MFs, the PL and OTD remain unchanged when the bud
vertical diameter grew to 0.75 ± 0.05 mm (stage 9, May 16th) (Fig.
[109]2d). Later, some limited growth may continue from stage 9 to
flowering for the carpel in MFs, which is often only one-quarter to
one-half the full length that in the HFs (Fig. [110]3j).
Fig. 3.
[111]Fig. 3
[112]Open in a new tab
Comparison of the sex differentiation between Hermaphrodite Flowers and
Male Flowers. In T. sinensis, male and bisexual flowers have the same
calyx initiation, petal initiation and stamen initiation pattern (from
stage 1 to 5), so we begin to show their differences from stage 5. a
The morphology of 5 petals and 5 stamens. a-b~g The morphology of HF:
a, b, c Two or three gynoecium primordia begin to form as a ridge (the
ring meristem) around the center of the flower; c, d at stage 7, the
gynoecium grows as a hollow tube; at stage 8, the gynoecium tube grows
taller and wider; during stage 10, ovule primordiu were observable; g,
h at flowering time, the gynoecium becomes ready for fertilization, and
while the stamens extend out of the petals, fertilization occurs. a-b~g
The morphology of MF: a, b at stage 6, the gynoecium primordium also
begins to form, but the central zone or gynoecium primordium on the
flower apex was fused as a small bulge; c~e Abnormal pistil primordium
grows rapidly, and finger-like. f, g The pistil was characterized by a
solid, stem-like structure. vp, vestigial pistil; oi, outer integument;
ii, inner integument; rm., ring meristem; pp., pistil primordium; st,
stamen; pi, pistil; cz, central zone
The SEM results showed that there were no obvious differences observed
between MF and HF before stage 5 (Fig. [113]3a). However, at stage 6,
in HFs, the gynoecium primordial begins to form as a ridge of raised
cells around the center of the flower apex (Fig. [114]3a, b, c). At
stage 7, the gynoecium grows as a hollow tube (Fig. [115]3c, d); At
stage 8, the gynoecium tube grows taller and wider, the gynoecium
continues to grow upward to form a continuous hollow cylinder, and
placenta primordia were observable; during stage 10, in the centre of
the receptacle, ovule primordium were observable (Fig. [116]3e, f); at
the flowering time (stages 11–14), the gynoecium becomes ready for
fertilization, and when the stamens extend out of the petal,
fertilization occurs (Fig. [117]3g, h).
In most MFs, during stage 6, the gynoecium primordium began to form
(Fig. [118]3b), but the central zone was fused as a hemispherical bulge
and no differentiation of the ring meristem is present (Fig. [119]3c,
d), which was obviously different from the ring meristem in HFs (Fig.
[120]3b, c). The pistil structure of the MF was characterized by the
formation of a solid, stem-like structure (Fig. [121]3e-g). The key
stage at which ovule abortion occurred was stage 6 (Fig. [122]3b, B).
Structure of the flower
In T. sinensis, the course of individual flower development can be
divided into five phases: early differentiation, calyx initiation,
petal initiation, stamen initiation, and carpel initiation. The
perianth was the first organ to be formed on the floral apex
(Fig. [123]4a, b). Sepals differentiated concurrently and curved to
cover the central zone (Figs. [124]3a and [125]4c). In late stage 5,
petal and stamen primordia are initiated (Fig. [126]4d). As the stamens
formed, there was a bud vertical diameter of 76.92 ± 1.98 μm, and the
floral apex became broadly concave (Fig. [127]4e). The central zone
(floral meristem) of the flower bulges out to form a platform on which
the gynoecium primordia will develop (Figs. [128]4f, g and [129]3b). At
stage 6, the sepals enclose the flower bud, and two or three carpel
primordia begin to form, with raised cells around the floral apex
constituting a ridge (ring meristem) (Fig. [130]4h, i). At stage 7, the
medial stamens become stalked at their bases, and the ring meristem
grows as a hollow tube (vertical section) (Fig. [131]4j~l); at stage 8,
the placenta primordium is observable (Fig. [132]4m~o); then, the ovule
forms from the placenta (Fig. [133]4p, q). The two or three carpels
gradually fused with each other in hermaphrodites and eventually
developed into a normal pistil with a stigma, a style, and an ovary
(Additional file [134]4: Figure S1).
Fig. 4.
[135]Fig. 4
[136]Open in a new tab
Developmental anatomy of the HF. a the flower primordium initiated. b
The calyx primordium initiated. c When the calyx covered the central
zone of the flower, the petal primordium initiated. d, e The cetral
zone of the flower becomes broad, then the stamen primordium intiated.
f, g The central zone (cz) of the flower apex raised as a platform. h,
i During stage 6, the sepals enclose the central zone of the flower bud
and the carpel primordia begin to form, where the raised cells surround
the floral apex as a ridge. j, k Gynoecium continues to grow upward to
form a continuous hollow cylinder. l The central zone becomes broad,
where the placenta primordia will be formed; m, n At stage 7, the
medial stamens become stalked at their bases, the gynoecium grows as a
hollow tube, and placenta primordia are observable. o Rapid growth of
the placental primordia and style, and the unilocular formed. p Period
of the formation of ovule primordium. q In the centre of the ridge,
ovules will form. Bar = 100 μm. Cap, calyx primordium; pep, petal
primordium; cz, central zone; sp., stamen primordium; rm., ring
meristem; st, stamen; p, petal; pep, petal primordium; pl, placenta
primordium
From stage 1 to stage 5, hermaphroditic individuals and males have the
same developmental process and morphological structures (Fig.
[137]4a-k). However, after stage 5, the carpel primordium of the centre
of the floral apex in MFs exhibits different patterns from that of HFs.
Longitudinal sections of MF from the same tree showed that there are
three types of pistils. Type I: Most of the MFs without ovules or
ovaries were characterized by the formation of a solid, stem-like
residual pistil (Fig. [138]5a-i). Type II: A tiny minority (1.6%) of MF
had one ovary but no ovule structure (Fig. [139]5j-n). Type III: MFs
have a well-developed pistil that is sterile (3.33%); the ovary has one
chamber that has one ovule that has double integuments and is
crassinucellate and anatropous (Fig. [140]5o-u); however, when the MF
falls off after maturation, termination of embryo sac development
occurs at the triad stage (Fig. [141]5v; Additional file [142]4: Figure
S1).
Fig. 5.
[143]Fig. 5
[144]Open in a new tab
Developmental anatomy of the MF. Three types of gynoecia in MF: Type I,
Gynoecia without ovary (a-i); Type II, ovules were arrested at stage 6
(j-n); Type III, the MF have a well-developed pistil that is sterile
(o-v), particularly. Bar=100 μm
A pollination analysis found that a few pollen grains from self or HF
could germinate on MF stigmas 6 h after pollination. Remarkably,
compared with HF (Fig. [145]6a-c), in MF, although pollen tubes could
grow in a spiral manner over the surface of the papilla cells, some of
them penetrated into stigmas or styles and then succeeded in reaching
the abnormal ovules 24 h after pollination (Fig. [146]6a-c). Further
observation showed that the papilla cells on MF stigmas (Fig. [147]6d,
e) were similar to HFs at the maturation stage (Fig. [148]6d, e). The
anther and pollen grain sections show that many sterile pollen grains
(spg) were observed in HF (Fig. [149]6f) but not in MF (Fig. [150]6f).
There are two types of pollen in hermaphroditic flowers (Fig. [151]6g),
while there is only one type of pollen in MF (Fig. [152]6g). In
contrast to the perforate tectum of the olivary pollen grains in MF
(Fig. [153]6h, i), the pollen grains in HF have a reticulate tectum
(Fig. [154]6h, i).
Fig. 6.
[155]Fig. 6
[156]Open in a new tab
Pollen germination, pollen tube elongation and papilla cells in MF
(a-e) and HFs (a-e). a-c The pollen tubes are screw-shaped and
penetrate into papilla cells in the three types of MF at 24 h after
self-pollination. a-c At 24 h after self-pollination, the pollen tubes
passed into the micropyle of the ovule of HF. d, e Papilla cells of
undivided stigma. d, Bilobed stigma. e, Trilobate stigma. f, Anther and
pollen grain sections of HF; sterile pollen grains (spg) were observed.
g, Globose pollen grains in HF. h, Globose pollen grain; i, Detail of
perforate tectum. f, Anther and pollen grain sections. g, Olivary
pollen grains. h, Olivary pollen grain in MF; i, Detail of reticulate
tectum
Illumina sequencing and sequence assembly
In this project, the Illumina HiSeq™ 4000 sequencing platform was used
for transcriptome sequencing, and 962,641,935 clean reads were
obtained, with a mean length of 150 base pairs (bp). The average GC
percentage was 44.69%. The clean reads were de novo assembled using
Trinity [[157]55] into 63,573 contigs, with a mean length of 1132 bp
and an N50 of 1892 bp. Then contigs were used for sequence clustering
with the software TGICL (TIGR Gene Indices clustering tools) to form
unigenes [[158]56]. The total number of unigenes was 52,945, with a
mean length of 973 bp and N50 of 1768 bp. The average depth of sequence
coverage was 82.42%. In all, 36,353 unigenes within 201–1000 bp, which
was 57.59% percent. Additionally, 7766 unigenes (14.67%) ranged from
1001 to 2000 bp in length, while 7408 of the unigenes (13.99%) were
over 2000 bp in length.
Differentially expressed genes (DEGs) and enrichment analysis
Using the series test of cluster (STC) [[159]57], a total of 9989 and
7326 unigenes were identified and transcribed in the M1-M2-M3 term
(TERM1)_vs_H1-H2-H3 term (TERM2), respectively (Additional file [160]4:
Figure S3). M1, M2 and M3 represent MF at stages 5, 6 and 10,
respectively, while H1, H2 and H3 represent bisexual flowers at stages
5, 6 and 10, respectively. In TERM1, the expression patterns of 13,761
genes were analysed, and eight model profiles were used to summarize
(Additional file [161]4: Figure S2); five expression patterns of genes
showed significant p-values (p < 0.05) marked with coloured boxes,
including profile 3, profile 4, profile 5, profile 6, and profile 7
(Additional file [162]4: Figure S3). In TERM2, the expression patterns
of 16,130 genes were analysed, and eight model profiles were used to
summarize (Additional file [163]4: Figure S2); three expression
patterns of genes showed significant p-values (p < 0.05) marked with
coloured boxes, including profile 1, profile 6, and profile 7
(Additional file [164]4: Figure S3). Each box represents a model
expression profile with the model profile number and p-value. To
further clarify the functional differences between MF and HF at stage
6, using a false discovery rate (FDR) ≤ 0.05 and the absolute value of
|log2 ratio| ≥ 1 as criteria, 9623, 14,347 and 12,049 unigenes
exhibited significantly different expression levels in M1_vs_ H1, M2
_vs_H2, and M3_vs_H3, respectively (Additional file [165]4: Figure
S4A-C).
To determine whether the DEGs between MF and HF are significantly
related to specific pathways or biological function, GO and KEGG
databases were used for enrichment analysis (Additional file [166]4:
Figure S4D-F). For TERM1, TERM2, M1_vs_H1 DEGs, the GO terms that
showed enrichment included “flavonoid biosynthetic process”
(GO:0009813), “secondary metabolic process” (GO:0019748), “hormone
biosynthetic process” (GO:0042446), “cytokinin biosynthetic process”
(GO:0009691), “glycoside biosynthetic process” (GO:0016138), and
“regulation of fertilization” (GO:0080154) (Additional file [167]3). In
contrast, for TERM1, TERM2, M2_vs_H2, the enriched classification terms
were “cell wall organization or biogenesis” (GO:0071554), “external
encapsulating structure organization” (GO:0045229), “glucose metabolic
process” (GO:0006006), “embryo development” (GO:0009790), “response to
inorganic substance” (GO:0010035), and “response to chemical”
(GO:0042221) (Additional file [168]3).
Pathway assignment was carried out by KEGG to analyse the biological
functions of DEGs. For the TERM1, TERM2, and M1_vs_H1 comparison (Fig.
[169]7a), pathways related to “Phenylpropanoid biosynthesis” (ko00940),
“Plant-pathogen interaction” (ko04626), “Plant hormone signal
transduction” (ko04075), “Glucosinolate biosynthesis” (ko00966),
“Diterpenoid biosynthesis” (ko00904), “Starch and sucrose metabolism”
(ko00500), “Zeatin biosynthesis” (ko00908), “Ribosome” (ko03010), and
“Flavonoid biosynthesis” (ko00941) were enriched (Additional file
[170]3). For the TERM1, TERM2, and M1_vs_H1 comparison (Fig. [171]7b),
the main pathways included “Glycolysis/Gluconeogenesis” (ko00010),
“Pentose and glucuronate interconversions” (ko00040), and “Tryptophan
metabolism” (ko00380), “Citrate cycle (TCA cycle)” (ko00020), “Starch
and sucrose metabolism” (ko00500), “Plant hormone signal transduction”
(ko04075), and “Zeatin biosynthesis” (ko00908) (Additional file
[172]3). The enrichment analysis revealed a total of 63 DEGs related to
plant hormone signal transduction. Thus, the process of “plant hormone
signal transduction” was concluded to be a putative pathway affecting
T. sinensis flower development (Additional file [173]3).
Fig. 7.
[174]Fig. 7
[175]Open in a new tab
Unique and shared DEGs among different stages for both HF and MF. a The
number of DEGs among TERM1, M1_vs_H1, TERM2. b The number of DEGs among
TERM1, M2_vs_H2, TERM2. Venn diagrams were drawn using a online tool
Venn Diagrams ([176]http://bioinformatics.psb.ugent.be/webtools/Venn/)
In brief, the DEGs were characterized using the GO and KEGG databases.
Based on recent research results [[177]25, [178]29, [179]31, [180]42,
[181]44], we will focus on gene regulatory network models related to
carpel and ovule development, “Plant hormone signal transduction”
(ko04075), “Zeatin biosynthesis” (ko00908), “hormone biosynthetic
process” (GO:0042446), and “cytokinin biosynthetic process”
(GO:0009691) (Additional files [182]2).
Gene regulatory network (GRN) models related to carpel and ovule development
The results of scanning electron microscopy (SEM) and light microscopy
(LM) show that most of the vestigial pistils of MF of T. sinensis
developed into finger-like structures without ovary or ovule. Pollen
germination and pollen tube growth in MF indicated that pollen tubes
can grow in the stigma/style and that they are functional, and
therefore, ovules are never present is the main factor responsible for
female-sterility of male flowers (Fig. [183]8a). In Arabidopsis, ovule
cell fate is controlled by ovule identity genes [[184]31, [185]44,
[186]58–[187]60]. The pathway of regulation of Arabidopsis ovule
development is shown in Fig. [188]7, and key regulators include AG
[[189]30, [190]61, [191]62], SPL [[192]63], INO [[193]64], ANT and BEL1
[[194]31, [195]65]. Carpels are the most complex structures within
flowers; a GRN (gene regulatory network model) underlies their
development in Arabidopsis [[196]31], and in T. sinensis, the key
regulators were identified, including WUS (Unigene0010686), AGAMOUS
(Unigene0020941), SUP (Uniene0047137), SPL (Unigene0006964), INO
(Unigene0025624), STM (Unigene16944), ANT1 (Unigene0002499), STK
(AGL11) (Unigene0007536), SHP1 (Unigene0013614), AP2 (Unigene0015664),
TSO1 (Ungene0002282), SPT (Unigene0034455), LUG (Unigene0041412), PI
(Unigene0020391), SEUSS (Unigene0033961) and BEL1 (Unigene0039261)
(Fig. [197]8b, c). In the MF-HF comparison, these 16 genes showed
significantly different expression levels. The results indicated that
these genes may influence the flower development.
Fig. 8.
[198]Fig. 8
[199]Open in a new tab
Expression of candidate genes involved in pistil development. b Flower
development in the MF-HF difference model. b Based on prior work in
Arabidopsis [[200]23, [201]58–[202]63], a proposed model of
gene-dependent pathways playing potential roles in regulating gynoecium
development. c Expression of candidate genes involved in ovule
development. d Relative expression levels of 20 selected genes at
different stages of floral bud development in MF and HF. Relative
expression levels were calculated by the 2-ΔΔCt method with actin as a
standard. 11-Apr., 21-Apr., and 4-May represent the pistils of flowers
at stage 5, stage 6, and stage 8, respectively. Different colours from
blue to red show the relative log2 (expression ratio)
To confirm that the unique-match genes from the Illumina sequencing and
bioinformatics analysis were indeed differentially expressed, a total
of 16 genes were selected from the DEGs related to flower development
(including carpel, ovule, pollen) for quantitative RT-PCR assays. For
16 genes, the qRT-PCR analysis revealed the same expression trends as
the RNA-Seq data, despite some quantitative differences. All genes had
different expression levels between MF and HF at stage 6 (Fig.
[203]8d). For example, at stage 6, Unigene0020941, a transcription
factor that regulated carpel development, was down-regulated (0.018 to
0.045) in MF relative to HF by qRT-PCR, the same trend as the one (39
to 338) from RNA-Seq.
Identification of genes involved in ethylene biosynthesis and signal
transduction
Ethylene regulates many aspects of plant development and responses to
biotic and abiotic stress. We identified ten DEGs involved in the
ethylene biosynthesis and signal transduction process (Fig. [204]9a):
ACO (Unigene0029908), ACO1 (Unigene0037768), ACS (Unigene0021193), ERS1
(Unigene0024060), CTR1 (Unigene0036282), TIR1 (Unigene0037222), ETR1/2
(Unigene0020715, Unigene0027393), EIN2/3/4 (Unigene0029458,
Unigene0050718, Unigene0003076) (Fig. [205]9b). To confirm the results
of RNA sequencing, ten DEGs related to carpel development were selected
for qRT-PCR analysis (Fig. [206]9c). The results showed similar
expression trends to the transcriptome analysis, and only 2/16 genes
(EIN2 and EIN3) had different results in qRT-PCR vs RNA-Seq, confirming
the reliability of the RNA-Seq data. 1-aminocyclopropane-1-carboxylic
acid (ACC) synthase (ACS) is the rate-limiting enzyme in ethylene
synthesis, and it is important in regulating ethylene biosynthesis.
Unigene0021193 was downregulated at stage 6 in the M-H comparison. Two
ACC oxidase homologous genes showed different expression patterns, but
ACO and ACO1 were both upregulated in stage 6 in the MF-HF comparison.
Fig. 9.
[207]Fig. 9
[208]Open in a new tab
Expression of candidate genes involved in ovule development. a Ethylene
biosynthesis and signal transduction pathway, which comes from a KEGG
pathway enrichment analysis [[209]76, [210]77]. b Expression of
candidate genes involved in ethylene biosynthesis and signal
transduction. Gene expression was measured by log2Ratio. 11-Apr,
21-Apr, 4-May, 11-May and 21-May represent the pistils of flowers when
their development stages were 5, 6, 10, 12 and 14, respectively.
However, here, we mainly discussed stage 6, which is the most important
development stage in sex differentiation in T. sinensis (for details,
see ‘Structure of the Flower’). Different colours from blue to red show
the size of the correlation coefficient. c Relative expression levels
of 10 selected genes at different stages of floral bud development in
MF and HF. Relative expression levels were calculated by the 2-ΔΔCt
method with actin as a standard. * represents the P < 0.05 level of
significant difference, and ** represents the P < 0.01 level of
significant difference in independent-samples t-tests
In Arabidopsis, the perception of ethylene is achieved by five members
of a family of ER membrane-bound receptors: ETR1, ETR2, ERS1, ERS2 and
EIN4, some of which have histidine kinase activity [[211]72]. From mRNA
sequencing, we obtained four ethylene receptors: EIN4 (Unigene0003076),
ETR2 (Unigene0027393), ETR1 (Unigene0020715), ERS1 (Unigene0024060). At
stage 6, the relative expression of Unigene0020715 displayed
significant differences, but Unigene0027393 and Unigene0024060 were
indistinctive. Unigene0020715, the repressor of ethylene signal
transduction, was up-regulated in MF relative to HF, but Unigene0027393
displayed indistinctive expression. Similarly, Unigene0003076 displayed
almost the same trend as ETR1 during stage 6, while Unigene0024060 was
indistinctive. ERF1/2 (Unigene0025860, Unigene0005688), the downstream
gene in the ethylene signal transduction pathway, displayed
differential expression in both the M1-M2 and H1-H2 comparisons, being
up-regulated in MF relative to HF at the key stage. Obviously, at stage
6, ACO (Unigene0029908), ACO1 (Unigene0037768), CTR1 (Unigene0036282),
and EIN2 (Unigene0029458) were up-regulated in the MF-HF comparison,
while ACS and EIN3 were down-regulated in MF relative to HF.
To verify the role of ethylene in sex differentiation of T. sisensis,
we sprayed ethephon at concentrations of 150 mg/L. We found that
ethephon had a significant effect on MF transformed into HF.
Identification of genes involved in Auxin and Cytokinin signal transduction
The third related pathway was the auxin transport and signalling
pathway. As the first plant hormone studied, auxin has wide-ranging
effects on growth and development throughout the plant, including
gynoecium and ovule morphogenesis [[212]73–[213]75]. In our study, we
identified 24 DEGs involved in auxin and cytokinin signal transduction
(Fig. [214]10a, b): IAA31 (Unigene0018133), IAA26 (Unigene0035359),
IAA15 (Unigene0001116), IAA19 (Unigene0025767), GH3.1 (Unigene0015426),
SAUR2 (Unigene0017779), SAUR14 (Unigene0050945), SAUR72
(Unigene0005319), SAUR50 (Unigene0010084), SAUR12 (Unigene0018267),
SAUR31 (Unigene0049968), SAUR53 (Unigene0019490), KUP3
(Unigene0005158), MIR397B (Unigene0025946), ARR10 (Unigene0049560),
AHP4 (Unigene0004154), AHP6 (Unigene0011561), AHP4 (Unigene0004325),
IAA1 (Unigene0001116), ARF17 (Unigene0029174), CALS5 (Unigene0019058),
PIN6 (Unigene0020253), PIN2 (Unigene0033811), PIN1 (Unigene0000533). To
confirm that the unique-match genes from the Illumina sequencing and
bioinformatics analysis were indeed differentially expressed, a total
of 9 genes were selected from the DEGs related to flower development
for quantitative RT-PCR assays. Among 9 of the 24 genes with different
expression levels between MF and HF (Fig. [215]10c) at stage 6, six
genes were up-regulated in MF relative to HF, and three genes were
down-regulated.
Fig. 10.
[216]Fig. 10
[217]Open in a new tab
a Auxin and cytokinin signal transduction. b Expression of candidate
genes involved in ethylene biosynthesis and signal transduction based
on KEGG pathway enrichment analysis. The gene abundances were
calculated and normalized to RPKM (Reads Per kb per Million reads). c
The relative expression levels of 9 selected genes at different stages
of floral bud development in M and H. Relative expression levels were
calculated by the 2-ΔΔCt method with actin as a standard. 11-Apr,
21-Apr, and 4-May represent the pistils of flowers when their
development stages were 5, 6, and 8, respectively. However, here, we
mainly discussed stage 6, which is the most important development stage
for sex differentiation in T. sinensis (for details, see ‘Structure of
the flower’)
Discussion
The male flowers of T. sinensis originated from ancestral hermaphroditic
flowers
A variety of hypothesized evolutionary pathways have been proposed for
the evolution of dioecy from hermaphroditism, and the production of
only bisexual flowers is the ancestral condition in angiosperms
[[218]3, [219]7, [220]16, [221]78–[222]83]. Theory predicts that at
least two mutations are necessary to evolve separate sexes [[223]11].
One must result in male sterility while the second must result in
female sterility [[224]27, [225]29]. Therefore, we examined whether the
androdioecy of T. sinensis resulted from transformation from female
flowers to restore male function in dioecious plants or was caused by
pistil sterility in bisexual flowers.
In T. sinensis, male flowers in each natural population possess the
characteristic of a vestigial pistil, and the receptive stigma is
covered by papillary cells. Although the degree of sterility of pistils
is different, the males have the potential to form ovules. In addition,
in natural populations, no pollen limitation on the reproduction of T.
sinensis [[226]84, [227]85]. Therefore, female flowers evolved from
female flowers are unlikely to restore their male function to ensure
reproduction [[228]15]. Obviously, it is also unreasonable for MF to
evolve a pistil if androdioecy arose from dioecy. In addition, as in
Osmanthus fragrans L., the MF has a vestigial pistil, also suggesting
that males evolved from female-sterile mutants among hermaphrodites
[[229]13].
The breeding system and other traits mapped on the phylogenetic tree of
Fraxinus L. (Oleaceae) show that dioecy has three separate origins; in
one instance, dioecy evolved from hermaphroditism via androdioecy
following the transition from insect to wind pollination [[230]45,
[231]86]. Dioecy may also evolve in conjunction with the evolution of
wind pollination in some previously insect-pollinated species [[232]87,
[233]88]. T. sinensis is a wind-pollinated and insect-pollinated
species [[234]52, [235]53], and this mixed pollination pattern may be a
possible evolutionary pathway from hermaphroditism to dioecy [[236]55].
A trend in tectum sculpture character state change from perforate to
reticulate appears to be common in angiosperms [[237]66, [238]67,
[239]89]. In T. sinensis, the tectum sculpture is reticulate in MF
pollen grains, while that of HF is perforate. Additionally, globose
pollen grains have the largest volume-to-surface area ratio, but with
minimal regulatory function, which is a relatively primitive character
[[240]68]. Olive pollen grains can be transformed from the globose
pollen grains of HF. Evolutionarily, this result also implied that
hermaphroditism is more ancestral [[241]69].
Key regulatory factors of pistil development may be associated with T.
sinensis sex differentiation
Admittedly, ovule primordia arise from the placenta; however, both
ovules and placenta are completely absent in most MFs of T. sinensis.
Transcriptome analyses showed that sixteen genes related to carpel and
ovule development are expressed differentially at the key stage of MF
pistil abortion and were verified by qRT-PCR. The molecular mechanisms
of ovule development in Arabidopsis, cotton, and rice are relatively
well known [[242]24–[243]26, [244]90]. The C-class gene AG (AGAMOUS)
belongs to the MADS-box family, which determines stamen and carpel
identity and plays a role in ovule development [[245]31, [246]44]. The
SHP1 (SHATTERPROOF 1), SHP2 (SHATTERPROOF 2), and STK (SEEDSTICK) genes
together with the SEP (SEPALLATA) genes have overlapping expression
patterns with AG in an AG-independent manner [[247]31]. ANT
(AINTEGUMENTA), HLL (HUELLENLOS), SIN2 (SHORT INTEGUMENT 2), INNER NO
OUTER (INO), and SUPERMAN (SUP) as regulators of ovule outgrowth
[[248]31, [249]76, [250]90, [251]91]. AG is necessary for the
determination and termination of floral meristem activity. In the ag
mutant, delayed floral meristem activity results in the stamens and
carpels being replaced by flowers [[252]61]. WUSCHEL (WUS), which
specifies stem cell identity, and AG play crucial roles in regulating
the initiation or termination of floral meristems [[253]92, [254]93].
Previous studies have indicated that SPL/NOZZLE is required for the
initiation of sporogenesis and plays a central role in nucellus
formation by antagonizing the expression of BELL, ANT and INNER NO
OUTER [[255]94]. Here, four genes are shown to be differentially
expressed in MF relative to HF at stage 6. Moreover, SPL represses the
expression of ANT and INO to control ovule primordium formation
[[256]95]. SPL was downregulated in the key stage of pistil abortion in
MF relative to HF, which indicated the role of SPL in promoting ovule
development. However, during ovule development, BEL1 and SPL are
required for cytokinin and auxin signalling for the correct patterning
of the ovule [[257]96].
In addition, the STM gene is considered to play a pivotal role in
sustaining stem cell function in the floral apical meristem (STM is
required for maintenance of WUS expression but independent of WUS
expression) and controls KNOX gene expression independently of the
transcriptional repressor AS1 [[258]97, [259]98]. More precisely, STM
plays a crucial role in preventing meristem cell differentiation by
inducing the production of cytokinins (CK) and the ARR transduction
pathway [[260]31]. The progressive loss of STM causes floral phenotypes
ranging from reduced formation of placental tissues and inhibited
carpel fusion to complete loss of carpel development [[261]99]. In T.
sinensis, STM also showed differential expression, which indicated the
role of STM in promoting ovule development.
In Arabidopsis, the ovule primordia arise from the placentas flanking
the medial ridges, while leunig (lug) aintegumenta (ant) double mutants
lack medial tissues, and the double mutant seu-3 ant-1 results in a
complete loss of ovule initiation, caused by severe defects in early
gynoecium development [[262]100]. CRABS CLAW (CRC) is a member of the
YABBY gene family, whose expression is largely limited to carpels and
nectaries, and plants expressing 35SCaMV::CRC produce variable flowers
with partial transformations of sepals into carpels and reductions in
flower organ number. Carpels are short and solid (ovule is completely
absent) and composed primarily of stylar tissues. The morphology of the
pistillode of MF in T. sinensis was similar to that of the Arabidopsis
35SCaMV::CRC pistil [[263]71].
Ethylene may act as an upstream factor for female sterility in male flowers
Ethylene, as an upstream factor for pistil or ovule formation, has been
verified in many species; for instance, tobacco, pomegranate, and
cucumber [[264]42, [265]101, [266]102].
The ETR1 protein acts as negative regulator of ethylene responses
[[267]103], has the function of promoting cell elongation. In T.
sinensis, ETR1 (Unigene0020715) was up-regulated in MF relative to HF.
ERF1 and ERF2, as ethylene response signal factors, were differentially
expressed between MF and HF. The down-regulation of ERF1/2 and
up-regulation of ETR1 may be related to the flower development. ACO and
ACS (ACS and ACO are rate-limiting enzyme) are a pivotal enzyme in
ethylene biosynthesis, they are showed significant expression levels in
pistil formation ant stage 6. Jointly, our findings indicated that
ethylene may act as one of the upstream regulators influencing the
development of pistil. Furthermore, ANT, an APETALAP2-like Gene, plays
fundamental roles in ovule formation and megasporogenesis [[268]38,
[269]104, [270]105], the down-regulation of ANT may be responded to
ethylene. However, the expression of ANT and INO (INNER NO OUTER) may
be regulated by AG and SPL.
Auxin and Cytokinin are required for ovule Primordia formation
Auxin maxima are fundamental for the formation of primordia, and auxin
action has been well described for lateral roots and primordial flowers
[[271]44, [272]106, [273]107]. The directionality of the auxin flux
depends principally on the polar localization of the PIN proteins. In
Arabidopsis, there are eight PIN proteins (PIN1–8), from which only
PIN1 and PIN3 are expressed in the pistil and ovules. For example,
polar auxin transport is mostly mediated by the PINFORMED1 (PIN1)
efflux carriers [[274]83]. Based on the phenotypes of pin mutant,
gynoecia display a range of phenotypes from almost normal structures
with two valves, style and stigma, to stalk-like gynoecia with an
elongated gynophore topped with a style and stigma [[275]108]. The
protein kinase PINOID (PID) positively regulates auxin efflux, which
can control an early general step in meristem development. In the pid-l
flower, the internal fourth whorl organ is stem-like rather than
carpelloid, although it is capped with a ring of stigmatic papillae.
Furthermore, the morphology of the vestigial pistil of MF T. sinensis
was similar to that of the Arabidopsis pid-l pistil [[276]109].
Similarly, PINOID is required for formation of the stigma and style in
rice [[277]110]. In addition, the crosstalk between CK and auxin plays
a pivotal role in the boundaries in the inflorescence meristem and
ovule primordia formation [[278]44]. Simultaneously, the Arabidopsis
Histidine Phosphotransfer Protein 6 (AHP6), an inhibitor, mediates an
auxin-cytokinin crosstalk that regulates the timing of organ initiation
at the shoot apical meristem. In this study, genes that were expressed
differently in MF relative to HF also supported this point.
Conclusion
In total, 29,585 out of 52,945 unigenes were annotated, and the KEGG
enrichment analysis will provide valuable information on the sex
differentiation of T. sinensis. First, ethylene is an upstream factor
in T. sinensis gynoecium development, and signal transduction of auxin
and cytokinin is fundamental to the initiation of pistil primordia;
MADS-box genes were also shown to be involved in the determination of
the fourth-whorl flower organs. In addition, auxin signalling pathways
were thought to be involved in pollen abortion (partial pollen abortion
in T. sinensis HF populations; Fig. [279]5) and were important for
megagametophyte development (male mutant-Type III in MF; Fig. [280]1).
Overall, a valuable new model of sex differentiation in T. sinensis was
summarized, and further study of the functions of sex-related genes
will be helpful in elucidating the sex-determination mechanism of the
androdioecious plant T. sinensis (Fig. [281]11).
Fig. 11.
[282]Fig. 11
[283]Open in a new tab
A proposed model of ethylene, auxin, identity genes, and CUC gene
family-dependent pathways for gynoecium development. All genes
displayed in this figure correspond to DEGs in this study; other
related genes without differences are not shown here. Different colours
of rectangles represent different pathways or factors; red: ethylene
synthesis and metabolism; blue: polar auxin transport; light blue: CUC
gene family; grey: MADS-box C genes. Arrows and bars indicate positive
and negative regulatory interactions, respectively. These regulatory
networks come from the relationships determined in Arabidopsis
[[284]31, [285]44, [286]60, [287]105, [288]118–[289]122]
Methods
Flower sample collection
T. sinensis in natural populations has 2 sexual types: male individuals
(MF, male flower) and hermaphroditic individuals (HF, hermaphroditic
flower) with winter fruits. Two types of flowers of T. sinensis were
sampled from the natural geographical distribution in Shennongjia
National Nature Reserve (31°38′N, 110°25′E; Hubei, China) and the
voucher specimens were deposited archived in the herbarium of
Structural Botany Laboratory (SBL), Northwest University. Prof. Liu
Wenzhe undertook the formal identification of the samples. Different
growth stages of flowers were collected for our study. The flowers used
for microscopic observation were fixed immediately in 2.5%
glutaraldehyde, and those used for transcriptome analysis were stored
in liquid nitrogen (LN) at − 80 °C.
A total of 18 samples were sequenced, comprising males and
hermaphrodites in three different development stages (with three
biological replicates for each stage). To describe the samples clearly,
we used M1, M2, and M3 to represent MF buds when their development was
at stage 5, stage 6, and stage 8, respectively (Fig. [290]1a-e). Stages
are numbered as described by Sanders et al. (1999) [[291]111]. Three
stages of samples were collected mainly according to the results of
anatomy (Figs. [292]3d, e and [293]4); these stages represent the key
stages that implied the formation of sexual dimorphism in T. sinensis.
Similarly, H1, H2, and H3 represented the HF buds when their
development was at stage 5, stage 6, and stage 8, respectively.
Photographs are included showing the general developmental stages of
the flowers (Fig. [294]1a-e).
RNA sequencing library construction, Illumina sequencing, and data processing
Total RNA (40 μg) was isolated from the M1, M2, M3 and H1, H2, and H3
samples using the TRIzol Kit (Promega, Beijing, China) and treated with
RNase-free DNase (Takara Biotech Incorporation, Otsu, Japan). RNA
quantity and purity were assessed using a Nanodrop 1000
spectrophotometer (Thermo Scientific, Waltham, MA, USA). RNA integrity
(> 8.0) was determined using a Bioanalyzer 2100 (Agilent) with an RNA
6000 Nanochip (Agilent). Methods of enrich Mrna, remove rRNA, RNA
fragment, random haxamer primed cDNA synthesis, size selection and PCR
amplification were performed as previously described [[295]43,
[296]112, [297]113], and sequenced using Illumina HiSeqTM 4000 (Gene
Denovo Biotechnology Co., Guangzhou, China).
De novo assembly and functional annotation
Raw data were processed with Perl scripts to ensure the quality of the
data used in further analyses. The adopted filtering criteria are as
Zhang et al. (2018) [[298]55, [299]114]. The obtained clean data after
filtering were used for statistical analyses of quality, including Q30,
data quantity and quality, and base content statistics. The software
Trinity-v2.6.5 was used for de novo assembly [[300]55]. TransDecoder
identifies candidate coding regions within transcript sequences, such
as those generated by de novo RNA-Seq transcripts assembled by Trinity.
Functional annotation of transcript assemblies was performed using the
included Trinotate pipeline ([301]http://trinotate.sourceforge.net/) to
identify open reading frames and assign best hits to UniProtKB (1e-03),
PFAM-A (1e-03), Gene Ontology (GO) and eggNOG categories.
Screening of differentially expressed genes
To identify DEGs across samples or groups, the edgeR package
([302]http://www.r-project.org/) was used. We identified genes with
greater than a two-fold change and a false discovery rate (FDR) < 0.05
in a comparison as significant DEGs. DEGs were then subjected to an
enrichment analysis of GO functions and KEGG pathways. Read count for
each gene in each sample was determined by HTSeq (v0.6.0;
[303]http://www-huber.embl.de/users/anders/HTSeq/doc/overview.html),
and RPKM (reads per kilobase per million mapped reads) was then
calculated to estimate the expression level of the genes in each sample
[[304]115]. The RPKM method is able to eliminate the influence of
different gene lengths and sequencing data amounts on the calculation
of gene expression. Therefore, the calculated gene expression can be
directly used for comparing differences in gene expression among
samples.
Differential expression analysis of the three stages in MF and HF was
performed using DESeq [[305]116]. DESeq2 v1.4.5 was used for
differential gene expression analysis, and a p-value could be assigned
to each gene and adjusted by Benjamini and Hochberg’s approach for
controlling the false discovery rate. Genes with q ≤ 0.05 and
|log2_ratio| ≥ 1 were identified as DEGs. Gene expression patterns in
all samples were compared with expression profiles for known
Arabidopsis, Zea and Oryza carpel development or sex-determination
genes [[306]31, [307]58–[308]60, [309]70, [310]78, [311]106, [312]115].
GO and pathway enrichment analysis of DEGs
The GO (Gene Ontology; [313]http://geneontology.org/) enrichment of
DEGs was implemented by the hypergeometric test, in which a p-value is
calculated and adjusted as the q-value, and the data background is
genes in the whole genome. GO terms with q < 0.05 were considered to be
significantly enriched. GO enrichment analysis could reveal the
biological functions of the DEGs [[314]112], and we focused on genes
that play roles in certain biological functions, i.e., those for which
pathway enrichment analysis (KEGG, Kyoto Encyclopedia of Genes and
Genomes) was mainly concerned with metabolic pathways or signal
transduction pathways. Pathway analyses were performed with the KEGG
pathway database using BLASTx with an E-value < 1 e-5 [[315]117]. The
calculating formula is the same as that in GO analysis.
[MATH: P=1−∑i=0
m−1Mi
N−M
n−i
Nn
:MATH]
Here N is the number of all genes that with KEGG annotation, n is the
number of DEGs in N, M is the number of all genes annotated to specific
pathways, and m is number of DEGs in M. The calculated p-value was gone
through FDR Correction, taking FDR ≤ 0.05 as a threshold. Pathways
meeting this condition were defined as significantly enriched pathways
in DEGs.
Real-time RT-PCR analysis
qRT-PCR was performed not only to analyze the expression levels of DEGs
but also to verify the RNA-Seq analysis, which included ten genes
(CTR1, ACS, ACO, ACO1, ETR1/2, ERS1, EIN2/3/4) related to ethylene
biosynthesis and signal transduction, nine genes involved in auxin and
cytokinin signal transduction, and sixteen genes related to carpel and
ovule formation. Primer sequences were designed using Primer premier
5.0, and Primer sequences and β - actin (as an internal control to
normalize the expression levels of genes) sequence of internal
reference gene are listed in Additional file [316]1. RNA was extracted
from T. sinensis of three developmental stages, and RNA extraction,
cDNA synthesis, and Real-time PCR method as described in related
references [[317]116, [318]123].