Accepted Manuscript Genome-wide transcriptome profiling provides insights into panicle development of rice (Oryza sativa L.)
Shanwen Ke, Xin-Jiang Liu, Xin Luan, Weifeng Yang, Haitao Zhu, Guifu Liu, Guiquan Zhang, Shaokui Wang PII: DOI: Reference:
S0378-1119(18)30765-0 doi:10.1016/j.gene.2018.06.105 GENE 43038
To appear in:
Gene
Received date: Revised date: Accepted date:
6 February 2018 26 April 2018 28 June 2018
Please cite this article as: Shanwen Ke, Xin-Jiang Liu, Xin Luan, Weifeng Yang, Haitao Zhu, Guifu Liu, Guiquan Zhang, Shaokui Wang , Genome-wide transcriptome profiling provides insights into panicle development of rice (Oryza sativa L.). Gene (2018), doi:10.1016/j.gene.2018.06.105
This is a PDF file of an unedited manuscript that has been accepted for publication. As a service to our customers we are providing this early version of the manuscript. The manuscript will undergo copyediting, typesetting, and review of the resulting proof before it is published in its final form. Please note that during the production process errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain.
ACCEPTED MANUSCRIPT Title: Genome-wide transcriptome profiling provides insights into panicle development of rice (Oryza sativa L.)
Author names: Shanwen Ke, Xin-Jiang Liu, Xin Luan, Weifeng Yang, Haitao Zhu, Guifu Liu,
Affiliation: Key
Laboratory
for
Conservation
and
Utilization
SC
State
RI
PT
Guiquan Zhang*, Shaokui Wang*
of
Subtropical
Agro-Bioresources, Guangdong Key Laboratory of Plant Molecular Breeding,
NU
College of Agriculture, South China Agricultural University, Guangzhou 510642,
MA
China.
Corresponding author: Corresponding author.
D
*
PT E
E-mail address:
[email protected] (G. Zhang),
[email protected] (S.
AC
CE
Wang).
1
ACCEPTED MANUSCRIPT Abstract Panicle architecture is an important component of agronomic trait in rice, which is also a key ingredient that could influence yield and quality of rice. In the panicle growth and development process, there are a series of complicated molecular and cellular events which are regulated by many interlinking genes. In this study, to
PT
explore the potential mechanism and identify genes and pathways involved in the formation of rice panicle, we compared the transcriptional profile of rice panicles
RI
(NIL-GW8 and NIL-gw8Amol) at three different stages of panicle development: In5 (formation of higher-order branches), In6 (differentiation of glumes) and In7
SC
(differentiation of floral organs). A range of 40.5 to 54.1 million clean reads was aligned to 31,209 genes in our RNA-Seq analysis. In addition, we investigated
NU
transcriptomic changes between the two rice lines during different stages. A total of 726, 1,121 and 2,584 differentially expressed genes (DEGs) were identified at stages
MA
1, 2 and 3, respectively. Based on an impact analysis of the DEGs, we hypothesize that MADS-box gene family, cytochrome P450 (CYP) and pentatricopeptide repeat
D
(PPR) protein and various transcription factors may be involved in regulation of
PT E
panicle development. Further, we also explored the functional properties of DEGs by Gene Ontology analysis, and the results showed that different numbers of DEGs genes were associated with 53 GO groups. In KEGG pathway enrichment analysis, many
CE
DEGs related to biosynthesis of secondary metabolites and plant hormone signal transduction, suggesting their important roles during panicle development. This study
AC
provides the first examination of changes in gene expression between different panicle development stages in rice. Our results of transcriptomic characterization provide important information to elucidate the complex molecular and cellular events about the panicle formation in rice or other cereal crops. Also, the findings will be helpful for the further identification of the genes related to panicle development.
Keywords Rice; Panicle; RNA-Seq; Transcriptome; GW8 2
ACCEPTED MANUSCRIPT Highlights Tanscriptome analysis was used to study the panicle differentiation of rice. Compare gene expression profile of two rice near isogenic lines during different panicle development stages. Obtain a list of candidate genes related to the formation of rice panicle involved in
RI
PT
OsSPL16.
Abbreviations
SC
COG, Cluster of orthologous groups of proteins; CYP, Cytochrome P450 monooxygenase; DEGs, Differentially expressed genes; FDR, False discovery rate;
NU
FPKM, Fragments per kilobase of transcript per million fragments mapped; GO, Gene ontology; KEGG, Kyoto Encyclopedia of Genes and Genomes; NIL, Near isogenic
MA
lines; Nr, Non-redundant protein database; PPR, pentatricopeptide repeat protein; qRT-PCR, Quantitative real-time PCR validation; QTL, Quantitative trait locus; SBP,
AC
CE
PT E
D
Squamosa promoter binding protein.
3
ACCEPTED MANUSCRIPT 1. Introduction Rice (Oryza sativa L.) is one of the most important cereal crops in the world. Benefit from the green revolution and the application of heterosis, the production of rice has achieved the double increase in the last 50 years, and even reached the increase of more than three times in some countries and regions (Xing and Zhang,
PT
2010). But it is predicted that the world food production must rise by 50% to meet the demand of the increasing population in 2050 (Ashikari et al., 2005). Panicle
RI
architecture is an important component of agronomic trait in rice, which determines yield of rice. The formation of panicle is a systemic process including the
SC
development of axillary meristems, the establishment of inflorescence and the formation of grains.
NU
Although researchers have studied the formation of grains for the whole last century, but the study of the molecular mechanism of grain size determination was
MA
only about 20 years, and the results mainly focus on model plants (Sundaresan, 2005). Grain affects not only the production but the appearance. Improvement of grain shape
D
is an effective way to increase yield and quality for rice. With the acceleration of
PT E
mutant research and the use of the advanced biological technologies, as well as the large-scale development of molecular markers, many quality trait genes associated with grain shape have been identified and cloned (Luo et al., 2013; Si et al., 2016).
CE
Some of them are negative regulators of grain length. For example, GS3, one major QTL for grain length, encodes a putative membrane protein that regulates grain
AC
growth by limiting cell proliferation (Fan et al., 2006; Mao et al., 2010). GL3.1 regulates grain size and yield by regulating the cyclin-T1;3 that influences the replication of chromosomes and promote cell division, leading to shorter rice grain (Qi et al., 2012). Studies have shown that GL3.1 influences the grain size and yield mainly through changing cell number of glume to regulates its length, on the other hand, it also affects the filling rate of rice (Zhang et al., 2012). GS2 functions as a positive regulator of grain length through regulating cell division and subsequent elongation (Duan et al., 2016). In addition, a few other genes have been studied as the regulators of grain width. GW2 is the first cloned gene which controls the grain width 4
ACCEPTED MANUSCRIPT and weight; it may be involved in the degradation of protein to promote the process of cell division, and then influence grain weight and yield by regulating the rice glume shell size (Song et al., 2007). The gene, qSW5, was discovered through gene mapping analysis; it regulates grain width by influencing the clever shell cell number (Shomura et al., 2008). Its allele GW5 was also cloned; it affect the grain width likely through
PT
participating in protein degradation pathways (Weng et al., 2008). GS3, GL3.1, GS2, GW2, qSW5 and GW5 reflect the negative relation of sequence. However,
RI
GS5 positively regulates rice grain weight; overexpression of GS5 increases grain size by promoting cell proliferation (Li et al., 2011). Another major QTL for rice weight
SC
and length (TGW6) promotes IAA-glucose hydrolyzed into free IAA and glucose, and participates in the carbohydrate transportation between source and sink (Ishimaru et
NU
al., 2013).
On the other hand, the number of branch and spikelets per panicle is one of the
MA
key factors affecting rice yield. The spikelet is the basal unit of inflorescence in grasses, and it is crucial for reproductive success and cereal yield. Although the
D
growth of plant organs, like panicles, often affected by the external environment
PT E
signals, there is still plenty of evidence to suggest that the morphology of plant organs is mainly affected by the internal signals which is controlled by many genes (Liu et al., 2002). Gn1a controls branch differentiation; the decreased expression level of Gn1a
CE
leads to the accumulation of cytokinins in inflorescence meristem, resulting in the increase of spikelets number per panicle (Ashikari et al., 2005). DEP1 is postulated to
AC
function as a negative regulator of branch differentiation, when in function-defect of DEP1, branch number and spikelets number increases significantly (Huang et al., 2009). PROG1 is also a down-regulator for the number of branch and spikelets, and it is mainly expressed in axillary bud organization (Tan et al., 2008). A few positive regulatory factors have also been found. For instance, FZP plays an important role in the forming process of spikelet meristems (SMs) from spikelet-pair meristems (SPMs) which may participate in the development of branch meristems (Komatsu et al., 2003). Besides, there are many positive regulatory factors in the process of the development of young panicles, like APO1, LAX1, LAX2, RCN1, RCN2 (Nakagawa et al., 2002; 5
ACCEPTED MANUSCRIPT Terao et al., 2010; Tabuchi et al., 2011). In rice, some SBP (Squamosa promoter Binding Protein) protein family members are tightly tied to the development of branch. Such as OsSPL14, its high-level expression leads to the increase of branch number and grain number (Jiao et al., 2010). OsLG1, encoding a SBP protein, regulates the panicle development by influencing the
PT
parenchyma cells of the gaps between the cobs and branches (Ishii et al., 2013). In our previous work, we have cloned a new gene GW8 which belongs to the family of SBP
RI
transcription factors. It positively regulates cell proliferation, improving not only the 1000-grain weight but also the rice quality (Wang et al., 2012). Based on the analysis
SC
of the single segment substitution lines (SSSLs), we found a GW8 allele from an indica rice (Amol3), named gw8Amol. Compared with the NIL-GW8, the grains
NU
showed high quality and slender shape without yield difference in NIL-gw8Amol because of the increase in the number of branch and grains per panicle (Fig. 1) (Wang
MA
et al., 2012). Our previous studies have shown that the expression of Gn1a, Lax1, FZP, OsCLV1, OsCLV2, OsCLV3, RCN1, RCN2, DEP1 and OsSPL14 have no significant
D
difference in panicle at the early stage of inflorescence development between
PT E
NIL-GW8 and NIL-gw8Amol plants (Wang et al., 2012). We suspect that GW8 regulates the panicle development and differentiation through some unknown metabolic pathways. It is interesting to identify and characterize the regulators of GW8 (ROGs)
CE
which is an essential gene affecting the branch differentiation and regulation of the grain number per panicle in rice.
AC
The completion of the sequencing of rice genome and the development of high-density molecular markers provide a new opportunity for the research of the genome of rice and other crops (Devos, 2005; Ashikari and Matsuoka, 2006; Xing and Zhang, 2010). In recent years, with the rapid development of next-generation sequencing technology, mRNA sequencing (RNA-seq) technology has become a useful means for the study of transcription profiles (Wang et al., 2010). In this study, a wide grain rice line (NIL-GW8) and a slender grain rice line (NIL-gw8Amol) were used to compare gene expression changes in panicles by deep RNA sequencing. We did the transcriptome analysis of three developmental stages of panicles in these two rice 6
ACCEPTED MANUSCRIPT lines. The identification and functional validation of ROGs may help us to better understand the working mechanism of GW8 in the regulation of cell proliferation in rice and eventually create new varieties with higher yield and better quality.
2. Materials and methods
PT
2.1. Plant materials Two indica rice (O. sativa L. ssp. indica) near-isogenic lines (NIL) NIL-GW8 and
RI
NIL-gw8Amol in the HJX74 genetic background were used for the transcriptome sequencing analysis in this research. NIL-gw8Amol contains a very short chromosome
SC
segment inherited from Amol3 (Wang et al., 2012), carrying an allele of GW8 (gw8Amol) which has a mutational OsmiR156 binding site. And NIL-gw8Amol produces
NU
slender grains and more branches and spikelets per panicle than NIL-GW8. The two rice lines were grown in Teaching and Experimental Farm in South China
differentiation period
of
MA
Agricultural University (113°36'E, 23°16'N). We sampled the panicles in the primary branch and secondary branch.
A total of
six
D
inflorescence samples (three stages: In5 (formation of higher-order branches), In6
PT E
(differentiation of glumes) and In7 (differentiation of floral organs) from the two NILs) were collected for RNA-seq according to Ikeda et al. (Ikeda et al., 2004). These collected samples were immediately frozen in liquid nitrogen and stored at -80℃
CE
until the total RNA extraction. The panicle organs which used in this study are very small, and the differentiation periods are relatively shorter with one to three days
AC
(Ikeda et al., 2004). For biological replicates, every sample of each developmental stage was a mixture of fifty panicles from at least ten individual plants. The growth and development status of these sampled panicles in a specific development period were highly consistent.
2.2. RNA isolation and cDNA library construction Total RNA was extracted from panicles and spikelets using TRIzol® RNA Purification Kit (Invitrogen) according to the manufacturer’s protocol, and electrophoresed in 1.0 % agarose gel. Nanodrop 2000 Spectrophotometer, Qubit® 2.0 7
ACCEPTED MANUSCRIPT Spectrophotometer
(Life
Technologies,
Carlsbad,
CA)
and
Aglient
2100
Bioanalyzer® (Agilent Technologies, Palo Alto, Calif.) were used to measure and test the purity, concentration and integrity of RNA samples, respectively. The enrichment of the poly (A)-containing mRNA was performed by using Oligo (dT) magnetic beads. The mRNA was fragmentated and used to synthesize double-stranded cDNA. The Agencourt AMPure XP beads (Beckmann Coul-ter, Pasadena, CA) were used to
PT
purify dsDNA. After that, the ends of purified ds-cDNA were repaired. Next, the ‘A’
RI
tail and paired-end adapters for sequencing were ligated to the end of cDNA. In the end, different sizes of fragments were selected using AMPure XP beads for PCR
SC
amplification to create the final RNA-seq library. When the cDNA library construction is finished, the concentration and insert size of the libraries were tested
NU
using Qubit2.0 and Agilent 2100, respectively. To ensure the quality of library, using the qPCR accurately quantify the effective concentrations of the libraries above. After
MA
complete the testing, the libraries were sequenced on Illumina HiSeq2500 Sequencer. These sequence data were deposited in an NCBI database and are accessible via
PT E
D
BioProject ID: PRJNA439157.
2.3. Quantitative real-time PCR validation The validation of the RNA-seq results was performed through ten genes with
CE
different expression profiles by quantitative RT-PCR. The read counts and the primer sequences of the ten genes were provided in Table S6 and Table S7, respectively.
AC
Total RNA was extracted references the RNA-seq experiments and first-strand cDNA was generated using a TransScript One-Step gDNA Removal and Cdna Synthesis SuperMix kit (TransGen Biotech). Quantitative RT-PCRs were carried out using the ABI Step One Plus Real-Time PCR System (Applied Biosystems) with the SYBR Premix Ex Taq™ RT-PCR kit (Takara). The rice gene Actin1 (LOC_Os03g50885) was used to standardize the relative expression levels.
2.4. Quality analysis and mapping reads to the reference genome The software Phred was used to evaluate the probability of incorrect base calling 8
ACCEPTED MANUSCRIPT according to phred quality score (Ewing et al., 1998). Before mapping reads to the reference genome, read quality control was conducted as follows: first, the sequencing adapters and primer sequences were removed; then, all the low-quality raw reads were removed to get the clean reads. Then all the clean reads were mapped to the rice reference
genome
(Os-Nipponbare-Reference-IRGSP-1.0,
PT
ftp://ftp.ensemblgenomes.org/ pub/plants/release-24/fasta/oryza_sativa/) (Kim et al.,
RI
2013) by software TopHat V2 (Langmead et al., 2009).
2.5. Analysis of gene expression
SC
Analysis of gene expression was performed based on Negative Binomial Distribution. The Cuffdiff components of Cufflinks software were used to do
NU
quantitative analysis of transcripts and calculate genes expression level according to reads mapped to specific genes. FPKM (Fragments Per Kilobase of transcript per
MA
Million fragments mapped) value was directly used for assessing the expression
D
abundance of each gene with high sensitivity across six orders (10-2 to 104).
PT E
2.6. Identification of differentially expressed genes We applied EBSeq package to identify differentially expressed genes (DEGs) among different panicle developmental stages of the two rice varieties. Genes with
CE
Fold Change ≥ 2 and FDR (false discovery rate) < 0.01 were selected as DEGs. The log2FC (fold changes) was used to indicate the significance of expression difference
AC
of each gene. The test corrections of P-values were performed using the Benjamini-Hochber False Discovery Rate method. The numbers of DEGs between sample pairs were counted and plotted in a Venn diagram. For Gene ontology (GO) enrichment analysis, BLAST2GO Ver. 2.7.2 was used for the detected DEGs according to biological process, molecular function and cellular component ontologies. GO terms with a FDR cut-off of 0.01 were considered as significantly enriched. Then the DEGs were subjected to search against COG (Cluster of Orthologous Groups of proteins) database to predict and classify their possible functions. After that, the pathway enrichment and annotation analysis of DEGs were 9
ACCEPTED MANUSCRIPT carried out by comparing to the KEGG(Kyoto Encyclopedia of Genes and Genomes) database using Cytoscape software. Finally, the DEGs significantly enriched metabolic pathways or signal transduction pathways were identified with ClueGO.
3. Results
PT
3.1. Illumina sequencing and alignment against the reference genome The sequencing of mRNA generated six transcriptome datasets from two NILs
RI
(NIL-GW8 and NIL-gw8Amol) at three different panicle developmental stages [W1: stage 1 (In5) in NIL-GW8; W2: stage 2 (In6) in NIL-GW8; W3: stage 3 (In7) in
SC
NIL-GW8; S1: stage 1 (In5) in NIL-gw8Amol; S2: stage 2 (In6) in NIL-gw8Amol; S3: stage 3 (In7) in NIL-gw8Amol]. W1, W2, W3, S1, S2 and S3 libraries generated around
NU
46.5, 53.1, 54.1, 40.5, 40.8 and 46.1 million high-quality clean reads, respectively, each of which was of a Q30 (an error probability of 1%) value of at least 89.48%
MA
(Table 1). All the reads were mapped to the rice reference genome by TopHat Ver. 2.0.11. A total of 227.1 million reads with a mean mapping rate of 80.88 % (ranging
D
from 79.88 to 82.35% for the six libraries) were mapped to the genome. Among of the
PT E
mapped reads, 77.55-80.01% were matched to unique genomic locations of the reference genome. 17.65-20.12% of the total reads were unmapped (Table 1).
CE
3.2. Global analysis of gene expression In the study, a total of 5,971 (W1), 5,899 (W2), 6,206 (W3), 5,737 (S1), 5,606 (S2)
AC
and 6,147 (S3) genes with different length (most of them were no less than 200 bp) were obtained from the six samples, respectively. Of the total 31,209 genes, 4,693 (15.03 %) genes longer than 3,000 bp accounts for the biggest proportion, and only 28 (0.08 %) genes with shorter than 200 bp (Fig. 2). More than half of those genes (51.45 %) were longer than 1,500 bp (Fig. 2). To get a better visualization of the gene expression in the six libraries, the Circos diagram was generated using the Circos software package. In each library, fragments per kilobase of transcript per million fragments mapped (FPKM) was used as the index of gene expression. As shown in Fig. 3, the expression of the genes had 10
ACCEPTED MANUSCRIPT apparent differences among different chromosomes. The genes with high expression were distributed mainly over Chr1, Chr2 and Chr3, while the genes on Chr11 and Chr12 showed low expression compared with other chromosomes in the six libraries. We obtained a large number of genes after removing partially overlapping sequences, and was able to get plenty of data for the analysis of rice panicle
characters.
Their
expression
in
the
three
panicle
PT
development. In the six libraries, the numbers of expressed genes displayed similar developmental
stages
RI
(2.0-millimetre-long panicle, 3.0-millimetre-long spikelet and 5.0-millimetre-long spikelet, respectively) were summarized in several Venn diagrams. Among these
SC
genes in NIL-GW8, besides the 24,133 genes expressed at all three panicle developmental stages, 563 genes were expressed in stage1 (W1) and stage2 (W2)
NU
simultaneously, 583 genes co-expressed in stage 2 and stage 3 (W3), and 563 genes co-expressed in stage1 and stage3. However, there were 313, 347 and 1,907genes
MA
specifically expressed in stage1, stage2 and stage3, respectively (Fig. 4). In addition, it also had some uniquely expressed genes in the two types of panicles from
D
NIL-GW8 and NIL-gw8Amol. The number of variety-specifically expressed genes was
PT E
560 (NIL-GW8) and 1,778 (NIL-gw8Amol) at stage 1, 710 (NIL-GW8) and 1,446 (NIL-gw8Amol) at stage 2, 960 (NIL-GW8) and 1,181 (NIL-gw8Amol) at stage 3, respectively (Fig. 4). These uniquely expressed genes may play an important role in
CE
the formation process of panicles.
AC
3.3. Identification of DEGs The growth and development of rice is a complicated procedure,there are many differentially expressed genes were identified between a rice line and its NIL (Moumeni et al, 2015; Tariq et al, 2018). Also, many genes involved in the development of rice panicle, especially the inflorescence and spikelet morphogenesis (Meng et al, 2017). The development period we chose involves the critical period of panicle morphogenesis, rapid changes in gene expression in rapidly differentiated organization (Ikeda et al., 2004; Itoh et al., 2005). On the other hand, every plant sample was from different individual plants, the individual differences were enriched 11
ACCEPTED MANUSCRIPT and amplified. So, we obtained a large number of differentially expressed genes (DEGs) in this study. Further, the differences in gene expression among the different developmental stages in NIL-GW8 and NIL-gw8Amol were analyzed. To identify DEGs, the false discovery rate (FDR) threshold was set to 0.01 and a threshold of 3-fold expression
PT
change was used as the selection criteria. All of the differentially expressed genes between two samples during different panicle developmental stages in NIL-GW8 and
RI
NIL-gw8Amol were listed in Table S1. As shown in Fig. 5, we found that a total of 858 (142 up-regulated genes, 716 down-regulated genes) DEGs in NIL-GW8 and 253 (79
SC
up-regulated genes, 174 down-regulated genes) DEGs in NIL-gw8Amol between stage 1 and stage 2, indicating that a lot more down-regulated genes were identified in
NU
NIL-GW8 than NIL-gw8Amol between stages 1 and 2. Between stage 2 and stage 3, 4,265 (2,835 up-regulated genes, 1,430 down-regulated genes) and 4,295 (2,933
MA
up-regulated genes, 1,362 down-regulated genes) DEGs were found in NIL-GW8 and NIL-gw8Amol, respectively. Compared to the early stage of the panicle development,
D
there were more significantly changed genes detected in the later period, suggesting
PT E
that the GW8 gene was specifically involved in the developmental processes at the stage of 2 and 3. There are a large amount of DEGs identified between the same stage of NIL-GW8 and NIL-gw8Amol. 726 (460 up-regulated genes, 266 down-regulated
CE
genes), 1,121 (1,082 up-regulated genes, 319 down-regulated genes) and 2,584 (1,160 up-regulated genes, 1,424 down-regulated genes) DEGs were identified for stages 1, 2,
AC
and 3, respectively (Fig. 5). The numbers of DEGs in stage 2 and stage 3 were remarkably greater than stage 1 indicating the involvement of complex developmental events in the development. Besides, we also noticed that the number of down-regulated genes was higher than up-regulated genes only in the comparison of W3-S3 (stage 3). As shown in Table 2, some DEGs with significant expression difference were detected between NIL-GW8 and NIL-gw8Amol at different stages. Among these genes, some of them were regarded as new genes which could not be mapped to the rice reference genome. To obtain the annotation information of the new genes, these genes 12
ACCEPTED MANUSCRIPT were annotated to Nr, Swiss-Prot, GO, COG, KOG, Pfam and KEGG databases, respectively. Some genes, such as Os02g0127600, which belong to down-regulated genes at stages 1 and 2, showed significantly higher expression at stage3, and some up-regulated genes at stages 1 showed significantly lower expression at stages 2 and 3
PT
(Table 2).
3.4. Genes related to the development of panicle
RI
In the identified DEGs, some of them are important for the development of panicle branches and spikelets. APO2/RFL, playing an important role in controlling and
flower
development
through
interaction
SC
inflorescence
with
APO1
(Ikeda-Kawakatsu et al., 2012), can modulate the differentiation of axillary meristem
NU
by regulating the expression of LAX1 (Deshpande et al., 2015). They all exhibited different expression trends during stage transition from 2 to 3 in NIL-GW8 and
MA
NIL-gw8Amol panicles. Unexpectedly, APO2/RFL, APO1 and LAX1 were all down-regulated with -7.677, -5.396, -6.604 in NIL-GW8 and with -7.782, -5.963,
D
-4.901 in NIL-gw8Amol, respectively. Among them, APO2/RFL was the highest
PT E
down-regulated gene in NIL-GW8. TAW1, a special regulator of rice inflorescence architecture, was highly expressed in shoot apical meristems, inflorescence meristem and branch meristem (Yoshida et al., 2013). Its expression was declined with -6.542
CE
and -6.477 from stage 2 to 3 in NIL-GW8 and NIL-gw8Amol. However, the other genes associated with panicle development, like Gn1a, LAX2, OsREL2, Hd1, DEP1 and
AC
PROG1 (Ashikari et al., 2005; Tan et al., 2008; Huang et al., 2009; Endo-Higashi and Izawa, 2011; Tabuchi et al., 2011; Yoshida et al., 2012), they showed no significant expression difference among different stages.
3.4.1. MADS-box genes Previous research about rice has shown that MADS-box genes involved in the development of panicle branches and grains, like OsMADS34 (Zhang et al., 2016), regulated by SPL14 (Wang et al., 2015). In this study, OsMADS34 was significant differentially expressed between stage 2 and 3 in NIL-GW8 and NIL-gw8Amol (Table 13
ACCEPTED MANUSCRIPT 3). In addition, we found thirteen up-regulated and seven down-regulated MADS-box genes (|log2FC| > 2) in the two rice lines. Interestingly, the elevated level of eleven up-regulated genes was higher in NIL-GW8 than in NIL-gw8Amol. Among them, OsMADS13, OsMADS3, MADS6 and OsMADS58, the coordinative effects of which determine the development of ovule and the floral meristem identities (Li et al., 2011),
PT
were significantly up-regulated in this study. OsMADS58, known to play a vital role in early rice stamen development (Chen et al., 2015), was up-regulated in
RI
NIL-gw8Amol with 3.341, but was not increased obviously in NIL-GW8. Conversely, OsMADS29, a gene closely related to embryo and endosperm development in rice
SC
(Nayar et al., 2013), was up-regulated in NIL-GW8 with 3.396, but was not increased obviously in NIL-gw8Amol. In addition, OsMADS1 and OsMADS17 were up-regulated
NU
and OsMDP1 was down-regulated only in NIL-GW8 from stage 1 to 2. Previous studies have shown that OsMADS1 associated with spikelet development (Wang et al.,
MA
2017), and OsMADS17 is the direct target gene of OsMADS1 (Hu et al., 2015). OsMDP1, has a negative regulatory role in brassinolide (BR) signaling (Duan et al.,
D
2006). We speculate that the different expression patterns of these MADS-box genes
3.4.2. CYP genes
PT E
were mainly because of the different function of GW8 in NIL-GW8 and NIL-gw8Amol.
CE
The expressions of a series of genes named cytochrome P450 (CYP) were differentially expressed in NIL-GW8 and NIL-gw8Amol during stage transition from 2
AC
to 3 (Table 4). CYP comprises one of the largest families in plant enzyme protein, participating in a suite of developmental and physiological processes, including reproduction (Bolwell et al., 1994). In this research task, fourteen CYP genes were up-regulated (log2FC > 3.0) between stages 2 and 3 in the two rice lines. It should be noted that two CYPs(Os03g0332100 and Os05g0482400)were involved in the gibberellin (GA) biosynthesis and signal transduction in rice (Luo et al., 2006; Magome et al., 2013). However, the other gene PLA1 (Os10g0403000), encoding a cytochrome P450 CYP78A11 which acts downstream of GA signaling pathway (Mimura et al., 2012) and participates in the development of primary branch (Miyoshi 14
ACCEPTED MANUSCRIPT et al., 2004), was highly down-regulated in the two rice lines from stages 2 to 3. Furthermore, Os04g0570500 was up-regulated with 5.777 only in NIL-GW8, Os11g0289700 was up-regulated with 5.655 only in NIL- gw8Amol. During stage conversion from 1 to 2, there were no differently expressed CYP genes in NIL-GW8, whereas only one CYP gene (Os01g0227400, log2FC 4.742) was up-regulated in
PT
NIL-gw8Amol. The identification of these differentially expressed cytochrome P450
RI
genes, indicating that they may be correlated with the panicle formation.
3.4.3. PPR genes
SC
Pentatricopeptide repeat (PPR) protein is one of the largest families of proteins which are involved in each step of organelles gene expression in plant (Binder et al.,
NU
2011; Dahan and Mireau, 2013; Shikanai and Fujii, 2013). In Arabidopsis, two PPR genes GPR23 and PPR4 have already been linked to the embryo development (Ding
MA
et al., 2006; Kocabek et al., 2006). In rice, a PPR protein sped1-D may take part in the regulation of peduncle and secondary branch elongation by blocking the function of
D
GID1L2, APO2/RFL and WOX3 (Jiang et al., 2014). Even though the expression of
PT E
sped1-D there was not significantly different in NIL-GW8 and NIL-gw8Amol between different stages, but APO2/RFL and WOX3 were significant down-regulated in the two rice lines from stages 2 to 3. The APO2/RFL was the highest down-regulated
CE
gene in NIL-GW8 with -7.677 during the stage transition from 2 to 3. In this study, seventeen PPR genes were identified as DEGs with log2FC >2.0 between different
AC
stages in the two NILs (Table 5). Some of them were up-regulated in NIL-GW8 and down-regulated in NIL-gw8Amol, like Os01g0589900 and Os07g0244400. Some were down-regulated in NIL-GW8 and up-regulated in NIL-gw8Amol, like Os02g0127600. Besides that, Os08g0402600 was up-regulated from stages 1 to 2, but down-regulated from stages 2 to 3 in NIL-gw8Amol. Moreover, we found that some PPR genes were identified as DEGs only in NIL-GW8 or NIL-gw8Amol, like Os02g0759500, Os02g0750400, Os01g0815900, Os07g0635800, Os02g0170000 and Os07g0260000. Although the physiological function of these PPR genes mentioned above was still unknown, we hypothesized that they may play a role in the panicle cell proliferation 15
ACCEPTED MANUSCRIPT and differentiation in development.
3.5. Gene Ontology (GO) enrichment analysis of DEGs Next, we conducted GO enrichment analysis to classify the functions of the DEGs in the process of panicle development. To explore the relevant biological functions of
PT
the DEGs, the significantly enriched GO terms of DEGs were characterized. The identified DEGs could be categorized into three main GO categories (biological
RI
process, cellular component, and molecular function) which contain 53 functional groups, a range of 1 to 3,428 DEGs in each group includes (Table S2). In NIL-GW8,
SC
the DEGs were significantly enriched in 44, 50 and 48 GO groups in W1-W2, W2-W3 and W1-W3, respectively (Fig. S1). In the three comparison groups above,
NU
cellular process (GO: 0009987), cell part (GO: 0044464) and DNA binding (GO: 0005488) were the most enriched GO groups in biological process, cellular
MA
component, and molecular function, respectively. In addition, we found a high proportion of genes from functional groups of single-organism process (GO:
D
0044699), cell (GO: 005623), catalytic activity (GO: 0003824), metabolic process
PT E
(GO: 0008152), organelle (GO: 0043226), catalytic activity (GO: 0003824), response to stimulus (GO: 0050896) and biological regulation (GO: 0065007). Interestingly,
Fig. S2).
CE
the most enriched GO groups was very similar to NIL-GW8 in NIL-gw8Amol (Table S2;
Moreover, to foster learning about the molecular mechanisms underlying the
AC
formation of spikelet in rice panicle, we put emphasis on GO functional classification analysis for up-regulated and down-regulated DEGs, by comparing the NIL-GW8 and NIL-gw8Amol at different stages (Table S3, Fig. 6). In comparison W1-S1 (stage 1), the 726 DEGs were allocated to 46 GO groups, 460 up-regulated genes and 266 down-regulated genes were classified into 46 and 42 groups, respectively. Among them, cell (GO: 0005623), cell part (GO: 0044464), organelle (GO: 0043226), cellular process (GO: 0042995), single-organism process (GO: 0044699), metabolic process (GO: 0008152) and response to stimulus (GO: 0050896) were the most enriched GO groups. As expected, the most enriched GO groups of W1-S1 (stage 1), W2-S2 (stage 16
ACCEPTED MANUSCRIPT 2) and W3-S3 (stage 3) were similar and all contained groups related to cellular activity. The GO groups “nutrient reservoir activity” (GO: 0045735), “extracellular matrix component” (GO: 0044420), “extracellular matrix” (GO: 0031012), “metallochaperone activity” (GO: 0016530) and “nucleoid” (GO: 0009295) were identified with a few DEGs at three stages. It might suggest that the genes relevant to these groups may not be the main factors which lead to the grain shape
PT
differences between NIL-GW8 and NIL-gw8Amol.
RI
To explore the function of GW8 and other DEGs related panicle development. We found GW8 and these genes were simultaneously assigned to single-organism process
SC
(GO: 0044699), developmental process (GO: 0032502), biological regulation (GO: 0065007), multicellular organismal process (GO: 0032501), reproductive process (GO:
NU
0022414), reproduction (GO: 0000003), DNA binding (GO: 0005488), organelle part (GO: 0044422), cell (GO: 0005623), cell part (GO: 0044464), organelle (GO:
MA
0043226).
D
3.6. Metabolic pathway by KEGG analysis of DEGs
PT E
In order to further understand the biological functions of the DEGs and summarize the biological pathways associated panicle development, the DEGs in each comparison were mapped to KEGG database (http://www.genome.ad.jp/kegg/). The
CE
DEGs significantly enriched pathways involved in metabolic or signal transduction were identified. Especially, the analysis of KEGG pathway were performed for the
AC
DEGs in different comparisons between NIL-GW8 and NIL-gw8Amol, allowing a rapid visualization and an overview over metabolism alterations about the biological pathways correlated with the formation of rice panicle. At stage 1, a total of 70 significant DEGs with pathway annotation were assigned to 57 KEGG pathways. Among the mapped pathways, nine pathways were significantly enriched (Pvalue < 0.05) in the comparison of W1-S1. The most enriched metabolic pathways were observed in taurine and hypotaurine metabolism (three genes, 4.29%), glycolysis / gluconeogenesis (nine genes, 12.86%), phenylpropanoid biosynthesis (six genes, 8.57%). The gene Os04g0404800 17
ACCEPTED MANUSCRIPT (AMP-binding enzyme) was detected in four pathways (Table S4; Fig. 7a). At stage 2 (W2-S2), 176 significant DEGs with pathway annotation were assigned to 75 KEGG pathways and 13 metabolism pathways were significantly enriched (Pvalue < 0.05). The carotenoid biosynthesis, plant hormone signal transduction, phenylpropanoid biosynthesis, phenylalanine metabolism, taurine and hypotaurine
PT
metabolism were identified (Pvalue < 0.001, Qvalue < 0.05) as the top five pathways. We calculated the gene counts in each pathway, there were 14, 27, 14, 10 and 5 DEGs
RI
in the five pathways, respectively. Compared with stage 1, the analysis showed taurine and hypotaurine metabolism and phenylpropanoid biosynthesis were
SC
significantly identified at stage 2 (Table S4; Fig. 7b).
At stage 3, 295 significant DEGs with pathway annotation were assigned to 100
NU
KEGG pathways. As a result, there were also a large amount of annotated DEGs were significantly assigned to thirty pathways (Pvalue < 0.05). Compared with W1-S1 and
MA
W2-S2, we found that more DEGs were detected in each pathway in W3-S3. Specific enrichment of genes was observed for pathways involved in phenylpropanoid
D
biosynthesis (25 genes, 25.47%), phenylalanine metabolism (16 genes, 5.42%) and
PT E
plant hormone signal transduction (27 genes, 9.15%) which have been identified as significant enrichment pathways in W1-S1 and W2-S2 as well (Table S4; Fig. 7c). The differences in metabolic and regulatory pathways indicate the different
AC
NIL-gw8Amol.
CE
molecular mechanisms underlying the panicle development between NIL-GW8 and
3.7. Analysis of DEGs identified as putative transcription factors Transcription factors play an important role in the process of plant gene expression and regulation by binding to the specific cis-acting element of target gene promoter region, as like GW8 which belongs to the SBP family. Given that transcription factors might have a major effect on the network of GW8, one of the objectives of our work was to identify transcription factors. We paid more attention to DEGs belonging to transcription factor families and summarized the family, gene description and GO annotation of these transcription factor genes in the comparisons 18
ACCEPTED MANUSCRIPT of W1-S1, W2-S2 and W3-S3 (Table S5; Fig. 7). The transcription factor database (http://drtf.cbi.pku.edu.cn/) was used to investigate transcription factors associated with panicle development. Together, a total of 345 DEGs encode known or putative transcription factors. As shown in Fig. 8, 46 putative transcription factors belonging to 15 transcription factor
PT
families were identified with significantly different expression levels in the W1-S1 comparison. Among of them, thirteen genes were characterized as AP2/ERF
RI
(APETALA2/ethylene responsive factor) gene family, constituting a large proportion of the identified transcription factors at stage 1. Based on FPKM values, twenty-five
SC
putative transcription factor genes were down-regulated and twenty-one were up-regulated.
NU
In the W2-S2 comparison (stage 2), 129 putative transcription factor genes belonging to 31 transcription factor families (e.g., CAMTA, WRKY, NAC and MYB)
MA
were differentially expressed between the panicles of the two rice lines. These genes included relative large number of AP2/ERF, bHLH (basic helix-loop-helix) and
D
WRKY transcription factors, suggesting that they may play important roles in the
PT E
genes network with GW8 correlation.
In the W3-S3 comparison (stage 3), 170 putative transcription factor genes belonging to 32 transcription factor families were significantly differentially
CE
expressed. Compared with stages 1 and 2, there were no DEGs were identified as the members of FAR1, GRF, MIKC, SRS, Trihelix and WOX transcription factor
AC
families. The ERF, bHLH and NAC families have the most genes which contained 22, 21 and 17 putative transcription factor DEGs, respectively. We hypothesized that these genes might be important regulators contribute to the structure difference of NIL-GW8 and NIL-gw8Amol panicles which induced by GW8 and gw8Amol. These
transcription
factor
genes
were
significantly
up-regulated
or
down-regulated, but the trend of expression of them varied at the three stages. For example, the expression of ZF-HD was only detected in W3-S3, and MIKC was only detected in W1-S1 and W2-S2. In short, a lot of transcription factor genes exhibited dynamic expression profiles during the panicle development of NIL-GW8 and 19
ACCEPTED MANUSCRIPT NIL-gw8Amol. Further analysis is needed to explore the involvement of these transcription factor families in the regulation of rice panicle formation.
3.8. Validation of RNA-Seq technology To validate the RNA-Seq results, ten genes were randomly selected for analysis
PT
by qRT-PCR referred to the same RNA extracts as for RNA-seq experiments (Fig. 9). The expression profiles of the ten genes determined by qRT-PCR were consistent
RI
with the RNA-Seq results, confirming a good correspondence between qRT-PCR and
SC
RNA-seq.
4. Discussion
NU
4.1. Biosynthesis of secondary metabolites may contribute to the formation of rice panicle
MA
According to the KEGG pathway analyses of DEGs among the three panicle development stages in NIL-GW8 and NIL-gw8Amol, we found that many DEGs were
D
significantly enriched in the pathways related to biosynthesis of secondary
PT E
metabolites, such as phenylpropanoid biosynthesis (ko00940), phenylalanine metabolism (ko00360), stilbenoid, diarylheptanoid and gingerol biosynthesis (ko00945) and flavonoid biosynthesis (ko00941). Among of them, the pathway of
CE
flavonoid biosynthesis is a branch of phenylpropanoid pathway which is associated with phenylalanine metabolism (Winkel-Shirley, 2001; Buer et al., 2010). We
AC
speculate that these secondary metabolites may contribute to the formation of rice panicle.
We all know that plant secondary metabolites constitute an unusually large group of structurally diversified compounds that can be produced in planta from various primary metabolites or their biosynthetic intermediates, constitutively or in response to different environmental stimuli (Piasecka et al., 2015), and were not generally regarded as indispensable components for the normal life activities. Decades of studies have suggested that the secondary metabolites were necessary in maintaining normal metabolism and completing life cycle of the plants and play an important role 20
ACCEPTED MANUSCRIPT in the life activities (Kutchan, 2001). For instance, gibberellin (GA3) and indole acetic acid (IAA) are involved in regulating many activities in plant cells; terpenes (chlorophyll and carotenoid) are involved in the processes of photosynthesis; flavonoids participates in auxin transport by regulating the activity of phosphatase and protein kinase in shoot and root apical meristem; carotenoids content influenced the
PT
drought tolerance in dendrobium moniliforme (Peer and Murphy, 2007; Pattanayak and Tripathy, 2016; Wu et al., 2016).
RI
The phenotypic difference of NIL-GW8 and NIL-gw8Amol is mainly because of the different gene function of GW8 and gw8Amol which belong to SPL (SBP-like)
SC
transcription factor family. In recent years, research found that many transcription factors participated in the regulation of plant secondary metabolites biosynthesis and
NU
accumulation (Appelhagen et al., 2011). Two Petunias transcription factor genes an2 and an11 have been shown to be effective in regulating gene expression in the
MA
flavonoids biosynthetic pathway (Quattrocchio et al., 1998). Overexpression of AtSPL9 and AtSPL10 remarkably negatively regulated the synthesis of anthocyanin,
D
and it show a negative relationship between the activity of miR156 and the contents of
PT E
flavonoids in Arabidopsis (Gou et al., 2011). Based on these findings, we focused on the expression of the genes related to secondary metabolites in our research. Just like what we thought, we found some
CE
genes were significantly enriched in the biosynthetic pathway of secondary metabolites. For example, Os04g0310900 (4-coumarate-CoA ligase) was enriched in
AC
phenylalanine metabolism pathway (Table 2). Despite the functions of secondary metabolites in the panicle are not entirely clear, these studies confirm that they do play important roles in the panicle formation. In summary, we speculate that GW8 regulate the growth and development of rice panicle by influencing the biosynthesis of secondary metabolites.
4.2. Genes related to plant hormone signal transduction are important for the formation of rice panicle In higher plants, the generation and establishment of meristem identity of 21
ACCEPTED MANUSCRIPT axillary meristem have been regard as the key factors of plant morphological and inflorescence architecture (Oikawa and Kyozuka, 2009). And he growth and development of axillary meristem are regulated by plant hormones (McSteen, 2009). In this study, a lot of DEGs in panicles between NIL-GW8 and NIL-gw8Amol were categorized into the pathway of plant hormone signal transduction, which was
PT
represented by many hormones-related genes, especially at stages 2 and 3, indicating that these genes are likely to be involved in the regulation of panicle development in
RI
rice by participating directly or indirectly in the plant hormone signal transduction. As a signaling molecule, auxin is one of most important hormones involved in a
SC
variety of development processes in plants. Studies in recent years have shown that auxin is involved in the modulating development of rice panicles and spikelets
NU
through direct induction of other growth regulating genes (Gao et al., 2016). And other research shows that auxin is necessary for the rice axillary meristem initiation
MA
and the second branch differentiation (Benjamins and Scheres, 2008; Mockaitis and Estelle, 2008; McSteen, 2009). Also, auxin is also involved in regulating the
D
inflorescence meristem development (Zhao, 2008). During the rice reproductive
PT E
process, the content of auxin in lateral meristem determined the number of primary branches and secondary branches (Komatsu et al., 2003). ASP1, a regulatory factor of auxin, controls the development of panicles in rice, and smaller number of branches
CE
and abnormal panicles were found in the mutant asp1 (Yoshida et al., 2012). In this research, many genes related to the auxin (IAA) synthesis (OsJRA2 and LC1), auxin
AC
transport (OsAUX1) and auxin response (OsSAURs) were significantly differentially expressed between the NIL-GW8 and NIL-gw8Amol panicles. It indicated that these genes might be involved in panicle developmental events. During the study many cytokinin-signalling and cytokinin-responsive genes were among the transcripts that were differentially expressed during panicle development of the two rice lines at different stages, such as a series of A-type and B-type response regulators (OsRR1, OsRR6, OsRR11, ORR1 and ORR4). Some theoretical and experimental evidence showed that cytokinin is tightly correlated with the panicle length and the number of branches and spikelets per panicle in rice (Ashikari et al., 22
ACCEPTED MANUSCRIPT 2005; Kurakawa et al., 2007; Li et al., 2013). The QTL GNP1, encoding GA20ox1, increases grain number and yield by increasing cytokinin activity in rice panicle meristems (Wu et al., 2016). Overexpression of OsRR6, a negative regulator of cytokinin signal transduction pathway, has resulted in insufficient development of panicles in rice (Hirose et al., 2007). In addition, a few other studies also
PT
demonstrated that A-type and B-type response regulator families have a vital role in the cytokinin signaling pathway (Zhao et al., 2015; Cho et al., 2016). Accordingly, we
RI
hypothesized that the expression variations of OsRRs and ORRs between NIL-GW8 and NIL-gw8Amol panicles may affect the signal transduction of cytokinin, resulting in
SC
the phenotypic differences of panicles.
In Arabidopsis, plant hormone jasmonic acid (JA) participates in the signal
NU
transduction pathway through activating the expression of JA-responsive genes by degrading the JAZ (jasmonate ZIM-domain) repressors (Chini et al., 2007). A
MA
previous study on spikelet in rice showed that JA was associated with the determination of rice spikelet morphogenesis, and the OsJAZ1 mutant exhibited
D
altered spikelet morphology with changed floral organ identity and number (Cai et al.,
PT E
2014). Decreased expression of OsCOI1 (Coronatine-insensitive 1), the receptor of JA, in rice RNAi lines resulted in long panicles compared with the wild type (Yang et al., 2012). In this study, we identified nine JAZ genes (OsJAZ3, OsJAZ4, OsJAZ6,
CE
OsJAZ7, OsJAZ9, OsJAZ10, OsJAZ11, OsJAZ12 and OsJAZ13) and a JAZ-interacting transcription factor (OsMYC2) in plant hormone signal transduction pathway.
AC
In addition to the above genes, according to the analysis of KEGG pathways, several genes mainly related to other hormone signal transduction pathways of plants were detected between the two rice lines, such as ABA-signalling genes (OsSAPK1, OsSAPK2 and OsSAPK9), ABA responsive element binding factor (OsABF1), ethylene-insensitive gene (OsEIL1) and ethylene receptor (ETR2). So, the expressions of these genes are likely to be affected by GW8 and play integral roles in the process of growth and development of panicles. In summary, our results provide important information for future functional studies of plant hormone signal transduction-related genes and a research foundation 23
ACCEPTED MANUSCRIPT for the mechanisms involved in rice panicle formation.
4.3. Difference in the regulation pattern of the formation of rice panicle between GW8 and GW7 Our previous research found that GW7 is a major QTL for both grain length and
PT
width in rice, increase of GW7 expression, longitudinal cell division of the grain was promoted and the transverse cell division was reduced, as a result of the elongation of
RI
grain (Wang et al., 2015). GW8 plays an important role in the development of panicle in rice, GW8 bound directly to the GW7 promoter and repressed its expression (Wang
SC
et al., 2012; Wang et al., 2015). But for this study, we found that there were no statistically significant differences in the expression levels of GW7 between NIL-GW8
NU
and NIL-gw8Amol at the three stages. So, we speculate that GW7 regulates the development of panicle in the latter stage which after the period of floral organs
MA
differentiation (In7). We compared the transcript levels of GW7 in more different panicle differentiation stages of rice. The result was as we’d expect, qRT-PCR
D
analysis identified differences in GW7 transcript abundance between NIL-GW8 and
PT E
NIL-gw8Amol, and the gene being transcribed more strongly in NIL-gw8Amol than in NIL-GW8 after the stage In7 (Fig. S3). Moreover, the expressions of GW7-related genes were analyzed among the
CE
different developmental stages in NIL-GW8 and NIL-gw8Amol. In the previous study, we have identified 18 candidate genes which interact with GW7 by yeast two-hybrid
AC
assay with GW7 as bait (Wang et al., 2015). For these 18 GW7-related genes, only three genes (Os01g0153800, Os02g0818000 and Os03g0135700) showed differences in expression between different stages in NIL-GW8 and NIL-gw8Amol with log2FC >2.0, but no significant difference between NIL-GW8 and NIL-gw8Amol at different stages in this study. Conclusion according to analysis on the expression differences of GW7 and its related genes, GW7 and GW7-related genes might not be the key regulation factors of spikelet formation as lateral branch of the inflorescence meristem. Our previous research has also shown that there are no differences in the number of primary and 24
ACCEPTED MANUSCRIPT secondary branches and spikelets per panicle between NIL-GW7TFA and NIL-gw7HJX74, but there exists a significant difference between NIL-GW8 and NIL-gw8Amol (Wang et al., 2012; Wang et al., 2015). Therefore, it would be assumed that the regulation pattern of GW7 different from GW8 in different stages of rice panicle development, although GW8 bound to the GW7 promoter and regulate its expression and they all
PT
regulate rice grain development.
RI
5. Conclusion
Our transcriptome analysis provided a comprehensive overview of the
SC
transcriptome changes that specifically occur in panicles of two NILs (NIL-GW8 and NIL-gw8Amol), which were characterized by distinct phenotypic differences of panicles.
NU
Based on the analysis of GO and KEGG, we found that gene expression profiles were dynamic during rice panicle development. Various types of stage-specifically
MA
expressed genes were identified from the three panicle developmental stages between NIL-GW8 and NIL-gw8Amol, and the identified DEGs were involved in diverse
D
pathways, including biosynthesis of secondary metabolites and plant hormone signal
PT E
transduction. The identification of DEGs including some MADS-box genes, cytochrome P450 (CYP) genes, pentatricopeptide repeat (PPR) protein genes and various transcription factors will enhance our understanding of the mechanism of
AC
CE
panicle formation affected by GW8.
25
ACCEPTED MANUSCRIPT Conflict of interests All authors declare that they have no conflict of interest.
Authors’ contributions S.K. performed the bioinformatics analysis and wrote the manuscript. X.J.L.
PT
conducted the sample collection and library construction. X.L. performed the RNA extraction. W.Y. performed the qRT-PCR. H.Z. participated in study coordination and
RI
project, and provided essential suggestions for this work. G.L. supervised and complemented the writing. G.Z. and S.W. conceived the study, participated in its
SC
design and coordination, and corrected the manuscript. All authors read and approved
NU
the final manuscript.
Acknowledgments
MA
This work was supported by the National Natural Science Foundation of China (91635302) and Foundation for Young Talents in Higher Education of Guangdong,
AC
CE
PT E
D
China (2014KQNCX037).
26
ACCEPTED MANUSCRIPT Figure legends Fig. 1. The phenotypic differences of the grain formed by NIL-GW8 and NIL-gw8Amol. Appearance of the grains (a). Grain length (b). Grain width (c). Bar = 5mm. Each column represents the mean ± s.e.m. (n = 120). Asterisks represent significant
PT
difference determined by Student's t-test at p-value < 0.01 (**).
Fig. 2. Length distribution of genes characterized from RNA-seq libraries in rice
SC
RI
developing panicles.
Fig. 3. The distribution of gene expression at three different panicle development
NU
stage in NIL-GW8 and NIL-gw8Amol. The twelve rice chromosomes are displayed by the first outer ring. Six concentric rings from outside to inside represent the
MA
corresponding six libraries (W1, W2, W3, S1, S2, and S3) in this figure, respectively. Colors represent (from green to red) the expression levels of all detected transcripts in
PT E
fragments mapped) value.
D
every library based on the FPKM (Fragments Per Kilobase of transcript per Million
Fig. 4. Venn diagram showing the genes expressed in each of six samples of of rice
CE
panicles development in NIL-GW8 and NIL-gw8Amol. (1) W1: stage 1 (In5) in NIL-GW8; (2) W2: stage 2 (In6) in NIL-GW8; (3) W3: stage 3 (In7) in NIL-GW8; (4) S1: stage 1 (In5) in NIL-gw8Amol; (5) S2: stage 2 (In6) in NIL-gw8Amol; (6) S3: stage 3
AC
(In7) in NIL-gw8Amol.
Fig. 5. Analysis of DEGs in each comparison among the different panicles development stages in NIL-GW8 and NIL-gw8Amol. The number of up-regulated and down-regulated genes in different comparisons is summarized.
Fig. 6. Functional analysis of differentially expressed genes (DEGs) based on RNA-Seq data in different comparisons. W1-S1 (a), W2-S2 (b), W3-S3 (c). 27
ACCEPTED MANUSCRIPT
Fig. 7. KEGG enrichment of DEGs at different panicle stages between NIL-GW8 and NIL-gw8Amol. The scatterplot of enriched GO terms between W1 and S1 (a); the scatterplot of enriched GO terms between W2 and S2 (b); the scatterplot of enriched
PT
GO terms between W3 and S3 (c).
Fig. 8. The number of DEGs belonging to different transcription factor families
RI
detected in the W1-S1, W2-S2, and W3-S3 comparisons. In total, we identified 345
SC
putative TF DEGs, which are classified into 38 families.
Fig. 9. qRT-PCR validation of the relative expression data of ten genes obtained in
NU
RNA-seq analysis. y axes represent the relative expression levels at different stages in NIL-GW8 and NIL-gw8Amol which are represented relative to W1 (assigned as a value
MA
of 1). OsActin1 was used as an internal control gene for transcript normalization.
PT E
Table legends
D
Error bars show standard deviations for three independent biological assays.
CE
Table 1 Summary of reads mapping to rice genome
AC
Table 2 Statistics of the main up-expressed and down-expressed genes
Table 3 Statistics of significant differently expressed MADS-box genes
Table 4 Statistics of significant differently expressed Cytochrome P450 genes
Table 5 Statistics of significant differently expressed PPR genes
28
PT E
D
MA
NU
SC
RI
PT
ACCEPTED MANUSCRIPT
CE
Fig. 1. The phenotypic differences of the grain formed by NIL-GW8 and NIL-gw8Amol. Appearance of the grains (a). Grain length (b). Grain width (c). Bar = 5mm. Each
AC
column represents the mean ± s.e.m. (n = 120). Asterisks represent significant difference determined by Student's t-test at p-value < 0.01 (**).
29
MA
NU
SC
RI
PT
ACCEPTED MANUSCRIPT
Fig. 2. Length distribution of genes characterized from RNA-seq libraries in rice
AC
CE
PT E
D
developing panicles.
30
PT E
D
MA
NU
SC
RI
PT
ACCEPTED MANUSCRIPT
Fig. 3. The distribution of gene expression at three different panicle development
CE
stage in NIL-GW8 and NIL-gw8Amol. The twelve rice chromosomes are displayed by the first outer ring. Six concentric rings from outside to inside represent the
AC
corresponding six libraries (W1, W2, W3, S1, S2, and S3) in this figure, respectively. Colors represent (from green to red) the expression levels of all detected transcripts in every library based on the FPKM (Fragments Per Kilobase of transcript per Million fragments mapped) value.
31
NU
SC
RI
PT
ACCEPTED MANUSCRIPT
MA
Fig. 4. Venn diagram showing the genes expressed in each of six samples of of rice panicles development in NIL-GW8 and NIL-gw8Amol. (1) W1: stage 1 (In5) in NIL-GW8; (2) W2: stage 2 (In6) in NIL-GW8; (3) W3: stage 3 (In7) in NIL-GW8; (4)
AC
CE
PT E
(In7) in NIL-gw8Amol.
D
S1: stage 1 (In5) in NIL-gw8Amol; (5) S2: stage 2 (In6) in NIL-gw8Amol; (6) S3: stage 3
32
PT E
D
MA
NU
SC
RI
PT
ACCEPTED MANUSCRIPT
CE
Fig. 5. Analysis of DEGs in each comparison among the different panicles development stages in NIL-GW8 and NIL-gw8Amol. The number of up-regulated and
AC
down-regulated genes in different comparisons is summarized.
33
AC
CE
PT E
D
MA
NU
SC
RI
PT
ACCEPTED MANUSCRIPT
Fig. 6. Functional analysis of differentially expressed genes (DEGs) based on RNA-Seq data in different comparisons. W1-S1 (a), W2-S2 (b), W3-S3 (c).
34
ACCEPTED MANUSCRIPT
T P
I R
C S U
N A
D E
M
T P E
C C
Fig. 7. KEGG enrichment of DEGs at different panicle stages between NIL-GW8 and NIL-gw8Amol. The scatterplot of enriched GO terms
A
between W1 and S1 (a); the scatterplot of enriched GO terms between W2 and S2 (b); the scatterplot of enriched GO terms between W3 and S3 (c).
35
ACCEPTED MANUSCRIPT
T P
I R
C S U
N A
D E
M
T P E
C C
Fig. 8. The number of DEGs belonging to different transcription factor families detected in the W1-S1, W2-S2, and W3-S3 comparisons. In total,
A
we identified 345 putative TF DEGs, which are classified into 38 families.
36
AC
CE
PT E
D
MA
NU
SC
RI
PT
ACCEPTED MANUSCRIPT
Fig. 9. qRT-PCR validation of the relative expression data of ten genes obtained in RNA-seq analysis. y axes represent the relative expression levels at different stages in NIL-GW8 and NIL-gw8Amol which are represented relative to W1 (assigned as a value of 1). OsActin1 was used as an internal control gene for transcript normalization. Error bars show standard deviations for three independent biological assays.
37
ACCEPTED MANUSCRIPT Table 1 Summary of reads mapping to rice genome. Total Reads
Mapped Reads
Mapped Ratio
Unique Mapped Reads
Unique Mapped Ratio
W1
46,479,240
37,345,492
80.35%
36,297,170
78.09%
W2
53,148,798
42,452,813
79.88%
41,216,953
77.55%
W3
54,120,416
43,417,595
80.22%
42,234,568
78.04%
S1
40,492,882
33,347,884
82.35%
32,398,176
80.01%
S2
40,822,488
33,533,503
82.14%
32,345,222
79.23%
S3
46,103,974
37,048,147
80.36%
35,813,220
77.68%
AC
CE
PT E
D
MA
NU
SC
RI
PT
Sample
38
ACCEPTED MANUSCRIPT Table 2 Statistics of the main up-expressed and down-expressed genes. Annotation
Regulated
NewGene_446
Hypothetical protein
NewGene_3734
Unknown protein
Os08g0527600
log2 FC W2 - S2
W3 - S3
up
9.17
9.10
6.04
up
3.78
5.24
7.53
Unnamed protein product
up
6.37
6.37
7.38
Os12g0210601
Hypothetical protein
up
7.23
7.36
7.12
Os08g0528100
Ubiquitin carboxyl-terminal
up
5.78
5.81
7.10
6.91
6.35
5.91
5.14
8.56
-
7.92
6.57
-
6.61
6.21
5.69
-5.43
-4.19
6.33
up
4.56
3.83
5.72
up/down
-6.78
-
5.62
up/down
6.13
-3.93
4.83
up/down
5.85
-6.17
-5.37
up
5.74
5.62
-
PT
W1 - S1
NU
Gene ID
Hypothetical protein
up
Os05g0469600
Pyruvate decarboxylase
up
Os03g0725200
Formin-like protein
up
Os10g0366400
Unknown protein
up
Os02g0127600
putative pentatricopeptide (PPR) repeat-containing protein hypothetical protein
Os02g0260500
Uncharacterized protein
Os01g0159000
Proline-rich receptor
Os01g0589900
putative pentatricopeptide (PPR)
up/down
MA
NewGene_3514
repeat-containing protein
SC
NewGene_3689
RI
hydrolase
Pectinesterase
Os03g0172100
Protein CHUP1, chloroplastic
up
5.59
5.44
Os05g0469701
Pyruvate decarboxylase
up
5.25
7.07
Os04g0310900
4-coumarate-CoA ligase
up
5.22
-
3.79
Os05g0162800
Transcription factor
up
-
6.45
4.87
PT E
D
Os02g0288100
Os06g0100550
Unknown protein
up
4.35
4.56
4.83
Hypothetical protein
down
-8.08
-6.85
-6.32
Hypothetical protein
down
-7.14
-5.19
-5.69
Unknown protein
down
-6.52
-6.21
-
Hypothetical protein
down
-6.30
-7.65
-5.95
Os09g0513900
Unnamed protein product
down
-5.95
-5.83
-4.38
Os12g0234000
Unknown protein
down
-5.91
-5.47
-2.05
Os01g0290800
ABC transporter B family
up/down
4.85
-5.72
1.55
down
-
4.91
-
Os10g0370100 Os07g0646200
AC
NewGene_3441
CE
NewGene_448
member NewGene_3073
Hypothetical protein
Os06g0361500
Unknown protein
down
-
-4.80
-4.30
Os03g0226200
Non-symbiotic hemoglobin
up/down
-
3.06
-7.50
Os12g0446500
metallopeptidase family M24
up/down
-
6.01
-7.09
down
-5.89
-
-5.56
down
-1.89
-2.30
-5.47
containing protein Os03g0229500
MATE efflux family protein
Os02g0586900
Glycine-rich cell wall structural protein (Precursor)
39
ACCEPTED MANUSCRIPT Table 3 Statistics of significant differently expressed MADS-box genes.
7.386 7.604 7.420 7.436 4.537 4.247 4.915 3.341 4.215 3.575 2.621 2.269 -3.252 -4.177 -2.024 -2.977 -3.903 -3.734 -7.073
PT
6.924 6.669 6.663 6.662 3.932 3.580 3.557 3.396 2.837 2.779 2.742 2.112 -2.467 -2.959 -3.153 -3.612 -4.549 -4.640 -5.823
AC
CE
PT E
D
S2-S3
RI
OsMADS13 OSMADS3 OsMADS7 OsMADS8 OsMADS17 OsMADS21 OsMADS1 OsMADS29 OsMADS58 OsMADS6 OsMADS4 OsMADS56 OsMADS26 OsMADS22 OsMDP1 OsMADS34 OsMADS50 OsMADS20 OsMADS32 OsMADS55
W2-W3
SC
Os12g0207000 Os01g0201700 Os08g0531700 Os09g0507200 Os04g0580700 Os01g0886200 Os03g0215400 Os02g0170300 Os05g0203800 Os02g0682200 Os05g0423400 Os10g0536100 Os08g0112700 Os02g0761000 Os03g0186600 Os03g0753100 Os03g0122600 Os12g0501700 Os01g0726400 Os06g0217300
log2 FC
NU
Gene Name
MA
Gene ID
40
ACCEPTED MANUSCRIPT Table 4 Statistics of significant differently expressed Cytochrome P450 genes. log2 FC
Gene Name
Os10g0525000 Os01g0804400 Os03g0332100 Os02g0186900 Os04g0570500 Os05g0482400 Os04g0573900 Os09g0441400 Os04g0407900 Os07g0635300 Os01g0227800 Os07g0603700 Os11g0289700 Os10g0486100 Os10g0403000
CYP714B2
S2-S3
6.510 6.179 5.856 5.782 5.777 4.952 4.704 4.548 4.035 4.022 3.993 3.935 3.003 -7.630
5.025 5.241 4.163 6.060 6.485 6.312 5.416 5.382 4.098 3.221 4.943 5.655 4.135 -6.811
SC
RI
EUI1
W2-W3
PT
Gene ID
NU
CYP78A13
AC
CE
PT E
D
MA
PLA1
41
ACCEPTED MANUSCRIPT Table 5 Statistics of significant differently expressed PPR genes. log2 FC
Gene ID
W1-W2
Os01g0589900 Os08g0402600 Os07g0244400 Os03g0363700 Os12g0552300 Os02g0759500 Os02g0750400 Os08g0340900 Os10g0116000 Os11g0433101 Os01g0815900 Os07g0635800 Os07g0621100 Os02g0127600 Os03g0283500 Os02g0170000 Os07g0260000
W2-W3
S1-S2
6.407
-6.071 4.894
-4.247 -2.258
PT
4.870 2.101 -4.793
S2-S3
RI
4.511 -4.545
NU
SC
4.740 4.144
MA
3.935 -3.599 -3.172
4.098 3.745 -3.604 -3.810 6.446 2.984
AC
CE
PT E
D
-3.148
42
ACCEPTED MANUSCRIPT References Appelhagen, I., Lu, G., Huep, G., Schmelzer, E., Weisshaar, B., Sagasser, M., 2011. TRANSPARENT TESTA1 interacts with R2R3-MYB factors and affects early and late steps of flavonoid biosynthesis in the endothelium of Arabidopsis thaliana seeds. Plant J. 67 (3), 406-419. https://doi.org/10.1111/j.1365-313X.2011.04603.x. Ashikari, M., Matsuoka, M., 2006. Identification, isolation and pyramiding of quantitative trait loci for rice breeding. Trends Plant Sci. 11 (7), 344-350. https://doi.org/10.1016/j.tplants.2006.05.008. Ashikari, M., Sakakibara, H., Lin, S.Y., Yamamoto, T., Takashi, T., Nishimura, A., Angeles, E.R., Qian, Q., Kitano, H., Matsuoka, M., 2005. Cytokinin oxidase regulates rice grain production.
PT
Science 309 (5735), 741-745. https://doi.org/10.1126/science.1113373.
Benjamins, R., Scheres, B., 2008. Auxin: The looping star in plant development. Annu. Rev. Plant Biol.
RI
59 (59), 443-465. https://doi.org/10.1146/annurev.arplant.58.032806.103805.
Binder, S., Holzle, A., Jonietz, C., 2011. RNA Processing and RNA Stability in Plant Mitochondria.
SC
Adv. Plant Biochem. Mol. Biol. 1, 107-130. https://doi.org/10.1007/978-0-387-89781-3_5. Bolwell, G.P., Bozak, K., Zimmerlin, A., 1994. Plant cytochrome P450. Phytochemistry 37, 1491-1506. https://doi.org/10.1016/S0031-9422(00)89567-9.
NU
Buer, C.S., Imin, N., Djordjevic, M.A., 2010. Flavonoids: New Roles for Old Molecules. J. Integr. Plant Biol. 52, 98-111. https://doi.org/10.1111/j.1744-7909.2010.00905.x. Cai, Q., Yuan, Z., Chen, M., Yin, C., Luo, Z., Zhao, X., Liang, W., Hu, J., Zhang, D., 2014. Jasmonic regulates
spikelet
development
in
rice.
Nat.
Commun.
5
(3),
319-333.
MA
acid
https://doi.org/10.1038/ncomms4476.
Chen, R., Shen, L., Wang, D., Wang, F., Zeng, H., Chen, Z., Peng, Y., Lin, Y., Tang, X., Deng, M., Yao, N., Luo, J., Xu, Z., Bai, S., 2015. A Gene Expression Profiling of Early Rice Stamen
D
Development that Reveals Inhibition of Photosynthetic Genes by OsMADS58. Mol. Plant. 8 (7), 1069-1089. https://doi.org/10.1016/j.molp.2015.02.004.
PT E
Chini, A., Fonseca, S., Fernandez, G., Adie, B., Chico, J.M., Lorenzo, O., Garcia-Casado, G., Lopez-Vidriero, I., Lozano, F.M., Ponce, M.R., Micol, J.L., Solano, R., 2007. The JAZ family of repressors is the missing link in jasmonate signalling. Nature 448 (7154), 666-671. https://doi.org/10.1038/nature06006.
CE
Cho, L., Yoon, J., Pasriga, R., An, G., 2016. Homodimerization of Ehd1 Is Required to Induce Flowering in Rice. Plant Physiol. 170 (4), 2159-2171. https://doi.org/10.1104/pp.15.01723. Dahan, J., Mireau, H., 2013. The Rf and Rf-like PPR in higher plants, a fast-evolving subclass of PPR
AC
genes. RNA Biol. 10 (9), 1469-1476. https://doi.org/10.4161/rna.25568. Deshpande, G.M., Ramakrishna, K., Chongloi, G.L., Vijayraghavan, U., 2015. Functions for rice RFL in vegetative axillary meristem specification and outgrowth. J. Exp. Bot. 66 (9), 2773-2784. https://doi.org/10.1093/jxb/erv092. Devos, K.M., 2005. Updating the 'Crop circle'. Curr. Opin. Plant Biol. 8 (2), 155-162. https://doi.org/10.1016/j.pbi.2005.01.005. Ding, Y.H., Liu, N.Y., Tang, Z.S., Liu, J., Yang, W.C., 2006. Arabidopsis GLUTAMINE-RICH PROTEIN23 is essential for early embryogenesis and encodes a novel nuclear PPR motif protein that interacts with RNA polymerase II subunit III. Plant Cell 18 (4), 815-830. https://doi.org/10.1105/tpc.105.039495. Duan, K., Li, L., Hu, P., Xu, S., Xu, Z., Xue, H., 2006. A brassinolide-suppressed rice MADS-box transcription factor, OsMDP1, has a negative regulatory role in BR signaling. Plant J. 47 (4), 43
ACCEPTED MANUSCRIPT 519-531. https://doi.org/10.1111/j.1365-313X.2006.02804.x. Duan, P., Ni, S., Wang, J., Zhang, B., Xu, R., Wang, Y., Chen, H., Zhu, X., Li, Y., 2016. Regulation of OsGRF4 by OsmiR396 controls grain size and yield in rice. Nat. Plants 2 (1), 15203. http://dx.doi.org/10.1038/nplants.2015.203. Endo-Higashi, N., Izawa, T., 2011. Time Genes Heading date 1 and Early heading date 1 Together Control
Panicle
Development
in
Rice.
Plant
Cell
Physiol.
52
(6),
1083-1094.
https://doi.org/10.1093/pcp/pcr059. Ewing, B., Hillier, L., Wendl, M.C., Green, P., 1998. Base-calling of automated sequencer traces using phred. I. Accuracy assessment. Genome Res. 8 (3), 175-185. https://doi.org/10.1101/gr.8.3.175.
PT
Fan, C.H., Xing, Y.Z., Mao, H.L., Lu, T.T., Han, B., Xu, C.G., Li, X.H., Zhang, Q.F., 2006. GS3, a major QTL for grain length and weight and minor QTL for grain width and thickness in rice,
RI
encodes a putative transmembrane protein. Theor. Appl. Genet. 112 (6), 1164-1171. https://doi.org/10.1007/s00122-006-0218-1.
SC
Gao, F., Wang, K., Liu, Y., Chen, Y., Chen, P., Shi, Z., Luo, J., Jiang, D., Fan, F., Zhu, Y., Li, S., 2016. Blocking miR396 increases rice yield by shaping inflorescence architecture. 15196. https://doi.org/10.1038/nplants.2015.196.
Nat. Plants 2 (1),
NU
Gou, J., Felippes, F.F., Liu, C., Weigel, D., Wang, J., 2011. Negative Regulation of Anthocyanin Biosynthesis in Arabidopsis by a miR156-Targeted SPL Transcription Factor. Plant Cell 23 (4), 1512-1522. https://doi.org/10.1105/tpc.111.084525.
MA
Hirose, N., Makita, N., Kojima, M., Kamada-Nobusada, T., Sakakibara, H., 2007. Overexpression of a type-A response regulator alters rice morphology and cytokinin metabolism. Plant Cell Physiol. 48 (3), 523-539. https://doi.org/10.1093/pcp/pcm022.
Hu, Y., Liang, W., Yin, C., Yang, X., Ping, B., Li, A., Jia, R., Chen, M., Luo, Z., Cai, Q., Zhao, X.,
D
Zhang, D., Yuan, Z., 2015. Interactions of OsMADS1 with Floral Homeotic Genes in Rice Flower Development. Mol. Plant 8 (9), 1366-1384. https://doi.org/10.1016/j.molp.2015.04.009.
PT E
Huang, X., Qian, Q., Liu, Z., Sun, H., He, S., Luo, D., Xia, G., Chu, C., Li, J., Fu, X., 2009. Natural variation at the DEP1 locus enhances grain yield in rice. Nat. Genet. 41 (4), 494-497. https://doi.org/10.1038/ng.352.
Ikeda, K., Sunohara, H. and Nagato, Y., 2004. Developmental course of inflorescence and spikelet in
CE
rice. Breed. Sci. 54 (2), 147-156. https://doi.org/10.1270/jsbbs.54.147. Ikeda-Kawakatsu, K., Maekawa, M., Izawa, T., Itoh, J., Nagato, Y., 2012. ABERRANT PANICLE ORGANIZATION 2/RFL, the rice ortholog of Arabidopsis LEAFY, suppresses the transition from
AC
inflorescence meristem to floral meristem through interaction with APO1. Plant J. 69 (1), 168-180. https://doi.org/10.1111/j.1365-313X.2011.04781.x. Ishii, T., Numaguchi, K., Miura, K., Yoshida, K., Pham, T.T., Than, M.H., Yamasaki, M., Komeda, N., Matsumoto, T., Terauchi, R., Ishikawa, R., Ashikari, M., 2013. OsLG1 regulates a closed panicle trait in domesticated rice. Nat. Genet. 45 (4), 462-465. https://doi.org/10.1038/ng.2567. Ishimaru, K., Hirotsu, N., Madoka, Y., Murakami, N., Hara, N., Onodera, H., Kashiwagi, T., Ujiie, K., Shimizu, B., Onishi, A., Miyagawa, H., Katoh, E., 2013. Loss of function of the IAA-glucose hydrolase gene TGW6 enhances rice grain weight and increases yield. Nat. Genet. 45 (6), 707-711. https://doi.org/10.1038/ng.2612. Itoh, J., Nonomura, K., Ikeda, K., Yamaki, S., Inukai, Y., Yamagishi, H., Kitano H., Nagato Y., 2005. Rice plant development: from zygote to spikelet. Plant Cell Physiol. 46(1), 23-47. https://doi.org/10.1093/pcp/pci501. 44
ACCEPTED MANUSCRIPT Jiang, G., Xiang, Y., Zhao, J., Yin, D., Zhao, X., Zhu, L., Zhai, W., 2014. Regulation of Inflorescence Branch Development in Rice Through a Novel Pathway Involving the Pentatricopeptide Repeat Protein sped1-D. Genetics 197 (4), 1395-1407. https://doi.org/10.1534/genetics.114.163931. Jiao, Y., Wang, Y., Xue, D., Wang, J., Yan, M., Liu, G., Dong, G., Zeng, D., Lu, Z., Zhu, X., Qian, Q., Li, J., 2010. Regulation of OsSPL14 by OsmiR156 defines ideal plant architecture in rice. Nat. Genet. 42 (6), 541-544. https://doi.org/10.1038/ng.591. Kim, D., Pertea, G., Trapnell, C., Pimentel, H., Kelley, R., Salzberg, S.L., 2013. TopHat2: accurate alignment of transcriptomes in the presence of insertions, deletions and gene fusions. Genome Biol. 14, R36. https://doi.org/10.1186/gb-2013-14-4-r36.
PT
Kocabek, T., Repkova, J., Dudova, M., Hoyerova, K., Vrba, L., 2006. Isolation and characterization of a novel semi-lethal Arabidopsis thaliana mutant of gene for pentatricopeptide (PPR)
RI
repeat-containing protein. Genetica 128 (1), 395-407. https://doi.org/10.1007/s10709-006-7518-x. Komatsu, M., Chujo, A., Nagato, Y., Shimamoto, K., Kyozuka, J., 2003. FRIZZY PANICLE is required
SC
to prevent the formation of axillary meristems and to establish floral meristem identity in rice spikelets. Development 130 (16), 3841-3850. https://doi.org/10.1242/dev.00564. Kurakawa, T., Ueda, N., Maekawa, M., Kobayashi, K., Kojima, M., Nagato, Y., Sakakibara, H.,
NU
Kyozuka, J., 2007. Direct control of shoot meristem activity by a cytokinin-activating enzyme. Nature 445 (7128), 652-655. https://doi.org/10.1038/nature05504. Kutchan, T.M., 2001. Ecological arsenal and developmental dispatcher. The paradigm of secondary
MA
metabolism. Plant Physiol. 125 (1), 58-60. https://doi.org/10.1104/pp.125.1.58. Langmead, B., Trapnell, C., Pop, M., Salzberg, S.L., 2009. Ultrafast and memory-efficient alignment of
short
DNA
sequences
to
the
human
genome.
Genome
Biol.
10
(3),
1-10.
https://doi.org/10.1186/gb-2009-10-3-r25.
D
Li, H., Liang, W., Hu, Y., Zhu, L., Yin, C., Xu, J., Dreni, L., Kater, M.M., Zhang, D., 2011. Rice MADS6 Interacts with the Floral Homeotic Genes SUPERWOMAN1, MADS3, MADS58, MADS13,
PT E
and DROOPING LEAF in Specifying Floral Organ Identities and Meristem Fate. Plant Cell 23 (7), 2536-2552. https://doi.org/10.1105/tpc.111.087262. Li, S., Zhao, B., Yuan, D., Duan, M., Qian, Q., Tang, L., Wang, B., Liu, X., Zhang, J., Wang, J., Sun, J., Liu, Z., Feng, Y., Yuan, L., Li, C., 2013. Rice zinc finger protein DST enhances grain production
CE
through controlling Gn1a/OsCKX2 expression. Proc. Natl. Acad. Sci. U. S. A. 110 (8), 3167-3172. https://doi.org/10.1073/pnas.1300359110. Li, Y., Fan, C., Xing, Y., Jiang, Y., Luo, L., Sun, L., Shao, D., Xu, C., Li, X., Xiao, J., He, Y., Zhang,
AC
Q., 2011. Natural variation in GS5 plays an important role in regulating grain size and yield in rice. Nat. Genet. 43 (12), 1266-1269. https://doi.org/10.1038/ng.977. Liu, J.P., Van Eck, J., Cong, B., Tanksley, S.D., 2002. A new class of regulatory genes underlying the cause of pear-shaped tomato fruit. Proc. Natl. Acad. Sci. U. S. A. 99 (20), 13302-13306. https://doi.org/10.1073/pnas.162485999. Luo, A.D., Qian, Q., Yin, H.F., Liu, X.Q., Yin, C.X., Lan, Y., Tang, J.Y., Tang, Z.S., Cao, S.Y., Wang, X.J., Xia, K., Fu, X.D., Luo, D., Chu, C.C., 2006. EUI1, encoding a putative cytochrome P450 monooxygenase, regulates internode elongation by modulating gibberellin responses in rice. Plant Cell Physiol. 47 (2), 181-191. https://doi.org/10.1093/pcp/pci233. Luo, J., Liu, H., Zhou, T., Gu, B., Huang, X., Shangguan, Y., Zhu, J., Li, Y., Zhao, Y., Wang, Y., Zhao, Q., Wang, A., Wang, Z., Sang, T., Wang, Z., Han, B., 2013. An-1 Encodes a Basic Helix-Loop-Helix Protein That Regulates Awn Development, Grain Size, and Grain Number in 45
ACCEPTED MANUSCRIPT Rice. Plant Cell 25 (9), 3360-3376. https://doi.org/10.1105/tpc.113.113589. Magome, H., Nomura, T., Hanada, A., Takeda-Kamiya, N., Ohnishi, T., Shinma, Y., Katsumata, T., Kawaide, H., Kamiya, Y., Yamaguchi, S., 2013. CYP714B1 and CYP714B2 encode gibberellin 13-oxidases that reduce gibberellin activity in rice. Proc. Natl. Acad. Sci. U. S. A. 110 (5), 1947-1952. https://doi.org/10.1073/pnas.1215788110. Mao, H., Sun, S., Yao, J., Wang, C., Yu, S., Xu, C., Li, X., Zhang, Q., 2010. Linking differential domain functions of the GS3 protein to natural variation of grain size in rice. Sci. U. S. A. 107 (45), 19579-19584. https://doi.org/10.1073/pnas.1014419107. McSteen, P., 2009. Hormonal Regulation of Branching in Grasses. Plant Physiol. 149 (1), 46-55.
PT
https://doi.org/10.1104/pp.108.129056.
Meng, Q., Li, X., Zhu, W., Yang, L., Liang, W., Dreni, L., Zhang, D., 2017. Regulatory network and
RI
genetic interactions established by OsMADS34 in rice inflorescence and spikelet morphogenesis. J. Integr. Plant Biol. 59(9), 693-707. https://doi.org/10.1111/jipb.12594.
SC
Mimura, M., Nagato, Y., Itoh, J., 2012. Rice PLASTOCHRON genes regulate leaf maturation downstream of the gibberellin signal transduction pathway. Planta 235 (5), 1081-1089. https://doi.org/10.1007/s00425-012-1639-5.
NU
Miyoshi, K., Ahn, B.O., Kawakatsu, T., Ito, Y., Itoh, J.I., Nagato, Y., Kurata, N., 2004. PLASTOCHRON1, a timekeeper of leaf initiation in rice, encodes cytochrome P450. Proc. Natl. Acad. Sci. U. S. A. 101 (4), 875-880. https://doi.org/10.1073/pnas.2636936100. Annu.
Rev.
Cell
MA
Mockaitis, K., Estelle, M., 2008. Auxin Receptors and Plant Development: A New Signaling Paradigm. Dev.
Biol.
24
(24),
55-80.
https://doi.org/10.1146/annurev.cellbio.23.090506.123214. Moumeni, A., Satoh, K., Venuprasad, R., Serraj, R., Kumar, A., Leung, H., Kikuchi S., 2015.
D
Transcriptional profiling of the leaves of near-isogenic rice lines with contrasting drought tolerance at the reproductive stage in response to water deficit. BMC Genomics, 16(1), 1110.
PT E
https://doi.org/10.1186/s12864-015-2335-1. Nakagawa, M., Shimamoto, K., Kyozuka, J., 2002. Overexpression of RCN1 and RCN2, rice TERMINAL FLOWER 1/CENTRORADIALIS homologs, confers delay of phase transition and altered
panicle
morphology
in
rice.
Plant
J.
29
(6),
743-750.
CE
https://doi.org/10.1046/j.1365-313X.2002.01255.x. Nayar, S., Sharma, R., Tyagi, A.K., Kapoor, S., 2013. Functional delineation of rice MADS29 reveals its role in embryo and endosperm development by affecting hormone homeostasis. J. Exp. Bot.
AC
64 (14), 4239-4253. https://doi.org/10.1093/jxb/ert231. Oikawa, T., Kyozuka, J., 2009. Two-Step Regulation of LAX PANICLE1 Protein Accumulation in Axillary
Meristem
Formation
in
Rice.
Plant
Cell
21
(4),
1095-1108.
https://doi.org/10.1105/tpc.108.065425. Pattanayak, G.K., Tripathy, B.C., 2016. Modulation of biosynthesis of photosynthetic pigments and light-harvesting complex in wild-type and gun5 mutant of Arabidopsis thaliana during impaired chloroplast
development.
Protoplasma
253
(3),
747-752.
https://doi.org/10.1007/s00709-016-0958-y. Peer, W.A. and Murphy, A.S., 2007. Flavonoids and auxin transport: modulators or regulators? Trends Plant Sci. 12 (12), 556-563. https://doi.org/10.1016/j.tplants.2007.10.003. Piasecka, A., Jedrzejczak-Rey, N., Bednarek, P., 2015. Secondary metabolites in plant innate immunity: conserved function of divergent chemicals.
New
Phytol.
206
(3),
948-964. 46
ACCEPTED MANUSCRIPT https://doi.org/10.1111/nph.13325. Qi, P., Lin, Y., Song, X., Shen, J., Huang, W., Shan, J., Zhu, M., Jiang, L., Gao, J., Lin, H., 2012. The novel quantitative trait locus GL3.1 controls rice grain size and yield by regulating Cyclin-T1;3. Cell Res. 22 (12), 1666-1680. https://doi.org/10.1038/cr.2012.151. Quattrocchio, F., Wing, J.F., Van Der Woude, K., Mol, J.N.M., Koes, R., 1998. Analysis of bHLH and MYB domain proteins: Species-specific regulatory differences are caused by divergent evolution of
target
anthocyanin
genes.
Plant
J.
13
(4),
475-488.
https://doi.org/10.1046/j.1365-313X.1998.00046.x. Shikanai, T., Fujii, S., 2013. Function of PPR proteins in plastid gene expression. RNA Biol. 10 (9),
PT
1446-1456. https://doi.org/10.4161/rna.25207.
Shomura, A., Izawa, T., Ebana, K., Ebitani, T., Kanegae, H., Konishi, S., Yano, M., 2008. Deletion in a
RI
gene associated with grain size increased yields during rice domestication. Nat. Genet. 40 (8), 1023-1028. https://doi.org/10.1038/ng.169.
SC
Si, L., Chen, J., Huang, X., Gong, H., Luo, J., Hou, Q., Zhou, T., Lu, T., Zhu, J., Shangguan, Y., Chen, E., Gong, C., Zhao, Q., Jing, Y., Zhao, Y., Li, Y., Cui, L., Fan, D., Lu, Y., Weng, Q., Wang, Y., Zhan, Q., Liu, K., Wei, X., An, K., An, G., Han, B., 2016. OsSPL13 controls grain size in
NU
cultivated rice. Nat. Genet. 48 (4), 447-456. https://doi.org/10.1038/ng.3518. Song, X., Huang, W., Shi, M., Zhu, M., Lin, H., 2007. A QTL for rice grain width and weight encodes a previously unknown RING-type E3 ubiquitin ligase. Nat. Genet. 39 (5), 623-630.
MA
https://doi.org/10.1038/ng2014.
Sundaresan, V., 2005. Control of seed size in plants. Proc. Natl. Acad. Sci. U. S. A. 102 (50), 17887-17888. https://doi.org/10.1073/pnas.0509021102. Tabuchi, H., Zhang, Y., Hattori, S., Omae, M., Shimizu-Sato, S., Oikawa, T., Qian, Q., Nishimura, M.,
D
Kitano, H., Xie, H., Fang, X., Yoshida, H., Kyozuka, J., Chen, F., Sato, Y., 2011. LAX PANICLE2 of Rice Encodes a Novel Nuclear Protein and Regulates the Formation of Axillary Meristems.
PT E
Plant Cell 23 (9), 3276-3287. https://doi.org/10.1105/tpc.111.088765. Tan, L., Li, X., Liu, F., Sun, X., Li, C., Zhu, Z., Fu, Y., Cai, H., Wang, X., Xie, D., Sun, C., 2008. Control of a key transition from prostrate to erect growth in rice domestication. Nat. Genet. 40 (11), 1360-1364. https://doi.org/10.1038/ng.197.
CE
Tariq, R., Wang, C., Qin, T., Xu, F., Tang, Y., & Gao, Y., Ji Z., Zhao K., 2018. Comparative transcriptome profiling of rice near-isogenic line carrying xa23 under infection of xanthomonas oryzae pv. oryzae. Int. J. Mol. Sci. 19(3), 717. https://doi.org/10.3390/ijms19030717.
AC
Terao, T., Nagata, K., Morino, K., Hirose, T., 2010. A gene controlling the number of primary rachis branches also controls the vascular bundle formation and hence is responsible to increase the harvest
index
and
grain
yield
in
rice.
Theor.
Appl.
Genet.
120
(5),
875-893.
https://doi.org/10.1007/s00122-009-1218-8. Wang, L., Li, P., Brutnell, T.P., 2010. Exploring plant transcriptomes using ultra high-throughput sequencing. Briefings Funct. Genomics 9 (2), 118-128. https://doi.org/10.1093/bfgp/elp057. Wang, L., Sun, S., Jin, J., Fu, D., Yang, X., Weng, X., Xu, C., Li, X., Xiao, J., Zhang, Q., 2015. Coordinated regulation of vegetative and reproductive branching in rice. Proc. Natl. Acad. Sci. U. S. A. 112 (50), 15504-15509. https://doi.org/10.1073/pnas.1521949112. Wang, L., Zeng, X., Zhuang, H., Shen, Y., Chen, H., Wang, Z., Long, J., Ling, Y., He, G., Li, Y., 2017. Ectopic expression of OsMADS1 caused dwarfism and spikelet alteration in rice. Plant Growth Regul. 81 (3), 433-442. https://doi.org/10.1007/s10725-016-0220-9. 47
ACCEPTED MANUSCRIPT Wang, S., Li S., Liu Q., Wu K., Zhang J., Wang S., Wang Y., Chen X., Zhang Y., Gao C., Wang F., Huang H., Fu X., 2015. The OsSPL16-GW7 regulatory module determines grain shape and simultaneously improves rice yield and grain quality. Nat. Genet. 47(8): 949-954. https://doi.org/10.1038/ng.3352. Wang, S., Wu, K., Yuan, Q., Liu, X., Liu, Z., Lin, X., Zeng, R., Zhu, H., Dong, G., Qian, Q., Zhang, G., Fu, X., 2012. Control of grain size, shape and quality by OsSPL16 in rice. Nat. Genet. 44 (8), 950-954. https://doi.org/10.1038/ng.2327. Weng, J., Gu, S., Wan, X., Gao, H., Guo, T., Su, N., Lei, C., Zhang, X., Cheng, Z., Guo, X., Wang, J., Jiang, L., Zhai, H., Wan, J., 2008. Isolation and initial characterization of GW5, a major QTL with
rice
grain
width
and
weight.
Cell
Res.
https://doi.org/10.1038/cr.2008.307.
18
(12),
PT
associated
1199-1209.
RI
Winkel-Shirley, B., 2001. Flavonoid biosynthesis. A colorful model for genetics, biochemistry, cell biology, and biotechnology. Plant Physiol. 126 (2), 485-493. https://doi.org/10.1104/pp.126.2.485.
SC
Wu, Y., Wang, Y., Mi, X., Shan, J., Li, X., Xu, J., Lin, H., 2016. The QTL GNP1 Encodes GA20ox1, Which Increases Grain Number and Yield by Increasing Cytokinin Activity in Rice Panicle Meristems. PLoS Genet. 12 (10): e1006386. https://doi.org/10.1371/journal.pgen.1006386.
NU
Xing, Y., Zhang, Q., 2010. Genetic and molecular bases of rice yield. Plant Biology 61 (61), 421-442. https://doi.org/10.1146/annurev-arplant-042809-112209.
Yang, D., Yao, J., Mei, C., Tong, X., Zeng, L., Li, Q., Xiao, L., Sun, T., Li, J., Deng, X., Lee, C.M.,
MA
Thomashow, M.F., Yang, Y., He, Z., He, S.Y., 2012. Plant hormone jasmonate prioritizes defense over growth by interfering with gibberellin signaling cascade. Proc. Natl. Acad. Sci. U. S. A. 109 (19), 1192-1200. https://doi.org/10.1073/pnas.1201616109. Yoshida, A., Ohmori, Y., Kitano, H., Taguchi-Shiobara, F., Hirano, H., 2012. ABERRANT SPIKELET regulation
of
meristem
D
AND PANICLE1, encoding a TOPLESS-related transcriptional co-repressor, is involved in the fate
in
rice.
Plant
J.
70
(2),
327-339.
PT E
https://doi.org/10.1111/j.1365-313X.2011.04872.x. Yoshida, A., Sasao, M., Yasuno, N., Takagi, K., Daimon, Y., Chen, R., Yamazaki, R., Tokunaga, H., Kitaguchi, Y., Sato, Y., Nagamura, Y., Ushijima, T., Kumamaru, T., Iida, S., Maekawa, M., Kyozuka, J., 2013. TAWAWA1, a regulator of rice inflorescence architecture, functions through the
CE
suppression of meristem phase transition. Proc. Natl. Acad. Sci. U. S. A. 110 (2), 767-772. https://doi.org/10.1073/pnas.1216151110. Zhang, X., Wang, J., Huang, J., Lan, H., Wang, C., Yin, C., Wu, Y., Tang, H., Qian, Q., Li, J., Zhang,
AC
H., 2012. Rare allele of OsPPKL1 associated with grain length causes extra-large grain and a significant yield increase in rice. Proc. Natl. Acad. Sci. U. S. A. 109 (52), 21534-21539. https://doi.org/10.1073/pnas.1219776110. Zhang, Y., Yu, H., Liu, J., Wang, W., Sun, J., Gao, Q., Zhang, Y., Ma, D., Wang, J., Xu, Z., Chen, W., 2016. Loss of function of OsMADS34 leads to large sterile lemma and low grain yield in rice (Oryza sativa L.). Mol. Breed. 36, 147. https://doi.org/10.1007/s11032-016-0578-4. Zhao, Y., 2008. The role of local biosynthesis of auxin and cytokinin in plant development. Curr. Opin. Plant Biol. 11 (1), 16-22. https://doi.org/10.1016/j.pbi.2007.10.008. Zhao, Y., Cheng, S., Song, Y., Huang, Y., Zhou, S., Liu, X., Zhou, D., 2015. The Interaction between Rice ERF3 and WOX11 Promotes Crown Root Development by Regulating Gene Expression Involved
in
Cytokinin
Signaling.
Plant
Cell
27
(9),
2469-2483.
https://doi.org/10.1105/tpc.15.00227. 48