Identification of abiotic stress-responsive genes: a genome-wide analysis of the cytokinin response regulator gene family in rice

ABSTRACT

Response regulators (RRs) are implicated in various developmental processes as well as environmental responses by acting as either positive or negative regulators, and are crucial components of cytokinin signaling in plants. We characterized 36 RRs in rice (Oryza sativa L.; Os) using in silico analysis of publicly available data. A comprehensive analysis of OsRR family members covered their physicochemical properties, chromosomal distribution, subcellular localization, phylogeny, gene structure, distribution of conserved motifs and domains, and gene duplication events. Gene Ontology analysis indicated that 22 OsRR genes contribute mainly to the cytokinin response and signal transduction. Predicted cis-elements in RR promoter sequences related to phytohormones and abiotic stresses indicated that RRs are involved in hormonal and environmental responses, supporting previous studies. MicroRNA (miRNA) target analysis showed that 148 miRNAs target 29 OsRR genes. In some cases, multiple RRs are targets of the same miRNA group, and may be controlled by common stimulus responses. Based on the analysis of publicly available gene expression data, OsRR4, OsRR6, OsRR9, OsRR10, OsRR22, OsPRR73 and OsPRR95 were found to be involved in responses to abiotic stresses. Using quantitative reverse transcription polymerase chain reaction we confirmed that six of these RRs, namely OsRR4, OsRR6, OsRR9, OsRR10, OsRR22 and OsPRR73, are involved in the response to salinity, osmotic, alkaline and wounding stresses, and can potentially be used as models to understand molecular mechanisms underlying stress responsiveness.

INTRODUCTION

Cytokinins are adenine derivative phytohormones (Sakakibara, 2006) that promote cell division and differentiation, shoot development, apical dominance, leaf senescence, chloroplast biogenesis, photomorphogenic development and stress responses (Mok and Mok, 2001; Tanaka et al., 2006; Roef and Van Onckelen, 2010; Jameson and Song, 2016). The cytokinin signal is perceived by a membrane-localized histidine kinase (HK), transmitted downstream via a histidine–aspartate phosphorelay, and ultimately triggers an intracellular molecular response. Cytokinin signaling was initially likened to the two-component system in prokaryotes; however, it was later elucidated as a multistep phosphorelay system in higher plants (Hwang and Sheen, 2001; Heyl and Schmulling, 2003; Qiang et al., 2023). After perception of the cytokinin signal by HK, the conserved C-terminal receiver domain of HK undergoes autophosphorylation. Subsequently, the phosphorylation relay targets histidine-containing phosphotransfer protein (HPt), which in turn facilitates the activation of response regulators (RRs) by HPt. Once phosphorylated, HPt shuttles to the nucleus and phosphorylates aspartate residues of the RRs. Certain subfamilies of RRs function as transcription factors that initiate transcription of their target genes (Stock et al., 2000; Grefen and Harter, 2004; Schaller et al., 2008).

RRs in Arabidopsis thaliana are classified into four sub-groups, type-A, type-B, type-C and pseudo RR (PRR), based on the structure of the core receiver domain and the sequence of the C-terminal domain (Huo et al., 2020). Type-A RRs contain a receiver domain at their N terminus and a short C-terminal extension. Type-B RRs have an N-terminal receiver domain along with a C-terminal extension containing a GARP/Myb-like DNA-binding domain, which is important for the binding and activity of transcription factors (Huo et al., 2020; Jiang et al., 2022). In addition, type-B RRs are very likely located in the nucleus and are transcription factors that directly bind to the promoter of target genes, including type-A RRs, to initiate transcription (Hwang and Sheen, 2001; Sakai et al., 2001; Hosoda et al., 2002; Kieber and Schaller, 2014). Type-C RRs have a similar structure to type-A RRs but are phylogenetically different and unresponsive to cytokinins (Schaller et al., 2008; Huo et al., 2020). Along with a unique CCT domain, PRRs have an atypical pseudo-receiver domain that lacks the conserved aspartate phosphorylation site and are involved in the control of the circadian rhythm (Pareek et al., 2006; Schaller et al., 2008; Wang et al., 2020).

Type-A RRs function as transcriptional repressors. They compete with the type-B RRs for phosphorylation, thereby negatively regulating cytokinin signaling by attenuating the strength of the cytokinin signal (To et al., 2007; Geng et al., 2022). On the other hand, type-B RRs are positive regulators that activate downstream cytokinin-responsive genes (Hutchison and Kieber, 2002). Type-A RRs are considered early cytokinin-responsive genes because their expression increases rapidly in the presence of exogenous cytokinin, whereas type-B RRs show the opposite phenomenon (D’Agostino et al., 2000).

The RRs are reported to be involved in plant developmental processes, transduction of cytokinin signals, and abiotic stress tolerance. For instance, overexpression of a type-A Arabidopsis RR, ARR8, represses shoot regeneration and greening of calli in Arabidopsis (Osakabe et al., 2002). ARR7 and ARR15 are also type-A RRs. It was shown that their loss-of-function mutants regulate embryogenesis by interacting with auxin (Müller and Sheen, 2008; Bielach et al., 2017). ARR16 was reported to be involved in root growth, shoot regeneration and leaf senescence (Ren et al., 2009). In maize (Zea mays L.; Zm), Zeng et al. (2021) demonstrated that the natural variation of ZmRR1 played an important role in imparting cold tolerance. In a recent study, Li et al. (2024) found that silencing of the melon (Cucumis melo L.; Cm) CmRR6 and CmPRR3 genes imparted higher cold tolerance in melon.

In rice (Oryza sativa L.; Os), overexpression of a type-B RR (OsRR22) resulted in reduced plant height (Shi et al., 2020). Hirose et al. (2007) found that overexpression of OsRR6 (a type-A RR) in rice plants resulted in poorly developed roots and reduced shoot regeneration. With respect to stress tolerance, the OsRR22 mutant hst1 acquired salinity tolerance (Takagi et al., 2015). OsRR6 overexpression in Arabidopsis led to increased drought and salinity tolerance (Bhaskar et al., 2021). In addition, the double mutant osrr9/osrr10 (two type-A RRs), generated by genome editing, showed higher salinity tolerance at the seedling stage in rice (Wang et al., 2019).

Thirty-six OsRRs were identified and named in previous studies (Pareek et al., 2006; Schaller et al., 2007). However, recent advances in multiple databases and analysis tools have made a genome-wide comprehensive characterization of the entire gene family possible. Therefore, we conducted a genome-wide study, aimed at characterizing the 36 OsRR genes, including evolutionary relationships, chromosome mapping and subcellular localization, gene structure, conserved motifs and domains, gene duplication, promoter analysis and microRNA (miRNA) analysis. This effort was undertaken to gain a deeper understanding of their evolution and potential roles. Furthermore, expression profiling using publicly available transcriptome data followed by validation through quantitative reverse transcription polymerase chain reaction (qPCR) was performed, and we found several genes that are involved in the response to abiotic stresses. Our findings provide detailed information on all OsRR genes, facilitating further functional analysis.

RESULTS

OsRR gene information and physicochemical properties of OsRR proteins

For our genome-wide analysis of RR genes we used 36 rice genes encoding OsRR proteins that were formerly identified and named (Pareek et al., 2006; Schaller et al., 2007). The basic gene information, synonyms used in different studies (Pareek et al., 2006; Schaller et al., 2007) and physicochemical properties of OsRR proteins are provided in Table 1. The lengths of the OsRR proteins range from 122 (OsRR8, OsRR12 and OsRR13) to 768 amino acids (OsPRR73), with sequences of their genes ranging from 366 to 2,304 bp. The molecular weights of the OsRR proteins range from 13.56 kDa (OsRR41) to 84.08 kDa (OsPRR73), and their theoretical isoelectric points (pIs) range from 4.21 (OsRR7) to 9.1 (OsRR30). Subcellular localization prediction by Plant-mPLoc indicated that OsRR proteins are primarily located in the nucleus except for OsRR18 and OsRR31, which are additionally located in the chloroplast and cytoplasm, respectively. However, because TargetP-2.0 and DeepLoc 2.0 predicted that OsRR18 and OsRR31 reside in the nucleus only, we consider that OsRRs are nucleus-located proteins. We next searched for putative nuclear localization signals (NLSs) on OsRR proteins. All OsRR proteins had a single NLS except OsRR26, OsRR27 and OsPRR73, which possessed two NLSs, and OsRR13, which had none despite being predicted to be located in the nucleus (Table 1). Based on their calculated instability index (IIp) (Gamage et al., 2019), OsRR8, OsRR41, OsRR42, OsRR17, OsRR18, OsRR24, OsRR28, OsRR29, OsRR30 and OsPRR11 are stable (instability index < 40), whereas the other OsRR proteins are predicted to be less stable (instability index > 40) (Table 1).

Table 1. Basic information for OsRR genes and physicochemical properties of OsRR proteins

Gene	Locus ID (MSU-ID)	Synonym(s)	cDNA Start to End	CDS (bp)	P. Loc.	NLS Loc. (AA)	PL (AA)	PMW (kDa)	pI	IIp	GRAVY
OsRR1	LOC_Os04g36070	OsRRA1a, Rr1b	22023831–22027680	696	Nu	173–177	232	25.27	6.58	59.56	0.138
OsRR2	LOC_Os02g35180	OsRRA2a, Rr2b	21128619–21131370	759	Nu	59–63	253	26.86	8.12	53.61	-0.014
OsRR3	LOC_Os02g58350	OsRRA3a, Rr3b	35690390–35689246	396	Nu	74–78	132	14.9	7.72	61.94	-0.222
OsRR4	LOC_Os01g72330	OsRRA4a , Rr4b	41948951–41952388	699	Nu	147–151	233	24.86	5.31	50.55	-0.317
OsRR5	LOC_Os04g44280	OsRRA5a , Rr5b	26227216–26228435	405	Nu	46–51	135	14.82	7.71	54.35	-0.149
OsRR6	LOC_Os04g57720	OsRRA6a, Rr6b	34375977–34377357	513	Nu	115–119	171	18	8.49	58.59	0.098
OsRR7	LOC_Os07g26720	OsRRA7a, Rr7b	15440888–15442180	621	Nu	153–158	207	22.65	4.21	55.90	-0.606
OsRR8	LOC_Os08g28900	OsRRA8a, Rr8b	17677909–17678691	366	Nu	70–75	122	13.68	6.27	38.35	0.107
OsRR9	LOC_Os11g04720	OsRRA9a, Rr9b	2013370–2015405	606	Nu	86–91	202	22.59	5.84	66.30	-0.599
OsRR10	LOC_Os12g04500	OsRRA10a, Rr10b	1915855–1917598	606	Nu	86–91	202	22.59	5.84	66.30	-0.599
OsRR11	LOC_Os02g42060	OsRRA11a, Rr11b	25287918–25286936	444	Nu	99–103	148	15.46	8.82	50.01	0.212
OsRR12	LOC_Os08g28950	OsRRA12a, Rr12b	17714142–17714938	366	Nu	70–75	122	13.63	5.89	43.14	0.180
OsRR13	LOC_Os08g26990	OsRRA13a, Rr13b	16489492–16490328	366	Nu	Not found	122	13.72	6.02	46.65	0.136
OsRR41	LOC_Os03g53100	OsRRA14a, Rr41b	30451517–30452854	384	Nu	73–77	128	13.56	5.17	30.45	0.053
OsRR42	LOC_Os04g13480	OsRRA15a, Rr42b	7521475–7522295	399	Nu	125–129	133	14.54	5.20	32.26	0.031
OsRR17	LOC_Os04g28120	OsRRA17a, Prr12b	16602470–16605138	759	Nu	79–83	253	27.77	4.95	26.92	-0.388
OsRR18	LOC_Os05g32890	OsRRA18a, Prr10b	19265596–19261579	1107	Nu	151–196	369	40.36	8.2	32.68	-0.204
OsRR21	LOC_Os03g12350	OsRRB1a, Rr21b	6512743–6518792	2076	Nu	173–187	692	73.81	6.06	46.47	-0.399
OsRR22	LOC_Os06g08440	OsRRB4a, Rr22b	4137247–4142663	2091	Nu	181–191	697	75.95	6.01	48.51	-0.406
OsRR23	LOC_Os02g55320	OsRRB5a, Rr23b	33864842–33859824	2067	Nu	179–190	689	74.38	6.45	49.99	-0.323
OsRR24	LOC_Os02g08500	OsRRB2a, Rr24b	4578363–4573805	1881	Nu	170–197	627	68.40	6.06	35.29	-0.602
OsRR25	LOC_Os06g43910	OsRRB3a, Rr25b	26450760–26454696	2085	Nu	126–179	695	76.62	5.65	46.09	-0.544
OsRR26	LOC_Os01g67770	OsRRB6a, Rr26b	39391101–39387229	1749	Nu	49–57; 165–174	583	64.81	5.04	52.76	-0.426
OsRR27	LOC_Os05g32880	OsRRA16a, Rr27b	19260019–19255490	1542	Nu	128–135; 153–159	514	54.59	5.45	42.12	-0.508
OsRR28	LOC_Os04g28160	OsRRA22a, Rr28b	16629916–16633083	1143	Nu	221–226	381	42.77	5.30	31.98	-0.499
OsRR29	LOC_Os04g28130	OsRRB7a, Rr29b	16614562–16618026	1173	Nu	167–171	391	43.33	5.57	37.14	-0.519
OsRR30	LOC_Os10g32600	Rr30b	17077979–17076097	504	Nu	137–141	168	18.68	9.10	32.67	-0.417
OsRR31	LOC_Os08g35650	Rr31b	22479169–22475055	1458	Nu	176–200	486	52.90	4.43	46.38	-0.169
OsRR32	LOC_Os08g17760	Rr32b	10879442–10875231	1713	Nu	123–147	571	60.70	4.42	49.37	-0.088
OsRR33	LOC_Os0g35670	OsRRA19a, Rr33b	22497789–22492747	1845	Nu	178–194	615	64.79	4.44	48.45	-0.080
OsPRR1	LOC_Os02g40510	OsPRR1a, Prr1b	24569293–24572560	1557	Nu	445–486	519	57.88	5.80	51.29	-0.665
OsPRR59	LOC_Os11g05930	OsPRR5a, Prr59b	2788995–2793735	2100	Nu	633–686	700	76.09	8.11	51.48	-0.635
OsPRR37	LOC_Os07g49460	OsPRR4a, Prr37b	29616704–29629223	2229	Nu	664–727	743	79.92	8.87	44.83	-0.848
OsPRR73	LOC_Os03g17570	OsPRR3a, Prr73b	9768690–9759478	2304	Nu	312–327; 705–742	768	84.08	5.98	48.01	-0.989
OsPRR95	LOC_Os09g36220	OsPRR2a, Prr95b	20889843–20885171	1872	Nu	562–616	624	69.62	6.62	58.17	-0.792
OsPRR11	LOC_Os04g28150	Prr11b	16624885–16627577	867	Nu	69–73	289	31.83	4.82	34.26	-0.522

^a Synonym used by Pareek et al. (2006); ^b synonym used by Schaller et al. (2007).

Proteins underlined in column IIp correspond to stable proteins.

cDNA, complementary DNA; CDS, coding sequence; bp, base pairs; P. Loc, protein location; Nu, nucleus; NLS Loc., location of nuclear localization signal; AA, amino acid; PL, protein length; PMW, protein molecular weight; pI, isoelectric point; II_p, protein instability index; GRAVY, grand average of hydropathy.

Phylogenetic relationships

In Arabidopsis, 12 type-B RRs (ARR1, 2, 10–14, 18–21 and 23), 12 type-A RRs (ARR3–9, 15–17, 22 and 24), and nine PRRs (APRR1–9) have been identified (D’Agostino and Kieber, 1999; D’Agostino et al., 2000; Nakamichi et al., 2004; Schaller et al., 2008; Ren et al., 2009). In maize, Chu et al. (2011) and Wang et al. (2022) identified 21 type-A RRs (1–21), seven type-B RRs (22–28) and nine PRR genes. To illustrate the evolutionary relationships across species, a circular phylogenetic tree was built based on the full-length protein alignment of 36 rice OsRRs, 37 maize ZmRRs and 33 Arabidopsis ARRs (Supplementary Fig. S1). The RR genes were classified into three groups, with Group I being the largest, consisting of 39 type-A and four type-C members (highlighted in red). Group II and Group III include 40 type-B and 23 PRR genes, respectively. Unlike the ARRs, OsRR proteins are closely related to ZmRRs, suggesting that RRs in rice and maize share similar functions (Qiang et al., 2023).

A rectangular phylogenetic tree was generated based on the alignment of 36 OsRR protein sequences using the neighbor-joining method with 1,000 bootstrap replications (Fig. 1A). OsRR genes are clustered into six groups, highlighted in different boxes of different colors. Clusters I, II and III consist of type-A, type-C and PRR genes, respectively, whereas type-B OsRR genes are grouped into three clusters, IV, V and VI (Fig. 1A). Cluster I is the largest, with 13 type-A members, whereas Cluster II and Cluster V are the smallest, with only two members.

Fig. 1. Phylogenetic tree, gene structure, conserved motif composition and conserved domains of OsRR genes and OsRR proteins. (A) The unrooted neighbor-joining phylogenetic tree was constructed following the maximum likelihood method with 1,000 bootstrap replications using MEGA 7, which distributed OsRR proteins into six groups highlighted in different background colors. (B) Exon–intron organization of OsRR genes. Exons, introns and untranslated regions are represented by pink boxes, black lines and green boxes, respectively. CDS, coding sequence. (C) Distribution of 20 conserved motifs in OsRR proteins, represented by different colors. Amino acid sequences of these motifs are summarized in Table 2. (D) Composition of conserved domains in OsRR family members.

Gene structure, conserved motif and protein domain analysis of the OsRR gene family

Investigation of gene structure identified that OsRR genes contain at least one intron. The number of exons ranges from two to eight. Exons of the type-A OsRR genes are smaller than those of the type-B genes and PRRs. Almost all information on mRNA and exon–intron structures was retrieved from the Phytozome database. However, we could not find any predicted untranslated sequences in OsRR18 or OsRR31 (Fig. 1B), or any evidence of their transcription from microarray and RNA-sequencing (RNA-seq) data, in the databases. The distribution of conserved motifs was explored using MEME Suite, and their composition was revealed using Pfam (Table 2). Sequence logos of the motifs are shown in Supplementary Fig. S2. The degree of conservation can be determined from the total height of each logo position. Almost all OsRR proteins have one or more copy of motifs 1, 3 and 6, which may constitute the signal transduction receiver domain. Motifs 4 and 7 are only found in type-B and pseudo-type OsRRs, respectively, and correspond to Myb-type DNA-binding and CCT domains, respectively (Fig. 1C; Table 2). The 13 types of conserved domains present across OsRR proteins are displayed in Fig. 1D. The receiver domains are distinct and separate in type-A, type-B, type-C and pseudo-type OsRRs. Myb-like DNA-binding and CCT domains are exclusive to type-B and pseudo-type OsRRs, respectively. However, OsRR31 and OsRR32 do not have any identifiable domains (Fig. 1D).

Table 2. Size and description of conserved motifs of OsRR proteins

Motif	Peptide sequence(s)	Length*	Description
1	YCMPEMTGYDLLKKIKESSELKDIPVVIMSSENVPSRISRC YWMPEMTGYDLLKKIKESSYCKHIPVIIMSSENIPTRIFRC	41	Response regulator receiver domain
2	GAEDFLLKPVRIEDL	15
3	LRVLAVDDDPVDRKVIEALLR	21	Signal transduction response regulator receiver domain
4	RVVWSVELHRKFVAAVNQLGIDKAVPKKILELMNVEKLTRE RVVWSVELHRKFVAAVNQLGIDKAVPKKILELMNVEYLTRE	41	Transcription factor LUX/BOA-like (Myb like)
5	CNYRVTTVDSGKKALEFL	18	Unknown
6	RENAEDFDLVJSDVHMPDMDGFKFL	25	Signal transduction response regulator receiver domain
7	SRREAALNKFRLKRKDRCFEKKVRYQSRKKLAEQRPRVRGQFVRQ CHREAALNKFRQKRKERCFEKKVRYQSRKKLAEQRPRVRGQFVRQ	45	Zinc finger protein CONSTANS-like (CCT like)
8	KNJWQHVWRKKL	12	Unknown
9	ELVGLEMDLPVIMLSADGETETVMKGVTH	29	Unknown
10	NVASHLQKYRLYLKR	15	Unknown
11	PMEDAGGEAADGSLNIGEGGMEIGWDLDLDDILMNNTNEFAFLDDLAWIE PMEDAGGEAADGSLNIEEGGMEIGWDLDLDDILMNNTNEFAFLDDLAWIE	50	Unknown
12	SEPBVNMIITD	11	Unknown
13	RCKLDYQQEQNKPSNADSDNSSNPTSCGSSDQTGRNSHKRKEVDEEILPE RCKEHYQQECHKPPNAESDHSHNPTTCGSSDQTGRNSHKRKEYDEEIEPH	50	Unknown
14	MTMNKGKAPMIELPFGLPVDDFLVGQTAYGGAGPSIG MTMNKGKAPMIELPFGMPVDDFLVGQTAYGGAGPSIG	37	Unknown
15	MEKIKDLLQGIGDESTCANEMNSFPENPKDGTKKKYYLMWTPHLQKKFLH MEKIKDMLQGIGDESTCANEMNSFPENPKDGTKKKYYLMWTPHLQKKFLH	50	Unknown
16	EHDDVYDAMRRALQYGTVFDESKYSSDPFSNEDE EYDDVYTAMRRALQYGTIFDESKYSSDPCSNEDE	34	Unknown
17	QSSWTKRAVEIDSPQQMSPDQPADPPDSTCAQVIHPKSEICSNRWLPC QSSWTKRAVEIDSPQQMSPDQPADPPDSTCAQVIHPKSEICSNRWLPC	48	Unknown
18	EFLNHLLLKATYIVRKPLDPAVMARLWRVVAWRRCCLEERI EFLNHLLLKATYIVRKPLDPTVMARLWRVVAWRRCCLEERI	41	Unknown
19	DGNDDDVVIIEEPQVHFKVVRSRGSRKRQ	29	Unknown
20	DAAAIYQYTNALSNNNAVGSLMVPPI	26	Unknown

* Number of amino acids.

Chromosome mapping and gene duplication

The 36 OsRR genes are unevenly scattered over the 12 chromosomes of rice (Fig. 2A), with chromosome 4 containing the highest number, eight OsRR genes. However, six OsRR genes are located on each of chromosomes 2 and 8. Chromosomes 9, 10 and 12 possess only one OsRR gene each (OsPRR95, OsRR30 and OsRR10, respectively) (Fig. 2A). Fourteen duplicated gene pairs were observed with 10 segmental duplications and four tandem duplications (Table 3). This suggests that segmental duplication has played a more important role in broadening the gene family (Fig. 2B; Table 3). Subsequently, we computed the ratios of nonsynonymous (Ka) and synonymous (Ks) substitutions (Ka/Ks ratios) to ascertain the evolutionary restrictions imposed on duplicated gene pairs. Ka/Ks ratios varied from 0.217 (OsRR12, OsRR8) to 1.049 (OsPRR95, OsRR27). Except for the OsPRR95, OsRR27 gene pair, the Ka/Ks ratios were less than 1, which suggests that these genes evolved under the action of negative or purifying selection, while OsPRR95, OsRR27 evolved under positive selection (Table 3).

Fig. 2. Chromosome mapping and gene duplication events of OsRR genes. (A) The physical location of 36 OsRR genes on 12 chromosomes of rice. Chromosome length corresponds to the vertical bar on the left side, measured in megabases (Mb). (B) Duplication events of OsRR genes. The measurement scale on the circle is in megabases.

Table 3. Gene duplication events of OsRR genes in rice

	Gene pair	Ka	Ks	Ka/Ks	Time (MYA)	Type of duplication	Purifying selection
1	OsRR17, OsRR18	0.380	0.644	0.591	49.55	Segmental	Yes
2	OsRR10, OsRR9	0.003	0.013	0.224	0.99	Segmental	Yes
3	OsRR28, OsRR29	0.074	0.108	0.686	8.30	Tandem	Yes
4	OsRR26, OsRR23	0.319	0.513	0.622	39.42	Segmental	Yes
5	OsRR4, OsRR6	0.126	0.458	0.274	35.22	Segmental	Yes
6	OsRR25, OsRR22	0.285	0.494	0.578	38.00	Tandem	Yes
7	OsRR7, OsRR24	0.275	0.408	0.675	31.35	Segmental	Yes
8	OsRR32, OsRR33	0.133	0.214	0.622	16.47	Tandem	Yes
9	OsRR2, OsRR1	0.108	0.281	0.384	21.61	Segmental	Yes
10	OsRR41, OsPRR11	0.243	0.383	0.634	29.45	Segmental	Yes
11	OsPRR95, OsRR27	0.497	0.474	1.049	36.43	Segmental	No (Positive)
12	OsPRR73, OsPRR37	0.111	0.289	0.383	22.22	Segmental	Yes
13	OsRR21, OsPRR1	0.438	0.518	0.845	39.86	Segmental	Yes
14	OsRR12, OsRR8	0.009	0.043	0.217	3.27	Tandem	Yes

Ka, Nonsynonymous substitutions; Ks, synonymous substitutions; MYA, million years ago.

Time = Ks/2λ, where λ= 6.5 x 10^-9

cis-elements in promoters

cis-acting elements are regulatory regions in a promoter where transcription factors bind and control gene expression. In silico analysis of cis-elements in the promoter regions (1,500 bp upstream) of OsRR genes identified cis-elements associated with responses to different phytohormones (gibberellins, abscisic acid (ABA), methyl jasmonate (MeJa), auxin and salicylic acid) and stresses (drought, low temperature, wound and defense stress), as shown in Fig. 3. cis-elements for ABA and MeJa were comparatively more abundant. The number of cis-elements was relatively higher in OsPRR1 (14 cis-elements) and OsRR5 (13 cis-elements). In contrast, only one element was found in OsRR3, OsRR22 and OsPRR11. This prediction result suggests that the transcription of OsRR genes is linked to hormonal crosstalk with abiotic stress responses.

Fig. 3. Phytohormone and stress-responsive cis-acting elements in OsRR promoters. PlantCARE was used to identify the cis-elements present in the 1,500-bp sequence immediately upstream of the transcription start site of each OsRR gene. TBtools was used to visualize the elements. Colored dots denote various elements as indicated to the right. The scale bar at the bottom corresponds to the location of elements within the promoter region. Some dots overlap due to the close location of the elements.

Gene Ontology (GO) analysis

GO is a framework designed to represent biological knowledge about gene products. It is divided into three primary categories: biological processes (BPs), molecular functions (MFs) and cellular components (CCs). To predict the function of OsRR genes we carried out a GO enrichment analysis (Table 4). GO enrichment analysis is a method for extracting groups of genes with highly frequent GO terms that is usually applied to differentially expressed gene sets from transcriptomics or proteomics. In this experiment, we used this method to determine the number of ‘child’ genes in each of the GO term ‘parent’ gene populations. The PANTHER and ShinyGO GO analysis tools considered OsRR9 and OsRR10 as the same protein, so the output is the result of GO analysis of 35 proteins in total. In biological processes, 31 out of 35 genes were assigned to the parent GO terms of cellular response to stimulus/intracellular signal transduction, and all 31 genes were assigned to the parent term of phosphorelay signal transduction system. A total of 74 genes were assigned to the parent term of phosphorelay signaling pathway, 31 of which were OsRR genes. In molecular functions, eight OsRR genes were assigned to the parent term of DNA-binding transcription factor activity; on the other hand, in cellular compartments, 33 OsRR genes were assigned to the parent term of nucleus. Supplementary Fig. S3 illustrates the interaction networks of the enriched GO terms, demonstrating their interrelations and functional interactions.

Table 4. GO analysis of OsRR genes

	GO term	No. of genes in group of		FDR	(Minus) log(FDR)	Fold enrichment
		the child term	the parent term
BPs	Cellular response to stimulus	31	1,898	3.68E-39	38.43	18.85
	Intracellular signal transduction	31	413	1.96E-59	58.71	86.62
	Signaling	31	1,111	2.45E-46	45.61	32.20
	Signal transduction	31	1,100	2.09E-46	45.68	32.52
	Phosphorelay signal transduction system	31	74	5.48E-85	84.26	483.45
	Hormone-mediated signaling pathway	22	396	1.37E-35	34.86	64.11
	Cytokinin-activated signaling pathway	22	43	7.54E-59	58.12	590.44
	Response to cytokinin	22	67	1.16E-53	52.94	378.94
	Cellular response to cytokinin stimulus	22	44	1.13E-58	57.95	577.02
	Cell communication	31	1,271	1.46E-44	43.83	28.15
MFs	Nucleic acid binding	9	3,649	0.0044399	2.35	2.85
	DNA binding	9	2,085	9.52E-05	4.02	4.98
	Transcription regulator activity	8	1,011	5.25E-06	5.28	9.13
	DNA-binding transcription factor activity	8	897	4.30E-06	5.37	10.29
CCs	Intracellular anatomical structure	18	10,370	8.76E-03	2.06	2.43
	Membrane-bound organelle	18	8,078	4.66E-04	3.33	3.12
	Intracellular organelle	18	8,820	1.26E-03	2.90	2.87
	Intracellular membrane-bound organelle	18	8,057	6.74E-04	3.17	3.13
	Nucleus	33	4,413	7.51E-08	7.12	5.80

BPs, biological processes; MFs, molecular functions; CCs, cellular components; FDR, false discovery rate.

Total number of genes in parent term = 43,659.

Prediction of miRNAs

miRNAs are short non-coding RNAs (21 bp long on average) that play important roles in plant development, signaling, reproduction and adaptation to various biotic and abiotic stresses (Hu et al., 2013; Yang et al., 2017). They typically function as negative regulators of gene expression by inhibiting translation of their target mRNA through processes such as cleavage or translational inhibition. In miRBase release 22.1 (http://www.mirbase.org/), accessed on 3 March, 2024, 713 mature published miRNAs were identified in rice. The interaction network between these miRNAs and their target OsRR genes is shown in Fig. 4A. Supplementary Table S1 shows miRNAs, their sequences and targeted OsRR genes, and their mode of inhibition of translation. Supplementary Table S2 shows the number of interacting miRNAs targeting OsRR genes. A total of 148 19–24-bp miRNAs target 29 OsRR genes. All OsRR genes are targets of miRNA except for OsRR5, OsRR9, OsRR10, OsRR13, OsRR17, OsRR42 and OsRR30. OsPRR37 is targeted by the highest number of miRNAs (19 miRNAs), followed by OsRR8 and OsRR12 (18 miRNAs each) (Fig. 4A, Supplementary Table S1 and S2). OsRR3, OsRR11, OsRR28 and OsRR31 are targeted by a single miRNA, namely osa-miR2863c, osa-miR3979-5p, osa-miR2275c and osa-miR1865-5p, respectively. Cleavage was the most prevalent inhibitory mechanism compared to translation inhibition (Fig. 4A, Supplementary Table S1).

Fig. 4. Prediction of miRNAs and protein–protein interactions among OsRRs. (A) Interaction network between putative miRNAs and their targeted OsRR genes. Yellow circles represent the miRNAs, while different colored ovals indicate the different OsRR genes. Connecting lines represent interactions between them. (B) Prediction of protein–protein interactions among two-component system proteins using STRING at high confidence (0.70). Each node represents a protein and its name is displayed on the right. Different colored lines represent sources of evidence for interaction. Disconnected nodes represent proteins for which no data are available to support interaction with any other of these proteins.

Prediction of protein–protein interactions

Cytokinin signaling takes places through a multistep phosphorelay two-component system (TCS) that involves interactions among signaling proteins to respond to signal stimuli. Signaling proteins responsible for the multistep phosphorelay TCS encoded in the rice genome comprise 11 histidine kinases (OsHKs), five histidine phosphotransferases (OsHPs) and 36 RRs (OsRRs) (Sharan et al., 2017). In this study, we used the STRING database to predict in silico protein–protein interactions among TCS proteins at high confidence (0.70) with the support of available data. Putative interactions between signal receiver OsHKs (OsHK1, OsHK2, OsHK3, OsHK4, OsHK5 and OsHK6) and OsHPs (OsAHP1, OsAHP2, OsPHP1, OsPHP2 and OsPHP3) were predicted and four OsHPs (OsAHP1, OsAHP2, OsPHP2 and OsPHP3) were predicted to interact with type-B RRs (OsRR21, OsRR22, OsRR23, OsRR24, OsRR26 and OsRR27). OsPRRs were predicted to interact only with their own group members (Fig. 4B).

Expression profiling of OsRR genes

Using publicly available RNA-seq and microarray data, we established expression profiles of OsRR genes under salinity, drought, alkaline and wound stress conditions, as described in Materials and Methods. The results show that OsRR10 had the highest expression in control shoots of the IR28 genotype at 6 h, followed by the CSR28 genotype in the same conditions. OsPRR95 showed moderately high expression in shoots of CSR28 and IR28 at 6 h of salinity treatment (Fig. 5A). In the case of drought stress, OsPRR95, OsPRR73 and OsPRR37 were upregulated in both ZH11- (susceptible) and SNAC1-overexpressing (tolerant) plants under normal and drought conditions (Fig. 5B). During alkaline stress treatment, OsRR21 and OsPRR95 showed the highest expression in both sensitive (CD) and tolerant (WD) cultivars. OsRR4 and OsRR6 produced slightly higher transcript levels in the sensitive cultivar, CD, exposed to alkali treatment (Fig. 5C). With respect to wound stress, OsRR11 and OsRR6 showed the highest expression at 2 h and 1 h after wounding in Nipponbare rice leaves, respectively (Fig. 5D).

Fig. 5. Transcription profiles of OsRR genes under different abiotic stresses. RNA-seq and microarray data were gathered from the publicly available NCBI GEO dataset. (A) Expression patterns of OsRR genes in two contrasting genotypes (CSR28, tolerant and IR28, susceptible), two treatments (control and 150 mM salinity), two tissues (root and shoot) and two sampling time points (6 h and 54 h). (B) Expression levels of OsRR genes in ZH11 wild type (ZH11; susceptible) and SNAC1-overexpressing (SNAC1; tolerant) plants under normal and moderate drought stress conditions. (C) Expression analysis of OsRR genes under alkaline stress. Transcriptome data were compared between an alkaline-tolerant [WD20342 (WD)] and -sensitive [Caidao (CD)] rice cultivar under control and alkaline stress (suffix T) conditions. (D) Expression profiles of OsRR genes under wounding stress. Wild-type rice leaves wounded for 30, 60, 120 or 240 min were compared with untreated control leaves. GEO accession numbers and computation of gene expression values are described in Materials and Methods. Heatmap Illustrator in TBtools was used to generate heatmaps, and the level of gene expression is represented by color scales at the side of each panel.

Among the OsRRs identified as stress response genes in the above analysis, OsRR4, OsRR6, OsRR9, OsRR10, OsRR22 and OsPRR73 were selected and their RNA-seq and microarray data were verified by qPCR analysis. OsRR9 was not seen to be induced under any stress conditions in microarray or RNA-seq data. However, as its coding sequence is 99% identical to OsRR10, we included it in the qPCR analysis to account for possible interference of OsRR9 transcripts with OsRR10 expression data. In addition, we included OsRR22 in the qPCR-tested genes because its mutant, hitomebore salt tolerant (hst1), is a well-known salinity tolerance mutant (Takagi et al., 2015), which was later confirmed using genome-edited rice (Zhang et al., 2019).

We subjected 14-day-old Nipponbare rice seedlings to salinity (250 mM NaCl), osmotic (25% polyethylene glycol (PEG)), alkaline (0.5% Na₂CO₃) and wounding stress for 0, 1, 6 and 12 h and selected OsRR4, OsRR6, OsRR9, OsRR10, OsRR22 and OsPRR73 to validate RNA-seq and microarray data by qPCR (Fig. 6).

Fig. 6. qPCR analysis of changes in OsRR gene expression under different abiotic stress conditions. The stress conditions from the left were 250 mM NaCl (for salinity stress), 25% PEG (for osmotic stress), 0.5% Na₂CO₃ (for alkaline stress) and wounding, for 0, 1, 6 and 12 h in 14-day-old Nipponbare rice seedlings. Statistical analysis using RStudio software (V1.3.1093, https://www.r-project.org/) included two-way analysis of variance and mean separation at P < 0.05 by Tukey’s HSD test. Error bars indicate mean ± standard error. Mean values accompanied by the same letter are not significantly different.

Under salinity stress, OsRR4 exhibited a significant upregulation at 6 h after treatment, which persisted to 12 h (Fig. 6, panel OsRR4, NaCl). OsRR6 started to produce significantly more transcripts at 6 h after osmotic stress treatment and continued to do so until the end of the treatment (12 h) (Fig. 6, panel OsRR6, PEG). Despite OsRR9 and OsRR10 sharing 99% identity in their coding sequences, they behaved differently under different abiotic stresses. Both OsRR9 and OsRR10 showed similar expression patterns in wounded plants, i.e., increased gene expression at the beginning of treatment (1 h) followed by a significant decline (Fig. 6, panel OsRR9 and OsRR10, wounding). However, during osmotic stress, OsRR10 produced significantly higher transcript levels during the entire treatment period (Fig. 6, panel OsRR10, PEG). OsRR22 responded to osmotic and alkaline stresses between 6 and 12 h of treatment (Fig. 6, panel OsRR22, PEG and Na₂CO₃). Expression of OsPRR73 was highly induced after 6 h of different stresses and significantly downregulated thereafter (Fig. 6, panel OsPRR73, NaCl, PEG and Na₂CO₃). However, it was significantly upregulated in wounded plants after 1 h of treatment (Fig. 6, panel OsPRR73, wounding). In summary, we observed substantial overlap between the transcriptome analysis and the qPCR data, with only a few deviations. Variations in the experimental setup and genotype are likely the major factors contributing to these exceptions.

DISCUSSION

To sense and acclimate to adverse environments, plants employ a complicated circuitry of signaling components (Zwack and Rashotte, 2015). Cytokinins are multifaceted phytohormones that carry out a variety of tasks related to plant growth and development as well as to responses to abiotic stresses. Cytokinins act by positively or negatively regulating gene transcription (Imamura et al., 1999; Ramírez-Carvajal et al., 2008). The whole genomes of many important crops have been sequenced and their sequences are publicly available, making it easier to conduct genome-wide in silico analysis of gene families to decipher features, structure, potential functions and other details, as well as to select candidate genes for functional analysis. In this study we provide a comprehensive description of the genome-wide characterization of 36 OsRR genes. This included analyses of phylogeny, gene structure, motif and domain compositions, chromosome localization, gene duplication, promoter analysis, GO analysis and miRNA analysis. Subsequently, transcriptome analyses under various abiotic stresses were carried out using publicly available databases, followed by validation by qPCR.

The physicochemical properties of OsRR proteins exhibit considerable variation. As pointed out by Zeng et al. (2017), subcellular localization of proteins is important because of the close relationship between their location and function. All of the OsRRs are found in the nucleus, where they can interact with other transcription factors (Table 1). Around 75% of OsRR proteins are predicted to be acidic, with a pI < 7.0. Most of the OsRRs (78%) are predicted to have a GRAVY (grand average of hydropathy) value < 0, showcasing their hydrophilic nature, in agreement with the findings of other studies (Geng et al., 2022; Qiang et al., 2023; Li et al., 2024) (Table 1). The structural diversity of genes and motifs helps to unravel the evolution of gene families (Sun et al., 2022). Exon–intron organization and distribution of conserved motifs and domains exhibit similarities within the same subfamily but vary across subfamilies of OsRR genes, suggesting functional diversification within each subfamily (Fig. 1) as previously proposed by Zhang et al. (2022) and Qiang et al. (2023).

Gene duplication contributes to expanding the gene family and consequently affects the function of the genes (Li et al., 2021). Usually, two types of gene duplication occur in nature: segmental and tandem. Segmental or whole-genome duplication takes place between chromosomes, while tandem duplications include the duplication of single or multiple genes on the same chromosome. The OsRR gene family has expanded through both segmental and tandem duplications. However, segmental duplication was found to be predominant (Fig. 2B, Table 3), which suggests that segmental duplication was the main driving force behind the evolution of the OsRR gene family. With the most recent duplication estimated to have happened about 0.99 million years ago, the OsRR gene family is an ancient gene family (Table 3). In population genetics, the Ka/Ks ratio serves as a powerful tool to assess natural selection acting on protein-coding genes. Ka/Ks ratios > 1, equal to 1 and < 1 indicate accelerated change through positive selection, neutrality and purifying selection, respectively. In our study, Ka/Ks ratios of all duplicated gene pairs (except OsPRR95, OsRR27) were < 1, indicating that these genes underwent evolution under purifying or stabilizing selection, and that their functions are relatively conserved (Table 3).

cis-elements are important for the regulation of gene transcription. We analyzed the 1,500 bp upstream of each of the OsRR genes using the PlantCARE web server to identify potential cis-elements. Different hormone-responsive and abiotic stress-related elements were identified. Among them, ABA, MeJA and drought-responsive elements were more abundant (Fig. 3). These elements contribute to the regulation of OsRR genes in growth, development and response to suboptimal environmental stimuli through phytohormone signaling, a viewpoint supported by Zhang et al. (2022). GO analysis revealed that cytokinin-activated signaling and DNA-binding transcription factor activity are the primary functions of the OsRR genes. So far, out of 74 rice genes involved in phosphorelay signaling pathways, 31 genes belong to the cytokinin response regulator family, suggesting that the OsRR genes constitute a major class of phosphorylation signaling pathways in rice cells. If additional phosphorelay signaling pathways are discovered in future, the proportion of OsRR genes will be changed. Of the 33 OsRR nuclear-localized proteins, eight function as DNA-binding transcription factors, while the other 25 OsRRs may have other functions or unknown DNA-binding domains (Table 4).

According to several studies (Nagasaki et al., 2007; Zhu et al., 2009; Xie et al., 2012; Fang et al., 2014), miRNAs participate in numerous biological processes in plants, such as leaf development, pattern formation, stress responses and sexual reproduction. In this study, a total of 29 OsRR genes were predicted to be the targets of 148 miRNAs (Fig. 4A, Supplementary Tables S1 and S2). Interestingly, some miRNAs target more than one gene, such as osa-miR414 (OsRR1, OsRR7, OsRR23, OsRR24, OsRR27 and OsPRR37), osa-miR1439 (OsRR4, OsRR27, OsRR32, OsRR33 and OsPRR37) and osa-miR2919 (OsRR1, OsRR21, OsRR24, OsRR32 and OsPRR37) (Fig. 4A, Supplementary Tables S1 and S2). Mangrauthia et al. (2017) reported that osa-miR1439 has a specific role in high-temperature stress tolerance. According to Peng et al. (2020), osa-miR414 delays flowering time by downregulating the expression of LOC_Os 05g51830. Additionally, Macovei and Tuteja (2012) found that osa-miR414 targets the PDH45 (Peas DNA helicase 45) gene under saline stress. Also, miR2919 was verified to play a role in drought stress (Kumar et al., 2023).

It is well known that type-B RRs function as transcription factors initiating the transcription of their target genes, including type-A RRs (Stock et al., 2000; Grefen and Harter, 2004; Schaller et al., 2008). Type-B RRs possess the GARP/Myb-like DNA-binding domain, which is important for the binding and activity of transcription factors. In our predicted protein–protein interaction model shown in Fig. 4B, we found interaction between HPs and some type-B RR proteins (OsRR21, OsRR22, OsRR23, OsRR24, OsRR26 and OsRR27), which have a Myb-like DNA-binding motif (motif 4) and/or domain (Fig. 1C; Fig. 1D). Type-A OsRRs do not show any interaction with HPs or other proteins because they are transcriptionally activated by type-B OsRRs.

From an analysis of publicly available transcriptome data, we found a diverse picture of expression profiles under different stresses, revealing the diverse functions of OsRR genes. After thorough consideration, we selected representatives of upregulated genes to study their expression by qPCR and verify the transcriptome data. According to qPCR analysis, all of the selected genes were upregulated in at least one of the abiotic stresses we examined. Significant upregulation observed at 6 h exposure to stress, and continuing to 12 h, was notable for OsRR4 (salinity), OsRR6 (osmotic) and OsRR22 (osmotic), indicating that they are not immediate stress-responsive genes (Fig. 6). Therefore, their expression may continue responding in later stages of stress exposure, consistent with the findings of Li et al. (2024) in melon. OsRR10 and OsPRR73 showed notable upregulation in response to multiple stresses, suggesting that they share some common pathways of different stresses responses via hormone crosstalk. Analysis of cis-elements in OsRR promoters revealed that OsRR9, OsRR10, OsRR22 and OsPRR73 possess drought-responsive elements and, as a result, are upregulated in drought stress (Fig. 3, Fig. 6). A correlation between cis-elements and the expression of RR genes was also reported in earlier studies (Zhang et al., 2022; Li et al., 2024).

CONCLUSION

Although cytokinins were originally studied as growth regulators, they were later discovered to be involved in stress tolerance in plants. In this study, we conducted a comprehensive characterization of 36 cytokinin response regulator genes in rice, laying the groundwork for a better understanding of their fundamental characteristics and roles in abiotic stress tolerance. Moreover, this study enabled us to pinpoint seven genes, namely OsRR4, OsRR6, OsRR9, OsRR10, OsRR22, OsPRR73 and OsPRR95, that should serve as models for functional analysis under various abiotic stresses and shed light on cytokinin-mediated stress response mechanisms.

MATERIALS AND METHODS

Sequence retrieval and database searching

Pareek et al. (2006) conducted a study in which they identified 32 response regulator genes in rice. Subsequently Schaller et al. (2007) expanded the number to 36 genes, renamed them, and classified them into four groups, namely type-A, type-B, type-C and PRRs, based on the structure of the conserved core receiver domain. In our study, the formerly identified genes were individually verified in Pfam (http://pfam.xfam.org/), HMMER (https://www.ebi.ac.uk/Tools/hmmer/) and the NCBI conserved domain database (https://www.ncbi.nlm.nih.gov/Structure/bwrpsb/bwrpsb.cgi) for analysis of the receiver domain (Pfam accession number: PF00072), Myb-like DNA-binding domain (PF00249) and CCT motif (PF06203). Gene ID, chromosome number, coding sequences, promoters and protein sequences were obtained from the Phytozome (https://phytozome-next.jgi.doe.gov/) (Goodstein et al., 2012) and NCBI (https://www.ncbi.nlm.nih.gov/) databases. All rice DNA sequences used in this study are based on the Nipponbare reference sequence, Os-Nipponbare-Reference-IRGSP-1.0 (International Rice Genome Project and Sasaki, 2005), latest update version released on 11 January, 2024. We retrieved information regarding protein molecular weight, pI and IIp from the Expasy online server (https://www.expasy.org/) (Duvaud et al., 2021). Subcellular localization of OsRR proteins was predicted using the Plant-mPLoc online tool (http://www.csbio.sjtu.edu.cn/bioinf/plant-multi/) (Chou and Shen, 2010) (accessed on 8 May, 2024), TargetP-2.0 (https://services.healthtech.dtu.dk/services/TargetP-2.0/) (Almagro Armenteros et al., 2019) and DeepLoc 2.0 (https://services.healthtech.dtu.dk/services/DeepLoc-2.0/) (Thumuluri et al., 2022). NLSs were predicted using the INSP (http://www.csbio.sjtu.edu.cn/bioinf/INSP/) (Guo et al., 2020) online tool (accessed on 8 May, 2024).

Multiple sequence alignment and phylogenetic relationships

TAIR (https://www.arabidopsis.org/) (Berardini et al., 2015) and maizeGDB (https://maizegdb.org/) (Woodhouse et al., 2021) databases were exploited to retrieve the RR protein sequences of Arabidopsis and maize, respectively. An unrooted tree was created by aligning the full-length RR protein sequences of rice, maize and Arabidopsis using the Clustal Omega (https://www.ebi.ac.uk/Tools/msa/clustalo/) multiple sequence alignment program. The alignment file was then submitted to iTOL (https://itol.embl.de/) (Letunic and Bork, 2024), which produced a colored circular phylogenetic tree. A rectangular phylogenetic tree of rice OsRR proteins was constructed using the maximum likelihood method in MEGA 7 software with 1,000 bootstrap replications (Sun et al., 2022).

Gene structure, conserved motifs and conserved domain prediction of OsRR genes

Gene Structure Display Server 2.0 (http://gsds.gao-lab.org/) (Hu et al., 2015) was used to predict the structure of OsRR genes. OsRR protein sequences were submitted to the MEME suite online tool (https://meme-suite.org/meme/tools/meme) (Bailey et al., 2015) to discover 20 motifs. The conserved domains of OsRR proteins were retrieved from NCBI’s Conserved Domain Database (CDD) (https://www.ncbi.nlm.nih.gov/cdd/) (Marchler-Bauer et al., 2015). Finally, TBtools (Chen et al., 2020) software was employed to draw the motif and domain composition from MEME.xml file output and CDD hit data, respectively.

Chromosome mapping and gene duplication events of OsRR genes

The location of OsRR genes on 12 rice chromosomes was drawn using the MapGene2Chrom web (http://mg2c.iask.in/mg2c_v2.0/) (Chao et al., 2015) server. The KaKs calculator (https://ngdc.cncb.ac.cn/biocode/tools/BT000001) (Zhang, 2022) was used to identify gene duplication events and to estimate divergence time. Nonsynonymous (Ka) and synonymous (Ks) substitutions were computed for this purpose. Divergence time was calculated using the formula: Time = Ks/2λ, whereλ = 6.5 x 10^-9 for rice (Rehman et al., 2022).

Promoter analysis

The 1.5-kb region immediately upstream of OsRR genes was analyzed to identify potential cis-acting elements using the PlantCARE (http://bioinformatics.psb.ugent.be/webtools/plantcare/html/) (Lescot et al., 2002) online server. Subsequently, TBtools software was used to visualize the elements.

GO analysis

GO analysis based on biological process, molecular function and cellular component categorization was done using ShinyGO 0.80 (http://bioinformatics.sdstate.edu/go/) (Ge et al., 2020) and PANTHER18.0 (https://pantherdb.org/) (Thomas et al., 2022) online tools. We set the top ten significantly enriched GO terms to be identified at a false discovery rate of < 0.05 during analysis. ShinyGO 0.80 was used to build interaction networks of enriched GO terms under each category.

miRNA prediction

To predict potential miRNAs targeting OsRR genes, we submitted the genomic sequences of OsRR genes to the psRNATarget server (https://www.zhaolab.org/psRNATarget/) (Dai et al., 2018) and set the E-value cut-off to 4.0 because a lower E-value indicates a higher degree of similarity between small RNAs and their target genes (Singh and Mukhopadhyay, 2021). To display the results, we created an interactive network using Cytoscape software (Shannon et al., 2003).

Prediction of protein–protein interactions

STRING V12.0 (https://string-db.org/) (Szklarczyk et al., 2023) was utilized to construct a protein–protein interaction network among two-component system members. We submitted the protein sequences of six OsHKs, five OsHPs and 36 OsRRs to the server and configured the evidence-based interaction network with the setting ‘network type’ in the full STRING network, ‘minimum required interaction score’ at a high confidence score (0.7), and ‘active interaction sources’ at Textmining-on/ Experiments-on/ Databases-on/ Co-expression-on/ Neighborhood-on/ Gene Fusion-on/ Co-occurrence-on, leaving other settings at default.

Transcript profiles of OsRR genes for different abiotic stresses

RNA-seq and microarray data from the NCBI GEO dataset (https://www.ncbi.nlm.nih.gov/gds) were used to generate expression profiles of OsRR genes under different abiotic stress conditions, including salinity, drought, alkaline and wounding stress. The GEO accession numbers of salinity, drought, alkaline and wounding stress experiments are GSE133480, GSE128495, GSE104928 and GSE77097, respectively. Gene expression values in salinity, drought and alkaline stress experiments were calculated as transcripts per kilobase million (TPM), fragments per kilobase per million reads (FPKM) and reads per kilobase per million mapped reads (RPKM), respectively. These values were further log2-transformed to generate heatmaps. For the wounding stress experiment, log10-normalized expression values of OsRR genes were retrieved from microarray data in the NCBI database. Heatmaps were generated using TBtools, following the Euclidean distance method for clustering.

Experimental conditions of these RNA-seq or microarray samples are as follows. In the salinity stress experiment, two tissues (root and shoot) were sampled from two contrasting genotypes (CSR28, tolerant and IR28, susceptible) and subjected to two treatments (control and 150 mM salinity) at two sampling time points (6 and 54 h). NAC family genes, encoding NAC (NAM, ATAF1/2 and CUC) transcription factors and SNAC1 (STRESS-RESPONSIVE NAC 1), are mainly induced in guard cells under drought stress. SNAC1-overexpressing (SNAC1-OE) rice showed salinity and drought tolerance in field trials (Hu et al., 2006). In the drought stress experiment, SNAC1-OE (tolerant) and ZH11 wild type (ZH11-WT; susceptible) rice seedlings at the four-leaf stage were grown under normal and moderate drought stress conditions (Li et al., 2019). SNAC1-OE lines showed a significantly higher survival rate after recovery from drought stress compared to ZH11-WT. In the alkaline stress experiment, RNA-seq data were compared between an alkaline-tolerant (WD20342) and an alkaline-sensitive (Caidao) cultivar under control and alkaline stress conditions (Li et al., 2018). To induce wounding stress, wild-type Nipponbare rice leaves were mechanically injured and sampled after 30 min, 1 h, 2 h and 4 h.

Plant material and stress treatments

Nipponbare rice seeds were sterilized with 2% sodium hypochlorite solution, placed on a sterile wetted filter paper in a petri dish, and allowed to germinate in a growth chamber for four days, 13 h light/11 h dark, at 28 °C, 50 μmol m⁻² s⁻¹ photosynthetic photon flux density (PPFD), with no humidity control. Germinated seeds were transplanted into magenta boxes filled with sterile synthetic cultivation soil and grown under a cycle of 13 h light at 27 °C and 11 h dark at 23 °C , 350 μmol m⁻² s⁻¹ PPFD and 60% relative humidity. On the 11^th day, healthy and uniform-looking seedlings were transferred into Yoshida solution (Yoshida et al., 1976) at pH 5.0. On the 14^th day, the seedlings were subjected to 250 mM NaCl (salinity stress), 25% PEG (osmotic stress), 0.5% Na₂CO₃ (alkaline stress) or mechanical wounding for 0, 1, 6 and 12 h with three biological replicates. Mechanical wounding was inflicted following the procedure described by Reymond et al. (2000). In brief, each leaf blade was randomly given eight cuts on both sides with a pair of scissors and also pricked six times at different points randomly using a fine needle.

qPCR analysis

Total RNA was extracted using a QIAGEN RNeasy Plant Mini Kit according to the manufacturer’s instructions. To remove genomic DNA and make first-strand cDNA, ReverTra Ace qPCR RT Master Mix (TOYOBO) was used. RT-PCR was carried out in a StepOne RT-PCR (Applied Biosystems) machine using PowerTrack SYBR Green Master Mix (Applied Biosystems). Each reaction was performed in a final volume of 10 μl, containing 4 μl of cDNA sample, 5 μl of SYBR Premix and 0.5 μl of 10 μM primers. The thermal profile was set as 95 °C for 20 s, followed by 40 cycles of 95 °C for 3 s and 60 °C for 30 s. The internal reference was rice Ubiquitin 5. The relative expression level was computed using the 2^–∆∆CT method (Livak and Schmittgen, 2001). Gene-specific primer details are shown in Table 5.

Table 5. Primers used for qPCR

Gene	Sequence
OsRR4_qF	GTATGATCTGCTCAAGAGGGTG
OsRR4_qR	AGCTTGACAGGTTTCAGGAAG
OsRR6_qF	GTCCCCAACGTCAACATGATC
OsRR6_qR	CACGTTCTCCGACGACATGAT
OsRR9_qF	CTCCTTGTAGTCTCTTCTGTTCTTTCAG
OsRR9/10_qR	CATGGAACGGAGCCTCTATAGCC
OsRR10_qF	GCTTTCTCCTTGTAGCCTCTTCC
OsRR22_qF	GCATGCTTGTCTCCCCTTCT
OsRR22_qR	AGAATTGCTTGCACCTCCGA
OsPRR73_qF	CTGGAGGTGGCAATGGAAGT
OsPRR73_qR	GGCGTGTAACTGAAAACAGGC
Ubiq5 qF	ACCACTTCGACCGCCACTACT
Ubiq5 qR	ACGCCTAAGCCTGCTGGTT

AUTHOR CONTRIBUTIONS

K. I. conceptualized and designed the research. S. R. S. and S. M. S. I. conducted the experiments. S. R. S. analyzed the data and wrote the manuscript. K. I. supervised the experiments, provided the required facilities and made the final draft of the manuscript.

FUNDING

This research was supported by a Research Grant from the Koizumi Foundation to S. R. S. It was also partially supported by POISE project no. 16817624 from the Japan Science and Technology Agency to K. I., and BRIDGE Fellowship program ID BR 0317838 from the Japan Society for the Promotion of Science to S. M. S. I.

CONFLICTS OF INTEREST

The authors declare no conflicts of interest.

ACKNOWLEDGMENTS

The authors express thanks to CCRF, Niigata University for providing research facilities. The authors are sincerely thankful to Dr. Toshiaki Mitsui and Dr. Murat Aycan for providing laboratory facilities and technical assistance.

REFERENCES

• Almagro Armenteros, J. J., Salvatore, M., Emanuelsson, O., Winther, O., von Heijne, G., Elofsson, A., and Nielsen, H. (2019) Detecting sequence signals in targeting peptides. Life Sci. Alliance 2, e201900429. DOI: 10.26508/lsa.201900429
• Bailey, T. L., Johnson, J., Grant, C. E., and Noble, W. S. (2015) The MEME suite. Nucleic Acids Res. 43, W39–W49. DOI: 10.1093/nar/gkv416
• Berardini, T. Z., Reiser, L., Li, D., Mezheritsky, Y., Muller, R., Strait, E., and Huala, E. (2015) The Arabidopsis information resource: making and mining the “gold standard” annotated reference plant genome. Genesis 53, 474–485. DOI: 10.1002/dvg.22877
• Bhaskar, A., Paul, L. K., Sharma, E., Jha, S., Jain, M., and Khurana, J. P. (2021) OsRR6, a type-A response regulator in rice, mediates cytokinin, light and stress responses when over-expressed in Arabidopsis. Plant Physiol. Biochem. 161, 98–112. DOI: 10.1016/j.plaphy.2021.01.047
• Bielach, A., Hrtyan, M., and Tognetti, V. B. (2017) Plants under stress: involvement of auxin and cytokinin. Int. J. Mol. Sci. 18, 1427. DOI: 10.3390/ijms18071427
• Chao, J.-T., Kong, Y.-Z., Wang, Q., Sun, Y.-H., Gong, D.-P., Lv, J., and Liu, G.-S. (2015) MapGene2Chrom, a tool to draw gene physical map based on Perl and SVG languages. Yi Chuan 37, 91–97 (in Chinese). DOI: 10.16288/j.yczz.2015.01.013
• Chen, C., Chen, H., Zhang, Y., Thomas, H. R., Frank, M. H., He, Y., and Xia, R. (2020) TBtools: an integrative toolkit developed for interactive analyses of big biological data. Mol. Plant 13, 1194–1202. DOI: 10.1016/j.molp.2020.06.009
• Chou, K.-C., and Shen, H.-B. (2010) Plant-mPLoc: a top-down strategy to augment the power for predicting plant protein subcellular localization. PLoS One 5, 11335. DOI: 10.1371/journal.pone.0011335
• Chu, Z. X., Ma, Q., Lin, Y. X., Tang, X. L., Zhou, Y. Q., Zhu, S. W., Fan, J., and Cheng, B. J. (2011) Genome-wide identification, classification, and analysis of two-component signal system genes in maize. Genet. Mol. Res. 10, 3316–3330. DOI: 10.4238/2011.December.8.3
• D’Agostino, I. B., Deruère, J., and Kieber, J. J. (2000) Characterization of the response of the Arabidopsis response regulator gene family to cytokinin. Plant Physiol. 124, 1706–1717. DOI: 10.1104/pp.124.4.1706
• D’Agostino, I. B., and Kieber, J. J. (1999) Phosphorelay signal transduction: the emerging family of plant response regulators. Trends Biochem. Sci. 24, 452–456. DOI: 10.1016/s0968-0004(99)01465-6
• Dai, X., Zhuang, Z., and Zhao, P. X. (2018) psRNATarget: a plant small RNA target analysis server (2017 release). Nucleic Acids Res. 46(W1), W49–W54. DOI: 10.1093/nar/gky316
• Duvaud, S., Gabella, C., Lisacek, F., Stockinger, H., Ioannidis, V., and Durinx, C. (2021) Expasy, the Swiss bioinformatics resource portal, as designed by its users. Nucleic Acids Res. 49(W1), W216–W227. DOI: 10.1093/nar/gkab225
• Fang, Y., Xie, K., and Xiong, L. (2014) Conserved miR164-targeted NAC genes negatively regulate drought resistance in rice. J. Exp. Bot. 65, 2119–2135. DOI: 10.1093/jxb/eru072
• Gamage, D. G., Gunaratne, A., Periyannan, G. R., and Russell, T. G. (2019) Applicability of instability index for in vitro protein stability prediction. Protein Pept. Lett. 26, 339–347. DOI: 10.2174/0929866526666190228144219
• Ge, S. X., Jung, D., and Yao, R. (2020) ShinyGO: a graphical gene-set enrichment tool for animals and plants. Bioinformatics 36, 2628–2629. DOI: 10.1093/bioinformatics/btz931
• Geng, X., Zhang, C., Wei, L., Lin, K., and Xu, Z.-F. (2022) Genome-wide identification and expression analysis of cytokinin response regulator (RR) genes in the woody plant Jatropha curcas and functional analysis of JcRR12 in Arabidopsis. Int. J. Mol. Sci. 23, 11388. DOI: 10.3390/ijms231911388
• Goodstein, D. M., Shu, S., Howson, R., Neupane, R., Hayes, R. D., Fazo, J., Mitros, T., Dirks, W., Hellsten, U., Putnam, N., et al. (2012) Phytozome: a comparative platform for green plant genomics. Nucleic Acids Res. 40, D1178–D1186. DOI: 10.1093/nar/gkr944
• Grefen, C., and Harter, K. (2004) Plant two-component systems: principles, functions, complexity and cross talk. Planta 219, 733–742. DOI: 10.1007/s00425-004-1316-4
• Guo, Y., Yang, Y., Huang, Y., and Shen, H.-B. (2020) Discovering nuclear targeting signal sequence through protein language learning and multivariate analysis. Anal Biochem. 591, 113565. DOI: 10.1016/j.ab.2019.113565
• Heyl, A., and Schmulling, T. (2003) Cytokinin signal perception and transduction. Curr. Opin. Plant Biol. 6, 480–488. DOI: 10.1016/s1369-5266(03)00087-6
• Hirose, N., Makita, N., Kojima, M., Kamada-Nobusada, T., and Sakakibara, H. (2007) Overexpression of a type-A response regulator alters rice morphology and cytokinin metabolism. Plant Cell Physiol. 48, 523–539. DOI: 10.1093/pcp/pcm022
• Hosoda, K., Imamura, A., Katoh, E., Hatta, T., Tachiki, M., Yamada, H., Mizuno, T., and Yamazaki, T. (2002) Molecular structure of the GARP family of plant Myb-related DNA binding motifs of the Arabidopsis response regulators. Plant Cell 14, 2015–2029. DOI: 10.1105/tpc.002733
• Hu, B., Jin, J., Guo, A.-Y., Zhang, H., Luo, J., and Gao, G. (2015) GSDS 2.0: an upgraded gene feature visualization server. Bioinformatics 31, 1296–1297. DOI: 10.1093/bioinformatics/btu817
• Hu, H., Dai, M., Yao, J., Xiao, B., Li, X., Zhang, Q., and Xiong, L. (2006) Overexpressing a NAM, ATAF, and CUC (NAC) transcription factor enhances drought resistance and salt tolerance in rice. Proc. Natl. Acad. Sci. USA 103, 12987–12992. DOI: 10.1073/pnas.0604882103
• Hu, J., Zhang, H., and Ding, Y. (2013) Identification of conserved microRNAs and their targets in the model legume Lotus japonicus. J. Biotechnol. 164, 520–524. DOI: 10.1016/j.jbiotec.2013.01.026
• Huo, R., Liu, Z., Yu, X., and Li, Z. (2020) The interaction network and signaling specificity of two-component system in Arabidopsis. Int. J. Mol. Sci. 21, 4898. DOI: 10.3390/ijms21144898
• Hutchison, C. E., and Kieber, J. J. (2002) Cytokinin signaling in Arabidopsis. Plant Cell 14(suppl1), S47–S59. DOI: 10.1105/tpc.010444
• Hwang, I., and Sheen, J. (2001) Two-component circuitry in Arabidopsis cytokinin signal transduction. Nature 413, 383–389. DOI: 10.1038/35096500
• Imamura, A., Hanaki, N., Nakamura, A., Suzuki, T., Taniguchi, M., Kiba, T., Ueguchi, C., Sugiyama, T., and Mizuno, T. (1999) Compilation and characterization of Arabidopsis thaliana response regulators implicated in His-Asp phosphorelay signal transduction. Plant Cell Physiol. 40, 733–742. DOI: 10.1093/oxfordjournals.pcp.a029600
• International Rice Genome Sequencing Project, and Sasaki, T. (2005) The map-based sequence of the rice genome. Nature 436, 793–800. DOI: 10.1038/nature03895
• Jameson, P. E., and Song, J. (2016) Cytokinin: a key driver of seed yield. J. Exp. Bot. 67, 593–606. DOI: 10.1093/jxb/erv461
• Jiang, J.-J., Li, N., Chen, W.-J., Wang, Y., Rong, H., Xie, T., and Wang, Y.-P. (2022) Genome-wide analysis of the type-B authentic response regulator gene family in Brassica napus. Genes 13, 1449. DOI: 10.3390/genes13081449
• Kieber, J. J., and Schaller, G. E. (2014) Cytokinins. Arabidopsis Book 12, e0168. DOI: 10.1199/tab.0168
• Kumar, D., Ramkumar, M. K., Dutta, B., Kumar, A., Pandey, R., Jain, P. K., Gaikwad, K., Mishra, D. C., Chaturvedi, K. K., Rai, A., et al. (2023) Integration of miRNA dynamics and drought tolerant QTLs in rice reveals the role of miR2919 in drought stress response. BMC Genomics 24, 526. DOI: 10.1186/s12864-023-09609-6
• Lescot, M., Déhais, P., Thijs, G., Marchal, K., Moreau, Y., Van de Peer, Y., Rouzé, P., and Rombauts, S. (2002) PlantCARE, a database of plant cis-acting regulatory elements and a portal to tools for in silico analysis of promoter sequences. Nucleic Acids Res. 30, 325–327. DOI: 10.1093/nar/30.1.325
• Letunic, I., and Bork, P. (2024) Interactive Tree of Life (iTOL) v6: recent updates to the phylogenetic tree display and annotation tool. Nucleic Acids Res. 52(W1), W78–W82. DOI: 10.1093/nar/gkae268
• Li, L., Zhang, X., Ding, F., Hou, J., Wang, J., Luo, R., Mao, W., Li, X., Zhu, H., Yang, L., et al. (2024) Genome-wide identification of the melon (Cucumis melo L.) response regulator gene family and functional analysis of CmRR6 and CmPRR3 in response to cold stress. J. Plant Physiol. 292, 154160. DOI: 10.1016/j.jplph.2023.154160
• Li, N., Liu, H., Sun, J., Zheng, H., Wang, J., Yang, L., Zhao, H., and Zou, D. (2018) Transcriptome analysis of two contrasting rice cultivars during alkaline stress. Sci. Rep. 8, 9586. DOI: 10.1038/s41598-018-27940-x
• Li, X., Chang, Y., Ma, S., Shen, J., Hu, H., and Xiong, L. (2019) Genome-wide identification of SNAC1-targeted genes involved in drought response in rice. Front. Plant Sci. 10, 982. DOI: 10.3389/fpls.2019.00982
• Li, X.-T., Feng, X.-Y., Zeng, Z., Liu, Y., and Shao, Z.-Q. (2021) Comparative analysis of HSF genes from Secale cereale and its Triticeae relatives reveal ancient and recent gene expansions. Front. Genet. 12, 801218. DOI: 10.3389/fgene.2021.801218
• Livak, K. J., and Schmittgen, T. D. (2001) Analysis of relative gene expression data using real-time quantitative PCR and the 2^−ΔΔCT method. Methods 25, 402–408. DOI: 10.1006/meth.2001.1262
• Macovei, A., and Tuteja, N. (2012) microRNAs targeting DEAD-box helicases are involved in salinity stress response in rice (Oryza sativa L.). BMC Plant Biol. 12, 183. DOI: 10.1186/1471-2229-12-183
• Mangrauthia, S. K., Bhogireddy, S., Agarwal, S., Prasanth, V. V., Voleti, S. R., Neelamraju, S., and Subrahmanyam, D. (2017) Genome-wide changes in microRNA expression during short and prolonged heat stress and recovery in contrasting rice cultivars. J. Exp. Bot. 68, 2399–2412. DOI: 10.1093/jxb/erx111
• Marchler-Bauer, A., Derbyshire, M. K., Gonzales, N. R., Lu, S., Chitsaz, F., Geer, L. Y., Geer, R. C., He, J., Gwadz, M., Hurwitz, D. I., et al. (2015) CDD: NCBI’s conserved domain database. Nucleic Acids Res. 43, D222–D226. DOI: 10.1093/nar/gku1221
• Mok, D. W., and Mok, M. C. (2001) Cytokinin metabolism and action. Annu. Rev. Plant Physiol. Plant Mol. Biol. 52, 89–118. DOI: 10.1146/annurev.arplant.52.1.89
• Müller, B., and Sheen, J. (2008) Cytokinin and auxin interaction in root stem-cell specification during early embryogenesis. Nature 453, 1094–1097. DOI: 10.1038/nature06943
• Nagasaki, H., Itoh, J.-i., Hayashi, K., Hibara, K.-i., Satoh-Nagasawa, N., Nosaka, M., Mukouhata, M., Ashikari, M., Kitano, H., Matsuoka, M., et al. (2007) The small interfering RNA production pathway is required for shoot meristem initiation in rice. Proc. Natl. Acad. Sci. USA 104, 14867–14871. DOI: 10.1073/pnas.0704339104
• Nakamichi, N., Ito, S., Oyama, T., Yamashino, T., Kondo, T., and Mizuno, T. (2004) Characterization of plant circadian rhythms by employing Arabidopsis cultured cells with bioluminescence reporters. Plant Cell Physiol. 45, 57–67. DOI: 10.1093/pcp/pch003
• Osakabe, Y., Miyata, S., Urao, T., Seki, M., Shinozaki, K., and Yamaguchi-Shinozaki, K. (2002) Overexpression of Arabidopsis response regulators, ARR4/ATRR1/IBC7 and ARR8/ATRR3, alters cytokinin responses differentially in the shoot and in callus formation. Biochem. Biophys. Res. Commun. 293, 806–815. DOI: 10.1016/S0006-291X(02)00286-3
• Pareek, A., Singh, A., Kumar, M., Kushwaha, H. R., Lynn, A. M., and Singla-Pareek, S. L. (2006) Whole-genome analysis of Oryza sativa reveals similar architecture of two-component signaling machinery with Arabidopsis. Plant Physiol. 142, 380–397. DOI: 10.1104/pp.106.086371
• Peng, Y., Zhang, X., Liu, Y., and Chen, X. (2020) Exploring heat-response mechanisms of microRNAs based on microarray data of rice post-meiosis panicle. Int. J. Genomics 2020, 7582612. DOI: 10.1155/2020/7582612
• Qiang, Y., He, X., Li, Z., Li, S., Zhang, J., Liu, T., Tursunniyaz, M., Wang, X., Liu, Z., and Fang, L. (2023) Genome-wide identification and expression analysis of the response regulator gene family in alfalfa (Medicago sativa L.) reveals their multifarious roles in stress response. Front. Plant Sci. 14, 1149880. DOI: 10.3389/fpls.2023.1149880
• Ramírez-Carvajal, G. A., Morse, A. M., and Davis, J. M. (2008) Transcript profiles of the cytokinin response regulator gene family in Populus imply diverse roles in plant development. New Phytol. 177, 77–89. DOI: 10.1111/j.1469-8137.2007.02240.x
• Rehman, O. U., Uzair, M., Chao, H., Fiaz, S., Khan, M. R., and Chen, M. (2022) Role of the type-B authentic response regulator gene family in fragrant rice under alkaline salt stress. Physiol. Plant. 174, e13696. DOI: 10.1111/ppl.13696
• Ren, B., Liang, Y., Deng, Y., Chen, Q., Zhang, J., Yang, X., and Zuo, J. (2009) Genome-wide comparative analysis of type-A Arabidopsis response regulator genes by overexpression studies reveals their diverse roles and regulatory mechanisms in cytokinin signaling. Cell Res. 19, 1178–1190. DOI: 10.1038/cr.2009.88
• Reymond, P., Weber, H., Damond, M., and Farmer, E. E. (2000) Differential gene expression in response to mechanical wounding and insect feeding in Arabidopsis. Plant Cell 12, 707–719. DOI: 10.1105/tpc.12.5.707
• Roef, L., and Van Onckelen, H. (2010) Cytokinin regulation of the cell division cycle. In Plant Hormones. (ed.: Davies, P. J.), pp. 241–261. Springer, Dordrecht. DOI: 10.1007/978-1-4020-2686-7_12
• Sakai, H., Honma, T., Aoyama, T., Sato, S., Kato, T., Tabata, S., and Oka, A. (2001) ARR1, a transcription factor for genes immediately responsive to cytokinins. Science 294, 1519–1521. DOI: 10.1126/science.1065201
• Sakakibara, H. (2006) Cytokinins: activity, biosynthesis, and translocation. Annu. Rev. Plant Biol. 57, 431–449. DOI: 10.1146/annurev.arplant.57.032905.105231
• Schaller, G. E., Doi, K., Hwang, I., Kieber, J. J., Khurana, J. P., Kurata, N., Mizuno, T., Pareek, A., Shiu, S.-H., Wu, P., et al. (2007) Nomenclature for two-component signaling elements of rice. Plant Physiol. 143, 555–557. DOI: 10.1104/pp.106.093666
• Schaller, G. E., Kieber, J. J., and Shiu, S.-H. (2008) Two-component signaling elements and histidyl-aspartyl phosphorelays. Arabidopsis Book 2008, e0112. DOI: 10.1199/tab.0112
• Shannon, P., Markiel, A., Ozier, O., Baliga, N. S., Wang, J. T., Ramage, D., Amin, N., Schwikowski, B., and Ideker, T. (2003) Cytoscape: a software environment for integrated models of biomolecular interaction networks. Genome Res. 13, 2498–2504. DOI: 10.1101/gr.1239303
• Sharan, A., Soni, P., Nongpiur, R. C., Singla-Pareek, S. L., and Pareek, A. (2017) Mapping the ‘Two-component system’ network in rice. Sci. Rep. 7, 9287. DOI: 10.1038/s41598-017-08076-w
• Shi, F., Wang, M., and An, Y. (2020) Overexpression of a B-type cytokinin response regulator (OsORR2) reduces plant height in rice. Plant Signal. Behav. 15, 1780405. DOI: 10.1080/15592324.2020.1780405
• Singh, P., and Mukhopadhyay, K. (2021) Comprehensive molecular dissection of TIFY Transcription factors reveal their dynamic responses to biotic and abiotic stress in wheat (Triticum aestivum L.). Sci. Rep. 11, 9739. DOI: 10.1038/s41598-021-87722-w
• Stock, A. M., Robinson, V. L., and Goudreau, P. N. (2000) Two-component signal transduction. Annu. Rev. Biochem. 69, 183–215. DOI: 10.1146/annurev.biochem.69.1.183
• Sun, L., Nasrullah, Ke, F., Nie, Z., Xu, J., Huang, X., Sun, J., and Wang, P. (2022) Genome-wide identification and transcript analysis during fruit ripening of ACS gene family in sweet orange (Citrus sinensis). Sci. Hortic. 294, 110786. DOI: 10.1016/j.scienta.2021.110786
• Szklarczyk, D., Kirsch, R., Koutrouli, M., Nastou, K., Mehryary, F., Hachilif, R., Gable, A. L., Fang, T., Doncheva, N. T., Pyysalo, S., et al. (2023) The STRING database in 2023: protein–protein association networks and functional enrichment analyses for any sequenced genome of interest. Nucleic Acids Res. 51, D638–D646. DOI: 10.1093/nar/gkac1000
• Takagi, H., Tamiru, M., Abe, A., Yoshida, K., Uemura, A., Yaegashi, H., Obara, T., Oikawa, K., Utsushi, H., Kanzaki, E., et al. (2015) MutMap accelerates breeding of a salt-tolerant rice cultivar. Nat. Biotechnol. 33, 445–449. DOI: 10.1038/nbt.3188
• Tanaka, M., Takei, K., Kojima, M., Sakakibara, H., and Mori, H. (2006) Auxin controls local cytokinin biosynthesis in the nodal stem in apical dominance. Plant J. 45, 1028–1036. DOI: 10.1111/j.1365-313X.2006.02656.x
• Thomas, P. D., Ebert, D., Muruganujan, A., Mushayahama, T., Albou, L.-P., and Mi, H. (2022) PANTHER: making genome-scale phylogenetics accessible to all. Protein Sci. 31, 8–22. DOI: 10.1002/pro.4218
• Thumuluri, V., Almagro Armenteros, J. J., Johansen, A. R., Nielsen, H., and Winther, O. (2022) DeepLoc 2.0: multi-label subcellular localization prediction using protein language models. Nucleic Acids Res. 50, W228–W234. DOI: 10.1093/nar/gkac278
• To, J. P. C., Deruere, J., Maxwell, B. B., Morris, V. F., Hutchison, C. E., Ferreira, F. J., Schaller, G. E., and Kieber, J. J. (2007) Cytokinin regulates type-A Arabidopsis response regulator activity and protein stability via two-component phosphorelay. Plant Cell 19, 3901–3914. DOI: 10.1105/tpc.107.052662
• Wang, C., Wang, L., Liu, Q., Zhang, Y., and Dong, K. (2022) Genome-wide identification and characterization of PRR gene family and their diurnal rhythmic expression profile in maize. Int. J. Genomics 2022, 6941607. DOI: 10.1155/2022/6941607
• Wang, J., Du, Z., Huo, X., Zhou, J., Chen, Y., Zhang, J., Pan, A., Wang, X., Wang, F., and Zhang, J. (2020) Genome-wide analysis of PRR gene family uncovers their roles in circadian rhythmic changes and response to drought stress in Gossypium hirsutum L. PeerJ 8, e9936. DOI: 10.7717/peerj.9936
• Wang, W.-C., Lin, T.-C., Kieber, J., and Tsai, Y.-C. (2019) Response regulators 9 and 10 negatively regulate salinity tolerance in rice. Plant Cell Physiol. 60, 2549–2563. DOI: 10.1093/pcp/pcz149
• Woodhouse, M. R., Cannon, E. K., Portwood, J. L., 2nd, Harper, L. C., Gardiner, J. M., Schaeffer, M. L., and Andorf, C. M. (2021) A pan-genomic approach to genome databases using maize as a model system. BMC Plant Biol. 21, 385. DOI: 10.1186/s12870-021-03173-5
• Xie, K., Shen, J., Hou, X., Yao, J., Li, X., Xiao, J., and Xiong, L. (2012) Gradual increase of miR156 regulates temporal expression changes of numerous genes during leaf development in rice. Plant Physiol. 158, 1382–1394. DOI: 10.1104/pp.111.190488
• Yang, Z.-B., Liu, G., Liu, J., Zhang, B., Meng, W., Müller, B., Hayashi, K.-i., Zhang, X., Zhao, Z., De Smet, I., et al. (2017) Synergistic action of auxin and cytokinin mediates aluminum-induced root growth inhibition in Arabidopsis. EMBO Rep. 18, 1213–1230. DOI: 10.15252/embr.201643806
• Yoshida, S., Forno, D.A., Cock, J. H., and Gomez, K. A. (1976) Routine procedure for growing rice plants in culture solution. In Laboratory Manual for Physiological Studies of Rice. (eds: Yoshida, S., Forno, D. A., Cock, J. H., and Gomez, K. A.), pp. 61–66. International Rice Research Institute, Los Baños, Philippines.
• Zeng, J., Zhu, X., Haider, M. S., Wang, X., Zhang, C., and Wang, C. (2017) Genome-wide identification and analysis of the type-B authentic response regulator gene family in peach (Prunus persica). Cytogenet. Genome Res. 151, 41–49. DOI: 10.1159/000458170
• Zeng, R., Li, Z., Shi, Y., Fu, D., Yin, P., Cheng, J., Jiang, C., and Yang, S. (2021) Natural variation in a type-A response regulator confers maize chilling tolerance. Nat. Commun. 12, 4713. DOI: 10.1038/s41467-021-25001-y
• Zhang, A., Liu, Y., Wang, F., Li, T., Chen, Z., Kong, D., Bi, J., Zhang, F., Luo, X., Wang, J., et al. (2019) Enhanced rice salinity tolerance via CRISPR/Cas9-763 targeted mutagenesis of the OsRR22 gene. Mol. Breed. 39, 47. DOI: 10.1007/s11032-019-0954-y
• Zhang, M., Wang, F., Wang, X., Feng, J., Yi, Q., Zhu, S., and Zhao, X. (2022) Mining key genes related to root morphogenesis through genome-wide identification and expression analysis of RR gene family in citrus. Front. Plant Sci. 13, 1068961. DOI: 10.3389/fpls.2022.1068961
• Zhang, Z. (2022) KaKs_Calculator 3.0: calculating selective pressure on coding and non-coding sequences. Genomics Proteomics Bioinformatics 20, 536–540. DOI: 10.1016/j.gpb.2021.12.002
• Zhu, Q.-H., Upadhyaya, N. M., Gubler, F., and Helliwell, C. A. (2009) Over-expression of miR172 causes loss of spikelet determinacy and floral organ abnormalities in rice (Oryza sativa). BMC Plant Biol. 9, 149. DOI: 10.1186/1471-2229-9-149
• Zwack, P. J., and Rashotte, A. M. (2015) Interactions between cytokinin signalling and abiotic stress responses. J. Exp. Bot. 66, 4863–4871. DOI: 10.1093/jxb/erv172