Distribution of BTB/POZ genes in O. rufipogon genome
The HMM profile and BlastP outputs yielded a total of 128 preliminary hits, which were confirmed by Pfam searches to verify the occurrence of bona fide BTB/POZ domains. Of the total preliminary hits, 110 genes were determined to contain bona fide BTB/POZ domains, while the rest were potential random homologies. Because the HMM profile search often misidentifies truncated domains, the 110 putative BTB/POZ protein-encoding gene loci were used as queries for another iteration of BlastP search across the O. rufipogon protein database. This extra iteration did not reveal additional gene copies, thus 110 putative BTB/POZ protein-encoding genes of O. rufipogon with genomic loci identifier (OrBTB1 to OrBTB110) comprised the final gene set from all chromosomes (Supplementary Table S2). In silico analysis of the physico-chemical attributes of putative BTB/POZ proteins revealed diverse properties across the gene family with molecular weights ranging from 18.49 kDa to 281.12 kDa and an average of 58.76 kDa (Supplementary Table S2), and isoelectric points (pI) ranging from 4.27 to 10.55 and average of 6.19.
BTB/POZ gene loci were distributed across all 12 chromosomes of O. rufipogon (Fig. 1A). Chromosomes-08 and -10 contained the largest number of copies with 24 and 23 loci, respectively, while chromosome-12 contained the lowest number with only three loci. Evident in the chromosomes with the largest number of loci, BTB/POZ genes occurred in small to large tandem arrays, indicating their origins through tandem duplication. For instance, a large cluster of 23 loci occurred in a region of chromosome-10 with coordinates of 13.12Mb and 13.74Mb (total cluster size of 0.62Mb). Smaller tandem arrays were located on chromosomes-11 and -02, with seven members (coordinates = 24.04Mb to 26.62Mb), and four members (coordinates = 11.25Mb to 11.33Mb), respectively. Two small tandem arrays were located on chromosome-08, with eight (coordinates = 6.79 Mb to 6.95Mb) and six (coordinates = 6.79Mb to 6.95Mb) members, respectively. Prediction of subcellular localization showed that most OrBTB proteins were localized in the cytoplasm (30 genes), chloroplast (29 genes), nucleus (23 genes), plasma membrane (17 genes), mitochondria (8 genes), and extracellular space (3 genes) (Supplementary Table S2; Supplementary Fig. S1).
Phylogeny of O. rufipogon BTB/POZ genes
The unrooted tree constructed by the Neighbor-Joining method showed that the 110 BTB/POZ-encoding genes of O. rufipogon formed three phylogenetically distinct sub-groups (Fig. 2). The largest group (clade-2) was comprised of 49 members, of which the majority (88% = 43 genes) had only the BTB/POZ domain, while a smaller minority of its members contained other additional domains such as MATH (5 genes) and PA/peptidase (1 gene). Eight of the BTB/POZ genes in clade-2 were part of a small tandem array on chromosome-08 (OrBTB8.6/ORUFI08G07460, OrBTB8.7/ORUFI08G07490, OrBTB8.8/ORUFI08G07520, OrBTB8.9/
ORUFI08G07540, OrBTB8.10/ ORUFI08G07550, OrBTB8.11/ORUFI08G07560, OrBTB8.12/ORUFI08G07610, OrBTB8.13/ORUFI08G07650) along the 6.7Mb to 6.9Mb genomic coordinates. As the general structures of these genes were not significantly similar to each other, they appeared to have been created through tandem duplication events followed by structural and functional diversification. Another tandem array on chromosome-08 comprised of OrBTB8.19/ ORUFI08G22900, OrBTB8.20/ORUFI08G22910, OrBTB8.21/ORUFI08G22920, OrBTB8.22/ORUFI08G22930, OrBTB8.23/ORUFI08G23440 and OrBTB8.24/ORUFI08G23450 in genome coordinates 23.2Mb to 23.6Mb clustered in the same sub-group in clade-2.
Clade-1 was comprised of 35 members, of which 21 also contained the MATH domain. Surprisingly, out of the 35 members in this clade, 22 were located in a tandem array on chromosome-10 (13.1Mb to 13.7 Mb), signifying that members of this sub-group were the outcomes of the largest duplication event in this protein family. The smallest clade (clade-3) consisted of 26 members, of which ten genes had an NPH3 domain (OrBTB 3.1/ ORUFI03G07610, OrBTB 3.2/ORUFI03G07690, OrBTB 3.3/ORUFI03G26720, OrBTB 3.4/ORUFI03G28320, OrBTB4.6/
ORUFI04G28350, OrBTB6.1/ ORUFI06G05610, OrBTB 8.4/ORUFI08G02180, OrBTB 9.4/ORUFI09G11290, OrBTB12.2/ORUFI12G19040, OrBTB12.2/
ORUFI12G19040), two had a BACK domain (OrBTB2.1/ORUFI02G11310, OrBTB6.4/ORUFI06G17340), two had an NPR1 domain (OrBTB1.1/
ORUFI01G35500, OrBTB3.5/ORUFI03G29880), and two had a TAZ domain (OrBTB1.2/ORUFI01G43500, OrBTB4.4/ORUFI04G17810). It appeared that this small clade was comprised of genes with much more diverse functions.
Structure of O. rufipogon BTB/POZ genes and their protein products
The majority (87 of 110) of the BTB/POZ proteins contained only one copy of the BTB/POZ domain, while the other members had multiple copies (2 to 7) (Fig. 3). The MATH domain was the most common combination with the BTB/POZ domain, with a total of 34 members (out of 110) having such a combination. Of these, 25 members had one MATH domain and three others had more than two copies of the MATH domain. Some BTB/POZ members combined with other functional domains such as Ank (13 members), NPH3 (10 members), DUF3420 (3 members), zf-TAZ (2 members), BACK (2 members), NPR1_like (2 members), Arm (1 member), F5_F8_ty (1 member), Methyltr (1 member), PA (1 member), and Peptidase (1 member). The combinatorial occurrence of these other domains suggests that the BTP/POZ gene family of O. rufipogon may have also undergone the process of domain shuffling.
Members of the BTB/POZ gene family of O. rufipogon were widely diverse in terms of intron-exon architectures, with exon numbers ranging from 1 to19 (Fig. 4). The majority of genes (70%, 77 out of 110 genes) had a relatively simple intron-exon architecture with one to four exons. Of these, 26 genes had a single exon. Complex intron-exon architectures (ten or more exons) were also evident across the gene family, although in a much smaller number of genes (10%; 11 out of 110). The most complex was evident in OrBTB5.2/ORUFI05G16760 with 19 exons.
Peptide motifs are highly conserved amino acid sequences with potential protein structural and functional implications. Analysis by MEME revealed ten significant peptide motifs among the BTB/POZ proteins, designated as motif-1 to motif-10 (Supplementary Fig. S2). Particularly noteworthy were motif-2 (ETFAAHRCVLAARSPVF) and motif-3 (IDDMEPAVFKALLHFIYTDSLP), occurring 84-times and 73-times, respectively. Motif-8 (YLRDDCFTIRCDVTVV) represented the least frequently occurring, only 57-times (Supplementary Fig. S3). Overall, trends in intron-exon architecture and conservation of peptide motifs suggest that apart from the important role of gene duplication, the evolution of the BTB/POZ protein family of O. rufipogon likely had significant contributions from domain-shuffling and exon-shuffling, which may have important implications to functional specialization.
Synteny of BTB/POZ tandem arrays across flowering plants
Evolutionary dynamics of the BTB/POZ genes of O. rufipogon were investigated to understand the role of gene duplication and selection on gene family expansion. Trends revealed by the analysis of amino acid sequence homology showed 24 BTB/POZ gene pairs representing co-paralogs (Supplementary Table S3). Of these, eight genes appeared to have originated from segmental duplication, which gave rise to four gene pairs (Fig. 1B). These results indicate that the continuous expansion of this large gene family in O. rufipogon was facilitated primarily by tandem duplication and secondarily by non-tandem duplication. It is noteworthy that the largest number of co-paralogs are located in the chromosome-10 (18 genes) and chromosome-8 (10 genes) clusters. Analysis of non-synonymous/synonymous ratios of the duplicated gene pairs indicated that these duplicated BTB/POZ genes have undergone purifying selection based on the Ka/Ks <1 (Hurst 2002). The average divergence time estimated for the duplicated genes was 54.2 MYA (Supplementary Table S3).
Synteny reflects the magnitude by which large gene families evolved independently in different lineages. Trends revealed from dual synteny analysis with four other reference genomes representing a wide diversity across flowering plants such as Oryza sativa ssp. japonica (domesticated from O. rufipogon; IRGSP-1.0), Arabidopsis thaliana (dicot lineage; TAIR10), Brachypodium distachyon (small-genome monocot distantly related to Oryza; Brachypodium_distachyon_v3.0), and Sorghum bicolor (large-genome diploid monocot as distant lineage to Oryza; Sorghum_bicolor_NCBIv3) indicated that colinear gene clusters indeed evolved from common ancestral gene copies (Fig. 5). As expected, the largest number of syntenic O. rufipogon BTB/POZ gene clusters were observed in its domesticated version Oryza sativa ssp. japonica with a total of 71 syntenic loci. Continuous breakage of syntenic blocks was evident across the comparative evolutionary panel with a significant decline in synteny with increasing phylogenetic distance, i.e., 56 conserved syntenic loci in B. distachyon and S. bicolor, and only 8 in the dicot lineage represented by A. thaliana (Supplementary Table S4).
It is noteworthy that all eight syntenic loci between O. rufipogon and A. thaliana were also conserved with the two monocot lineages (Sorghum, Brachypodium). Members of these syntenic loci across the monocot and dicot lineages may provide an important gene set for developing Conserved Ortholog Set (COS) markers in comparative plant genomics studies, and particularly for the mining of novel BTB/POZ alleles from O. rufipogon. Most importantly, dual synteny analysis revealed that five BTB/POZ orthologs of O. rufipogon are conserved in Brachypodium (OrBTB4.7/ORUFI04G29880, OrBTB8.8/ORUFI08G07520, OrBTB8.20/ ORUFI08G22910, OrBTB10.3/ORUFI10G11060, OrBTB11.6/ ORUFI11G21740) and six in Sorghum (OrBTB4.7/ORUFI04G29880, OrBTB6.6/ ORUFI06G25940, OrBTB8.8/ORUFI08G07520, OrBTB8.22/ORUFI08G22930, OrBTB10.2/ORUFI10G11050, OrBTB11.4/ORUFI11G21310), but not in its domesticated species. This loss of orthologs in O. sativa suggests putative gene copies that have been left behind in the wild gene pool during domestication.
Upstream sequences of O. rufipogon BTB/POZ genes
A total of 102 unique classes of sequence motifs were identified across the -2,000-bp region of the BTB/POZ genes (Fig. 6). On average, genes contained more than 50 classes of putative cis-elements. The potential functional diversity represented by these putative cis-elements (PlantCare database) suggests multiple pathways by which the BTB/POZ gene family members are regulated. Out of the 110 BTB/POZ family members in O. rufipogon, the OrBTB4.4/ORUFI04G17810 had the largest number of hits (n = 93) for known cis-elements, while the OrBTB11.7/
ORUFI11G21760 had the lowest (n = 29) (Supplementary Table S5). The vast majority of putative cis-elements appeared to be relevant to abiotic stress response mechanisms and hormonal signalling. This was further reiterated by the enrichment of GO keywords such as abiotic stress, jasmonic acid, abscisic acid, wound response, development, and combinatorial gene network regulation. Putative cis-elements of MYB-type transcription factors were the most significantly enriched (i.e., 484 MYB, 104 MYB-binding sites, 93 MYB-like sequences, 73 MYB recognition sites), suggesting key functions of MYB transcription factors in the cellular and biological processes by which BTB/POZ genes are either directly or indirectly involved (Ambawat et al. 2013; Park et al. 2010; Abe et al. 1997).
A strong association of BTB/POZ genes with general stress response pathways was implied by the enrichment of different classes of transcription factors involved in stress-mediated regulation, including 387 STRE (stress response element), 379 ABRE (ABA-responsive elements), 83 ABRE4, 83 ABRE3a, and 89 DRE (dehydration-responsive element). This enrichment of stress-responsive cis-elements is consistent with the reported functions of transcription factors in dehydration and osmotic ( Min et al. 2021; Muthurajan et al. 2021; De los Reyes et al. 2015; Yoshida et al. 2015; Hossain et al. 2010), oxidative (Yue et al. 2019), heat (Sharma et al. 2021), cold (Wang et al. 2021; Yun et al. 2010), and sugar starvation (Chen et al. 2019) stresses. These observations further support the potential contributions of BTB/POZ diversification to the inherent stress tolerance potential of the wild species O. rufipogon (Kitazumi et al. 2018; Atwell et al. 2014).
Another class of significantly enriched cis-elements among the BTB/POZ genes appeared to be targets of MYC transcription factors (399), particularly those involved in jasmonic acid (JA) signaling and wound-inducible responses during herbivory, pathogen invasion, and physical wounding. This trend is particularly evident from the enrichment of WOUND RESPONSIVE ELEMENTS (WRE), which occurred 133 times across the BTB/POZ gene family (Luo et al. 2019; Boter et al. 2004). The G-box, which is known to regulate senescence-induced gene expression, occurred in 440 upstream locations of BTB/POZ genes (Liu et al. 2016).
Novelty of O. rufipogon BTB/POZ genes based on exon-intron structure
We compared the intron-exon structures across orthologous O. rufipogon (OrBTB) and O. sativa ssp. japonica (OsBTB) genes that are located in three syntenic blocks located on chromosome-01 (four copies), chromosome-08 (four copies), and chromosome-10 (10 copies) in either one-to-one or one-to-few Clusters of Orthologous Groups (COGs) from available data in Ensembl Plants (https://plants.ensembl.org/index.html; Fig. 7A, Supplementary Table S6). In the COG from chromosome-01, we observed OrBTB1.3/ORUFI01G44480 (5 exons) and OrBTB1.4/ORUFI01G46370 (13 exons) as having the same exon-intron structures as their respective orthologous genes in O. sativa ssp. japonica, i.e., Os01t0908200 and Os01t0932600. In contrast, the respective orthologs of OrBTB1.2/
ORUFI01G43500 (7 exons) and OrBTB1.5/ORUFI01G47340 (2 exons) in O. sativa ssp. japonica, i.e., Os01t0893400 (5 exons) and Os01t0948900 (1 exon) appeared to have diversified through the loss of exons in the domesticated genome.
Similarly, in chromosome-08 COG, the O. sativa ssp. japonica gene Os08t0516200 had maintained the same intron-exon structure (1 exon) as its ortholog in O. rufipogon,OrBTB8.19/ORUFI08G22900 (1 exon). However, the three other O. sativa genes in this cluster, i.e.,Os08t0516500 (2 exons), Os08t0522700 (1 exon), Os08t0523400 (1 exon), appeared to have drastic exon reduction relative to their respective O. rufipogon orthologs which maintained 9 exons (OrBTB8.21/
ORUFI08G22920), 7 exons (OrBTB8.23/ORUFI08G23440), and 16 exons (OrBTB8.24/ORUFI08G23450).
We further observed a high degree of variability in intron-exon conservation in the chromosome-10 syntenic block between O. rufipogon and O. sativa ssp. japonica. Of the ten COGs in this block, only two orthologous pairs, i.e., O. sativa-Os10t0428900 with O. rufipogon-OrBTB10.8/ORUFI10G11150 (2 exons) and O. sativa-Os10t0434000 with O. rufipogon-OrBTB10.14 /ORUFI10G11410 (1 exon) had conserved exon-intron structures. Two other COGs in this cluster, i.e., OrBTB10.7/ORUFI10G11140 (1 exon) to Os10t0427300 (2 exons), and OrBTB10.18/ ORUFI10G11500 (1 exon) to Os10t0435400 (4 exons) appeared to have diversified through gain of exons in the domesticated genome. Six other COGs in this cluster showed indications of exon losses during domestication, i.e.,
OrBTB10.1/ORUFI10G11040 (2 exons) to Os10t0423300 (1 exon), OrBTB10.4/ ORUFI10G11080 (8 exons) to Os10t0425400 (1 exon), OrBTB10.11/
ORUFI10G11190 (3 exons) to Os10t0429200 (1 exon), OrBTB10.13/ ORUFI10G11230 (17 exons) to Os10t0430600 (7 exons), OrBTB10.20/ ORUFI10G11550 (4 exons) to Os10t0435900 (1 exon), and OrBTB10.22/ ORUFI10G11820 (4 exons) to Os10t0439333 (1 exon).
We also compared the gene body of the orthologous set of OrBTB genes that were lost in Oryza sativa ssp japonica but present in Brachypodium distachyon and Sorghum bicolor (Fig. 7B, C). The general trend indicated that both high and low degrees of conservation in BTB/POZ gene structures existed after the divergence of these three monocot species. For example, the respective orthologs of OrBTB10.3/ ORUFI10G11060 (1 exon) and OrBTB11.6/ORUFI11G21740 (2 exons) in B. distachyon (KQJ96779, KQJ87842) maintained a conserved exon-intron structure. Similarly, the respective orthologs of OrBTB4.7/ORUFI04G29880, OrBTB6.6/ ORUFI06G25940, OrBTB8.22/ORUFI08G22930 in S. bicolor (EES11589, KXG20566, KXG25577) all had maintained their conserved single exon structures after speciation. There were also examples of exon losses among BTB/POZ genes as a consequence of speciation of O. rufipogon, B. distachyon, and S. bicolor. For example, loss or gain of exons were evident between OrBTB4.7/ ORUFI04G29880 (1 exon), OrBTB8.8/ORUFI08G07520 (5 exons), and OrBTB8.20/ ORUFI08G22910 (2 exons) and their respective orthologs in S. bicolor with 2 (KQJ85084), 1 (KQJ95656), and 1 (PNT68446) exons, respectively. The respective orthologs of OrBTB8.8/ORUFI08G07520 (5 exons) and OrBTB10.2/ORUFI10G11050 (7 exons) had a higher number of exons than their Sorghum counterparts OQU80115 and EER94211, each with only a single exon.
These observations indicated a generally higher similarity of the COGs of O. rufipogon, B. distachyon, and S. bicolor amongst each other in terms of exon-intron structure compared to COGs in the direct descendant of O. rufipogon – O. sativa ssp. japonica. These observations further suggest that differences between O. rufipogon (progenitor) and O. sativa ssp. japonica were due to domestication. Many BTB/POZ genes/alleles found in the wild have been conserved during monocot speciation but eliminated or altered by more intense selection during the domestication of O. sativa ssp. japonica. These results further support the hypothesis that the progenitor species of the modern-day cultivated rice is a potential source of novel alleles for allelic enrichment of the domesticated germplasm.
Novelty of O. rufipogon BTB/POZ genes based on coding sequence variation
Orthologous BTB/POZ genes showed disparate dynamics in terms of mRNA and polypeptide length, with the rare exception of the perfectly homologous O. rufipogonOrBTB10.14/ORUFI10G11410 and O. sativa ssp. japonica Os10t0434000, which encodes a highly conserved CDS of 1113-bp and polypeptide product of 370 amino acid residues (Supplementary Table S6). For the majority of orthologous gene pairs across the syntenic blocks between O. rufipogon and O. sativa, amino acid sequence alignment showed a variable percentage of identity (pid). For example, in the chromosome-10 cluster, Os10t0435400 (49.6%) and Os10t0425400 (64.2%) showed some of the lowest pid under coverage percentage (cov) of 91.6% and 38.8%, respectively, compared to their O. rufipogon ortholog OrBTB10.18/
ORUFI10G11500 and OrBTB10.4/ORUFI10G11080 (Supplementary Fig. S4). In comparison with Brachypodium orthologs (5 orthologous loci), the pid ranged from 45.8 (KQJ96779 versus OrBTB10.3/ORUFI10G11060, with a cov of 89.9) to 24.4 (KQJ87842 versus OrBTB11.6/ORUFI11G21740, with a cov of 95.3). In comparison with Sorghum orthologs (6 orthologous loci), the pid ranged from 44.5 (EES11589 versus OrBTB4.7/ORUFI04G29880, with a cov of 86.8 to 27.2 (OQU80115 versus OrBTB8.8/ ORUFI08G07520 with a cov of 91.7) (Supplementary Fig. S5).
Upstream regulatory signatures of orthologous BTB/POZ genes
The regulatory sequences of orthologous O. rufipogon and O. sativa ssp. japonica BTB/POZ genes also exhibited significant variation. Trends in the COGs of a given syntenic block showed differential patterns of enrichment for various types of cis-elements. For example, in chromosome-01 COGs, the O. rufipogon orthologs had a total of 176 putative cis-elements while the O. sativa ssp. japonica orthologs had 216 putative cis-elements (Fig. 8A, Supplementary Table S7). The ten most predominant classes of cis-elements among O. sativa ssp. japonica orthologs include the MYB (n = 35), ABRE (n = 24), G-box (n = 23), ARE (n = 15), STRE (n = 15), MYC (n = 13), TGACG-motif (n = 10), as-1 (n = 10), MBS (n = 6), and DRE (n = 5). In O. rufipogon, the ten most predominant cis-elements were the MYB (n = 26), ABRE (n = 19), STRE (n = 17), G-box (n = 17), ARE (n = 12), MYC (n = 8), TGACG-motif (n = 1), as-1 (n = 6), MBS (n = 5), and WUN (n = 5) (Fig. 8A).
Similarly, dissection of COGs in chromosome-08 cluster showed a total of 191 putative cis-elements in O. sativa ssp. japonica and 163 in O. rufipogon. The ten most predominant classes in O. sativa ssp. japonica orthologs are MYB (n = 37), ABRE (n = 25), G-box (n = 17), TGACG-motif (n = 16), as-1 (n = 16), MYC (n = 12), ARE (n =7), MBS (n = 7), WRE3 (n = 6), and TCA (n = 6). Among O. rufipogon orthologs, the most predominant classes were MYB (n = 27), MYC (n= 17), as-1 (n = 13), TGACG (n = 13), WRE3 (n = 10), ABRE (n = 10), STRE (n = 8), G-box (n = 7), TCA (n = 5), and ARE (n = 6). COGs in chromosome-10 clusters showed generally similar trends of differential cis-element enrichment between the O. sativa ssp. japonica and O. rufipogon orthologs.
Inter-generic comparison (O. rufipogon vs. S. bicolor vs. B. distachyon) of upstream sequences of orthologous OrBTB genes lost in Oryza sativa ssp japonica relative to Brachypodium distachyon and Sorghum bicolor revealed that O. rufipogon orthologs have higher enrichment of cis-elements (n = 389) than S. bicolor (n = 285) and B. distachyon (n = 322) (Fig. 8B, Supplementary Table S7). This suggested that based on intra-generic history (i.e., domestication = O. rufipogon vs. O. sativa ssp. japonica) and inter-generic history (i.e., speciation = Oryza vs. Sorghum vs. Brachypodium), the regulatory mechanisms operating on orthologous BTB/POZ genes have diverged considerably due to speciation and domestication. Similar to the trends revealed by gene sequence and structure comparison, novel alleles that appeared to have distinct regulatory mechanisms are the potential features of O. rufipogon BTB/POZ genes relative to O. sativa ssp. japonica orthologs.
Potential role of miRNAs in the regulation of O. rufipogon BTB/POZ genes
Prediction of putative miRNAs with possible roles in regulating BTB/POZ genes suggested that 95 members of O. rufipogon BTB/POZ gene family are likely targets of 392 miRNAs belonging to more than 200 families (Supplementary Table S8). This initial set of candidate miRNAs was further reduced by false-positive filtration at a cut-off expectation score of 3.5, revealing 52 high-confidence candidate miRNA-target sites on BTB/POZ genes for 31 miRNA species (Supplementary Fig. S6). Of the 36 BTB/POZ genes that appeared to be regulated post-transcriptionally by miRNAs, four (OrBTB4.5/ORUFI04G27660, OrBTB4.7/ORUFI04G29880, OrBTB8.23/ORUFI08G23440, OrBTB10.11/ORUFI10G11190) had the highest number of potential miRNA target sites, with three sites each. Another set of eight BTB/POZ genes (OrBTB1.3/ORUFI01G44480, OrBTB3.7/ ORUFI03G38580, OrBTB4.1/ORUFI04G02660, OrBTB6.8/ORUFI06G26390, OrBTB7.6/ ORUFI07G26650, OrBTB8.13/ORUFI08G07650, OrBTB8.2/ORUFI08G02030, OrBTB8.21/ORUFI08G22920) had two putative target sites. The rest of the 24 BTB/POZ genes contained single putative miRNA-target sites.
Our analysis revealed that miR5075 (phytohormone synthesis), miR5809 (regulation of finger transcription factors and subtilisin), and miR2927 (plant growth) appeared to be the most prominent post-transcriptional regulators of BTB/POZ genes with eight, six, and four targets, respectively. We also observed that all mRNA cleavage targets of miRNA5075 in OrBTB1.1/ORUFI01G35500, OrBTB1.5/
ORUFI01G47340, OrBTB3.2/ORUFI03G07690, and OrBTB6.8/ORUFI06G26390 clustered in the same phylogenetic group (Fig. 2), consistent with the highly conserved nature of miRNA target sites.
Gene Ontology and protein-protein interaction analysis
Gene Ontology (GO) enrichment across the 110 members of the BTB/POZ gene family of O. rufipogon indicated eight high-level GO terms for biological process, i.e., protein binding, histone acetyltransferase activity, peptide-lysine-N-acetyltransferase activity, peptide N-acetyltransferase activity, transcription coregulator activity, N-acetyltransferase activity, N-acyltransferase activity, and acetyltransferase activity (Supplementary Fig. S7). Enrichment trends for molecular function indicated that 13 genes were distinctly involved in protein ubiquitination, protein modification by small protein conjugation, and protein modification by small protein conjugation or removal. Enrichment trends for cellular components revealed three BTB/POZ genes (OrBTB1.2/ORUFI01G43500,OrBTB2.6/ORUFI02G24050, OrBTB4.4/ORUFI04G17810) involved in mechanisms of host-pathogen interaction (Supplementary Fig. S7, Supplementary Table S9).
Prediction and modeling of protein-protein interaction among the BTB/POZ protein family members of O. rufipogon using the STRING database and Markov Cluster Algorithm (MCL) revealed three major interaction clusters (Fig. 9). The largest (Cluster-1) consisted of 87 genes (nodes) with 88 interactions (edges) (Supplementary Table S10). The average local clustering coefficient in this cluster was 0.985, with OrBTB 5.1/ORUFI05G13500 appearing to be the most likely central hub of the network. As the second-largest, Cluster-2 was comprised of three nodes (OrBTB4.4/ ORUFI04G17810, OrBTB11.2/ORUFI11G02560, OrBTB2.6/
ORUFI02G24050) and two edges, with an average local clustering coefficient of 0.667. The OrBTB11.2/ORUFI11G02560 connected the other two nodes of this cluster. Cluster-3 was the smallest, having only two nodes with an average local clustering coefficient of 1.0 (Supplementary Table S10).
Spatio-temporal expression of BTB/POZ genes
Publicly available transcriptome datasets of O. rufipogon included various tissue/organ (i.e., tiller bases, roots, leaf pulvini, leaf sheaths, nodes, culms, panicles >5 cm) and environmental response (i.e., Fe deficiency, salt stress, cold stress) RNA-Seq libraries. Mining these datasets to profile the range of spatio-temporal expression across the BTB/POZ gene family showed substantial expression differences among the members.
Hierarchical clustering of spatio-temporal expression showed that 40 of the 110 gene family members are expressed in all tissues/organs and 13 had non-detectable expression in all tissues/organs (Fig. 10). These observations implied varied functions of gene family members in cellular housekeeping processes as well tissue/organ-specific processes. For instance, 66 (tiller base), 57 (roots), 63 (leaf pulvini), 59 (leaf sheaths), 61 (nodes), 58 (culms), 62 (panicles, >5 cm), 63 (panicles, 1-5 cm), 86 (panicles, <1cm), and 53 (leaf blades) members of the gene family had significantly higher tissue/organ-limited expression.
The majority of BTB/POZ genes exhibited specific patterns of upregulation or downregulation in certain tissues/organs, with the exception of OrBTB1.2/
ORUFI01G43500, OrBTB8.3/ ORUFI08G02060, OrBTB8.24/ORUFI08G23450, OrBTB10.6/ORUFI10G11100, OrBTB4.3/ORUFI04G14290, and OrBTB10.2/
ORUFI10G11050. For instance, OrBTB1.2 was downregulated in the nodes, suggesting its potential role in vegetative organ development and growth.
Certain BTB/POZ gene family members were clearly regulated in response to environmental stressors. For example, under Fe (iron) deficiency, 51 gene family members had constitutive expression while 60 others were differentially expressed in the roots (Fig. 11A). In particular, OrBTB1.2/ORUFI01G43500 was upregulated by Fe deficiency while OrBTB3.7/ORUFI03G38580 and OrBTB 6.8/ORUFI06G26390 were downregulated. BTB/POZ genes also appeared to be involved in low-temperature response. In particular, OrBTB 1.3/ORUFI01G44480 and OrBTB 4.4/ORUFI04G17810 showed continuous upregulation under cold stress (Fig. 11B).
Salinity stress appeared to regulate the expression of 62 BTB/POZ genes in the leaves and 46 others in the roots (Fig. 11C). Notably, OrBTB10.22/
ORUFI10G11820, OrBTB8.13/ORUFI08G07650, OrBTB4.7/ORUFI04G29880, and OrBTB10.3/ORUFI10G11060 were downregulated in the leaves, while OrBTB8.23/ORUFI08G23440, OrBTB9.3/ ORUFI09G06920, OrBTB6.6/ ORUFI06G25940 were downregulated in the roots. These observations indicated that the functional specification of different BTB/POZ gene family members as defined by the integration of various intrinsic and extrinsic signals, consistent with the trends revealed by cis-element and protein-protein interaction analyses.
Direct involvement of BTB/POZ genes in various stress response mechanisms provided the impetus for further validation of expression changes under different abiotic stress regimes (Fig. 12 A, B, C) by qRT-PCR. The expression of nine representative genes, selected based on gene ontology, was profiled in O. rufipogon during exposure to nutrient and abiotic stress conditions. Results indicated that OrBTB3.5/ORUFI03G29880 was consistently downregulated with progressive exposure to heat stress, while other gene family members (OrBTB4.1/
ORUFI04G02660,OrBTB7.1/ORUFI07G00140, OrBTB8.1/ ORUFI08G00320, OrBTB11.2/ORUFI11G02560, OrBTB11.3/ORUFI11G12230) were consistently upregulated after prolonged exposure to heat stress (Fig. 12D).
Salinity stress rapidly but transiently induced several BTB/POZ genes, as evident from significant increases in transcript abundances as early as 24 hours after the onset of stress. However, elevated transcript levels were not sustained during prolonged exposure suggesting that these genes were responsive to osmotic shock but were not likely adaptive in nature (Fig. 12E). Additionally, a few other BTB/POZ genes, particularly OrBTB11.3/ORUFI11G12230, OrBTB3.3/ORUFI03G26720, OrBTB3.5/ORUFI03G29880, and OrBTB8.1/ORUFI08G00320 had relatively late but stable upregulation, suggesting potential roles in adaptive mechanisms. Under nitrogen deficiency, OrBTB3.3/ORUFI03G26720, OrBTB3.5/ORUFI03G29880, OrBTB8.1/ORUFI08G00320 and OrBTB4.1/ORUFI04G02660 were upregulated, and OrBTB3.3/ORUFI03G26720, OrBTB3.5/ORUFI03G29880, OrBTB8.1/
ORUFI08G00320 were downregulated. Iron (Fe) toxicity stress was a potent inducer of BTB/POZ genes as all nine members of the expression cohort showed a significant upregulation (Fig. 12F). This suggests that O. rufipogon has the potential for iron toxicity tolerance.