Transcriptionally and post-transcriptionally regulated microRNAs in heat stress response in barley

Summary Selected barley miRNAs and their targets are regulated upon heat stress. Splicing of introns carrying miRNAs was induced by heat and correlated with the accumulation of mature miRNA.


Introduction
Plants as sessile organisms have developed exceptional abilities to cope with multiple environmental stimuli. The exposure to various stressful environmental conditions requires the induction of proper molecular responses operating at the transcriptional, post-transcriptional, translational, and post-translational levels (Banival et al., 2004;Sunkar, 2010). The importance of microRNAs (miRNAs) as key regulators of gene expression at the post-transcriptional level has been widely documented (Mallory and Vaucheret, 2006;Sunkar, 2010;Bielewicz et al., 2013). First identified in Caenorhabditis stem . The miRNA/miRNA* duplex is diced out from the precursor by a DCL1 enzyme (Kurihara and Watanabe, 2004), which, accompanied mainly by a SERRATE (SE) and a HYPONASTIC LEAVES 1 (HYL1), forms a microprocessing complex (Kurihara et al., 2006;Yang et al., 2006) that ensures efficiency and accuracy of precursor processing to mature miRNA (Dong et al., 2008). Subsequently, the miRNA/miRNA* duplex is exported to the cytoplasm in a HASTY (HST)-dependent manner (Park et al., 2005), where the miRNA guide strand is selectively loaded into an ARGONAUTE1 (AGO1)-containing RNAinduced silencing complex (RISC) (Vaucheret et al., 2004;Baumberger and Baulcombe, 2005). MiRNAs regulate gene expression in a sequence-specific manner (Baulcombe, 2004;Voinnet, 2009) by targeting mRNAs for cleavage or translational repression (Bartel, 2004;Brodersen et al., 2008;Beauclair et al., 2010).
MiRNAs have two major physiological roles in plants: (i) inducing cell differentiation in response to an endogenous stimulus during development and (ii) activation of an adaptive response to a particular exogenous stress. A number of miRNAs have been demonstrated to function in biotic and abiotic stress responses in plants (Sunkar and Zhu, 2004;Khraiwesh et al., 2012;Kruszka et al., 2012;Bielewicz et al., 2013). The role of miRNAs in plants infected by pathogenic bacteria, viruses, nematodes, and fungi has been widely reported (Ruiz-Ferrer and Voinnet, 2009;Katiyar-Agarwal and Jin, 2010). Plant miRNAs have also been reported to play an important role in response to cold (Sunkar and Zhu, 2004;Zhang et al., 2009;Lv et al., 2010), drought (Bartels and Sunkar, 2005;Liu et al., 2008;Jian et al., 2010;Pieczynski et al., 2013), salt (Liu et al., 2008;Zhao et al., 2009), or oxidative stresses induced by heavy metals, salinity, and nutrient deficiency (Sunkar et al., 2006;Jagadeeswaran et al., 2009).
Heat stress is one of the major abiotic stresses that can induce severe plant damage (Baniwal et al., 2004;Wahid et al., 2007). Given the increasing amount of evidence for worldwide climate change, the possible consequences of a temperature increase are a concern not only for climatologists but also for farmers and crop breeders. Heat stress is very often accompanied by drought, and these two combined stressful conditions are the major limitations to food production worldwide, especially in areas that use rainfed agriculture. Observed in recent decades, global warming is predicted to have a general negative effect on plant growth, development, and reproduction mostly due to the destructive effect of high temperature (Bita and Gerats, 2013). The involvement of miRNAs in response to heat has been previously reported in plants such as wheat (Xin et al., 2010), Brassica rapa , Arabidopsis thaliana (Guan et al., 2013), and rice (Mittal et al., 2012). Despite the importance of miRNAs, little is known about the mode of regulation of miRNA expression. The importance of the intron within the pri-miRNA 163 and its splicing affecting miRNA biogenesis in A. thaliana was shown (Bielewicz et al., 2013;Schwab et al., 2013). No data are available for miRNA biogenesis regulation in monocotyledonous plants.
Barley (Hordeum vulgare) grain production ranks fourth in global cereal production (FAO, http://www.fao.org), indicating the economic importance of this monocotyledonous crop plant. For spring barley, an air temperature ~20-25 ºC is considered to be optimum for growth and development. An increase in temperature above these values often results in the reduction of agricultural productivity (Acevedo et al., 2002;Hossain et al., 2012). Hence the importance of determining and understanding the biological factors and processes involved in heat response in plants.
It was considered to be of interest to find out whether barley miRNAs might function in response to high temperature and what might be the mechanisms regulating the miRNA biogenesis upon heat stress. Here it is shown that four mature miRNAs-miR160a, 166a, 167h, and 5175aare up-regulated under heat stress in barley. The studies also revealed that the level of their pri-miRNAs was affected in the heat stress conditions. Surprisingly, the splicing of the intron-containing precursors pri-miR160a and pri-miR5175a was induced by heat and correlated with the accumulation of mature miRNAs, suggesting the post-trancriptional regulation of miRNA precursor processing. Furthermore, conserved (HD-Zip transcription factors and auxin response factors) as well as novel target genes (HOX9, ACC oxidase, and Nek5-like kinase) of the heat-responsive barley miRNAs were experimentally identified using the 5′ rapid amplification of cDNA ends (RACE) approach and the degradome data. The level of target genes is regulated through a slicing mechanism, and it was experimentally confirmed that the level of the putative target mRNAs was down-regulated in heat stress conditions. A miRNA-target gene network in response to heat stress in barley is discussed. The results of the present study shed new light on the regulation of miRNA precursors and mature miRNA biogenesis, and indicate the important role of miRNAs in response to heat stress in barley.

Plant material and stress conditions
Spring barley plants, cultivar Rolap (Devaux et al., 1992), were grown in a Conviron environmental chamber (Conviron, Winnipeg, Manitoba, Canada) with a 16 h day/8 h night photoperiod, and 800 μmol light conditions. Plants were grown at 22 ºC day/15 ºC night temperature in 250 ml pots containing field soil, and were watered to maintain optimal growth conditions of 70% SWC (soil water content). Two-week-old plants at the three-leaf stage (code 13 of Zadoks system; Zadoks et al., 1974) were subjected to heat treatment at 35.5 ºC. Control plants were grown at a constant temperature of 22 °C. Plants were collected after 3, 6, and 24 h of heat stress in three biological replicates. For developmental analyses of expression of selected miRNAs, the whole plants from five growth stages in three biological replicates for each growth stage were used as previously described (Kruszka et al., 2013).

RNA isolation and northern blots
A 100 mg aliquot of tissue was used for isolation of total RNA enriched in small RNA using a protocol as previously described Kruszka et al., 2013). The quality and quantity of RNA were measured with a NanoDrop ND-1000 spectrophotometer, and RNA integrity was estimated on agarose gels. RNA electrophoresis, blot transfer, and hybridization were performed as previously reported (Kruszka et al., 2013). DNA oligo probes (Sigma) were 5′ labelled with [γ-32 P]ATP (6000 Ci mmol -1 ; Hartmann Analytic GmbH, Germany). A DNA probe complementary to U6 small nuclear RNA (snRNA) was used, and the U6 hybridization signal was taken as a loading control. The blots were exposed for 1 week to a phosphorimaging screen (Fujifilm) and scanned with a Fujifilm FLA5100 reader (Fujifilm Co., Ltd, Tokyo, Japan). Blots were quantified with Multi Gauge V2.2 software.

Quantitative real-time PCR (RT-qPCR)
cDNA was prepared and RT-qPCR amplification was performed as previously reported (Kruszka et al., 2013). Expression levels were calculated with the relative quantification method (2 -ΔCt ) as fold change value and presented in the form log 10 2 -ΔCt . The R 2 values of analysed data (≥0.997) were calculated with LinRegPCR software (Ramakers et al., 2003). Since the pri-miRNA expression levels are lower than that of the reference gene, the expression profiles are shown in a positive data range without changing the actual values by shifting the zero value of the graph's y-axis to the basal expression level of the whole experiment . Primers designed and used for the validation of the levels of splicing isoforms were complementary to the exon-intron or exon-exon junctions. Primers are listed in Supplementary Table S1 at JXB online.

Target mRNA prediction and cleavage site analysis
The target identification procedure was based on three publicly available programs: miRanda, TAPIR, and psRNAtarget (John et al., 2004;Bonnet et al., 2010;Dai and Zhao, 2011). Full-length cDNA sequences from H. vulgare were downloaded from the NCBI and used as the input. Results of the prediction software were subsequently scored based on mismatches present in the predicted sRNA:target duplex alignments where: (i) G:U substitution was counted as 0.5 point; and (ii) other substitutions were given 1 point according to Jones- Rhoades and Bartel (2004) and Schwab et al. (2005). No insertions/deletions were allowed. Duplexes having a total mismatch score ≥4 and substitution score in positions 1:12 of sRNA >2.5 points were discarded. For further analysis, only duplexes predicted by at least two programs with MFE energy equal to or lower than -15 kcal mol -1 were selected. Protein domain prediction analysis using the HMMER3 tool was performed on putative targets. Protein sequences corresponding to identified cDNAs were extracted using the BLASTx algorithm on the H. vulgare protein data set from the NCBI (Altschul et al., 1997;Eddy, 2011). Proteins corresponding to targets with similarity equal to or lower than E-value 0.001 with at least 50% coverage of the known protein were selected. Finally, all putative targets were ranked based on the lowest MFE energy of the duplex and selected according to the position on the list.
To examine cleavage sites of predicted target mRNAs, 5′RACE experiments were conducted with the SMARTer RACE cDNA Amplification Kit (Clontech, Mountain View, CA, USA) according to the manufacturer's protocol. PCR products were cloned into the pGEM T-Easy vector (Promega) and sequenced (Faculty's Laboratory of Molecular Biology Techniques, Adam Mickiewicz University in Poznan, Poland).

Construction and analysis of degradome libraries
A library comprised of 26-27 nt cDNA tags derived from the 5′ ends of polyadenylated RNAs was prepared from poly(A)-enriched RNA from 47-and 68-day-old barley plants. Library construction was similar to previously described methods (Addo-Quaye et al., 2009;German et al., 2009) except that the RNA and DNA adaptors were modified for compatibility to the TruSeq sequencing system of Illumina. The degradome preparation and data analysis will be described in more detail in another publication. The degradome libraries were sequenced by Fasteris SA (Switzerland).

Mature microRNAs and their precursors are induced by heat stress
To determine the structures of miRNA genes and their transcripts, eight barley cDNA nucleotide sequences for microRNA precursors (160a, 166a, 167h, 530-5p, 1120b, 1432-5p, 5175a, and 5203) deposited in GenBank (http://www. ncbi.nlm.nih.gov/) (Benson et al., 2012) (see the Materials and methods) were selected for which hairpin structures, carrying conserved or newly estimated miRNA homologues, were computationally predicted. First the organization of all barley MIR genes was determined. The length and structure for each barley MIR gene determined, the nucleotide sequence, and the position of its miRNA/miRNA* within the gene, as well as plant orthologues of barley genes are presented in Table 1. The schematic organization of barley MIR genes as well as the miRNA-containing hairpin structures of the barley pre-miRNAs studied, and the orthologues of barley genes in rice, wheat, or Brachypodium distachyon are presented in Supplementary Fig. S1 at JXB online. Moreover, the expression profiles at the levels of pri-miRNAs and mature miRNAs in five developmental stages are also shown ( Supplementary  Fig. S2).
The regulation of expression of selected heat-responsive barley miRNAs was studied here. To address this issue, the structure of miRNA genes has to be known. In total the structures of only 16 barley MIR genes are available (Kruszka et al., 2013;this study). The accumulation of the mature miRNAs encoded by these genes was examined by northern hybridization ( Fig. 1; Supplementary Fig. 3 at JXB online). An increase in the mature miRNA level was observed for four mature miRNAs: miR166, 167, 160, and 5175 (Fig. 1). For miR166, a single band of 21 nt in length was detected. Its accumulation was 1.6-fold higher in plants subjected to 6 h and 24 h heat stress as compared with control plants (Fig. 1A). Northern hybridization with a specific probe for miR167 revealed the presence of two mature miRNAs, 21 nt and 22 nt long, which both were induced in heat-stressed plants (Fig. 1B). The expression level of both miR167 molecules was slightly increased in the case of 3 h and 6 h heat-treated plants (1.2-fold higher at both time points than in control plants), and was substantially elevated up to 1.6-fold in plants exposed to longer heat treatment. In the case of the mature miR160, its level was increased 1.4-and 1.8-fold in 6 h and 24 h heat-treated plants, respectively (Fig. 1C). Similarly to the above-described barley miRNAs, the expression level of miR5175 was up-regulated during heat treatment (Fig. 1D). The level of mature miR5175 was particularly elevated in plants exposed to high temperature for a longer period of time (at least 24 h). Using the same hybridization probe that was complementary to miR5175, a 24 nt long molecule was also detected. However, the level of accumulation of this molecule did not change in plants subjected to heat stress. The observed accumulation of mature miRNAs under heat stress raised a question regarding the simultaneous alterations in the level of their primary transcripts as well as putative splicing isoforms. The level of barley pri-miRNAs in heat stress conditions was tested using semi-quantitative RT-PCR. In the case of four miRNA precursors, pri-miR166a, 167h, 160a, and 5175a, considerable changes were also detected in their expression patterns in heat-treated plants, as compared with control plants (Fig. 2, left panels). The results obtained by RT-PCR were confirmed using RT-qPCR (Fig. 2, right  panels).
The intronless pri-miR166a accumulated during heat stress and its level increased gradually up to 13-fold after 24 h stress as compared with control plants ( Fig. 2A, right panel). For the other intronless miR167h primary transcript, the major accumulation was observed in plants stressed for 24 h, and its level was 2.6-fold higher than in control conditions (Fig. 2B, right panel). The accumulation of pri-miR166a and pri-miR167h which is accompanied by an increase in their mature miRNAs upon heat stress suggests the transcriptional regulation of the miRNA level.
Interestingly, the high temperature influenced not only the level of the primary transcripts of barley miRNAs but also the level of their spliced isoforms. miR160a is encoded within the second intron of a non-coding transcript containing three introns altogether ( Supplementary Fig. S1 at JXB online). Pri-miR160a unspliced transcript (isoform A) accumulated considerably in heat-stressed plants compared with control plants (Fig. 2C, left and middle panels). Additionally, two splicing isoforms were observed in heat-treated plants which were not detected in plants grown in control conditions. In the first splicing isoform (isoform B), the last intron was removed, while in the third isoform (isoform C) the last intron as well as the middle, miRNA-containing intron, were spliced. The splicing isoform C, in which intron 2, containing miR160a, was already removed, revealed the most upregulation in stressed plants, as shown by RT-PCR (Fig. 2C, right panel). Using an RT-qPCR pair of primers specific for the exon 2-exon 3 junction, and exon 3, strong accumulation of this splicing isoform C was confirmed; its level was 8-to 34-fold higher as compared with control plants. Taken together, the mature miR160 up-regulation and accumulation of the unspliced transcript as well as its spliced isoforms in response to heat stress indicate the transcriptional and post-transcriptional regulation of the mature miR160 level. miR5175a is encoded within the 10th intron of the RNA polymerase II phosphatase-like/MIR5175a gene. A similar effect of high temperature on the splicing efficiency of an miRNA-carrying intron was observed for the pri-miR5175a transcript (Fig. 2D, left panel). The level of miR5175a precursor containing intron 10 comprising a miRNA/miRNA* (isoform A) was slightly decreased in plants after 3 h and 6 h of heat stress, while in plants after 24 h of stress its level was similar to that of control plants. This observation was confirmed by RT-qPCR using primers specific for intron 10 in which miR5175a is located (Fig. 2D, middle panel). Simultaneously, a substantial accumulation of the splicing isoform lacking intron 10 (isoform B) in stress-treated plants was observed. RT-qPCR amplification with a specific pair of primers anchored in the exon 10-exon 11 junction and exon 11 confirmed strong accumulation of the spliced isoform whose level was 4-fold higher as compared with control plants (Fig. 2D, right panel). These results implicate the post-transcriptional regulation of miR5175a expression.
The described results show an unusual observation that heat-induced splicing and splicing efficiency of the Fig. 1. The level of mature barley miRNA166, 160, 167, and 5175 is affected by heat. The mature miRNAs miR166, 167, 160, and 5175 (A-D, respectively) were detected by northern hybridization in control and heat conditions. The level of mature miRNAs was analysed at different time points: 3, 6, and 24 h. U6 was used as a loading control. The level of miRNA in heat stress was quantified relative to that in control conditions at the respective time points tested. miRNA-bearing introns studied correlate with mature miRNA accumulation.

Elevated expression of barley heat-regulated miRNAs correlates with down-regulation of the expression of their target genes
To understand the function of barley miRNAs in response to heat stress, the level of putative miRNA-regulated target mRNAs was investigated. The expression level of target mRNAs was estimated by RT-qPCR. The 5′RACE experiments as well as the degradome sequencing were carried out to verify the predicted cleavage sites in target mRNAs.
In A. thaliana, miR166a is known to target mRNAs coding for HD-Zip transcription factors including PHAVOLUTA (PHV), which regulates auxiliary meristem initiation and leaf development . A barley orthologue of the Arabidopsis PHV gene (GenBank accession no. AK364215.1) was analysed as a potential target candidate. The expression level of PHV mRNA decreased nearly 2-fold in heat-stressed plants at all time points tested as compared with control plants. Using 5′RACE and degradome data, the cleavage site in the PHV transcript guided by miR166a was confirmed (Figs 3A, 4A). The REVOLUTA gene (REV; GenBank accession no. AK362009) is the other member of the HD-Zip transcription factor family predicted to be a Fig. 2. Barley pri-miRNA precursors are induced by heat stress. The panels on the left represent RT-PCR analysis of the expression of pri-miR166a, 167h, 160a, and 5175a (A-D, respectively) in control and heat stress conditions; primer positions are marked by black triangles on the pri-miRNA schematic structures; miRNA and miRNA* are marked in red and blue, respectively. Ubiquitin (UBQ) amplification was used as a loading control. The charts on the right show the real-time PCR measurements of the expression level of pri-miRNA precursors during heat stress. The expression studies were performed at different time points: 3, 6, and 24 h. The level of pri-miRNA in control conditions was assumed to be '1', and the levels of pri-miRNA during heat stress were quantified relative to this standard; bars represent the means of three independent biological samples ±SD. M, GeneRuler 100 bp Plus or 1 kb Plus DNA ladders.
putative target mRNA of barley miR166a. The expression level of the REV transcript was also reduced under heat stress, mostly after 24 h exposure to heat. The degradome data confirmed the predicted slicing site (Figs 3B, 4B). A homeoboxleucine zipper protein HOX9-like gene (GenBank accession no. AK354023.1) was also identified as a novel miR166a target gene in barley. The HOX9 protein is another member of the HD-Zip transcription factor family which was reported to regulate rice embryogenesis (Nagasaki et al., 2007). The expression level of HOX9 mRNA was decreased especially after 6 h and 24 h of heat stress, and the predicted slicing site guided by miR166a was confirmed by the degradome data (Figs 3C, 4C). mRNAs encoding the auxin response transcription factor ARF17 (GenBank accession no. AK363500.1) and ARF13 (GenBank accession no. AK364144) were confirmed as conserved miR160 targets in barley. In Arabidopsis thaliana miR160 was reported to target ARF genes . The level of barley ARF17 and ARF13 transcripts was considerably down-regulated to a 4-fold lower level in heat stress conditions than in control plants, especially in plants exposed to elevated temperature for 24 h (Fig. 3D, E). The cleavage site guided by miR160a in ARF17 mRNA was verified by the 5′RACE approach. In addition, the predicted cleavage sites in both ARF17 and ARF13 transcripts were confirmed by the degradome data (Fig. 4D, E).
Two potential mRNA targets of miR167h were predictedan auxin response transcription factor ARF8 gene (GenBank accession no. AK368822.1) and a serine/threonine-protein kinase Nek5-like gene (GenBank accession no. AK363207) (Fig. 3F, G). In Arabidopsis and rice, the members of the auxin response transcription factor family were previously reported as miR167 family targets (Jones- Rhoades and Bartel, 2004). In the present studies it was shown that barley ARF8 mRNA is strongly down-regulated under heat stress ( Fig. 3F). However, it was not possible to map the miR167hguided cleavage site of the ARF8 transcript precisely as all examined 5′RACE clones were at least 35 nt shorter than the expected 3′ cleavage products. The lack of precisely sliced 3′ products was most probably due to rapid degradation of the cleaved fragments of the ARF8 transcript. However, it was determined that the strong peak from the degradome library sequencing data corresponded precisely to the predicted cleavage site (Fig. 4F). A Nek5-like gene transcript was identified as a novel target gene of barley miR167h (Fig. 3G). The mRNA level of the NEK5-like kinase gene decreased considerably in plants exposed to 24 h heat stress. The cleavage site in the NEK5 transcript that is recognized by miR167h was also verified using the 5′RACE approach; however, it was not possible to detect the slicing site using the degradome data.
In the case of miR5175a, a novel putative target mRNA was predicted-an enzyme involved in ethylene biosynthesis ACC-like oxidase (GenBank accession no. AK370060.1) (Raz and Ecker, 1999). The expression level of ACC-like oxidase mRNA was down-regulated mainly in 24 h heatstressed plants (Fig. 3H). However, it was not possible to detect slicing products either in 5′RACE or in the degradome sequencing data.
The results show the conservation of the target mRNAs between barley and dicot plants and reveal the existence of new additional, as yet unrecognized target mRNAs.

Discussion
Barley MIR genes show a diverse genomic organization, as previously reported (Kruszka et al., 2013). Most of the MIR genes in this study contain at least one intron. Interestingly, some members of the same family were encoded by genes that differ noticeably in structure, as was revealed for the barley MIR166 family. A similar observation has been described for members of the MIR160 family in Arabidopsis (Szarzynska et al., 2009). These observations may indicate that miRNAs encoded by the different genes representing the same MIR gene family might differ in their expression regulation pattern at the post-transcriptional level. The intron-containing pri-miRNA transcripts investigated herein undergo a complex splicing pattern. This phenomenon, previously described for both Hordeum and Arabidopsis miRNA precursors (Szarzynska et al., 2009;Kruszka et al., 2013), may suggest a potential regulatory role of splicing in mature miRNA biogenesis. This suggestion has been recently confirmed for Arabidopsis miRNA precursors (Bielewicz et al., 2013;Schwab et al., 2013). The presence of an intron in pri-miRNA transcripts can promote mature miRNA accumulation, through a mechanism that probably operates at the level of miRNA processing or stability.
Remarkably, not only were intron-containing pri-miRNAs undergoing splicing detected, but, more importantly, the splicing events were observed to be strongly affected in barley plants subjected to heat stress. The findings suggest the post-transcriptional regulation of the expression of these two barley miRNAs. Surprisingly, substantial accumulation of miR160a and miR5175a splicing isoforms was observed, mostly lacking introns which contain the miRNA-bearing hairpin structures. Previous studies on pre-mRNA splicing have indicated the shutdown of pre-mRNA splicing in response to heat stress (Yost and Lindquist, 1986;Bond, 1988;Shin et al., 2004). Moreover, in plants, it has been previously shown that insertion of an artificial secondary structure (stem-loop) in an intron strongly inhibits splicing in vivo in tobacco protoplasts (Goodall and Filipowicz, 1991). A low level of barley pri-miRNAs splicing in control plants was observed, which might suggest the inhibitory role of stemloop structures in intron splicing. However, the presence of such secondary structure might also modulate the splicing in response to different conditions such as heat stress. RNA structures are highly sensitive to the concentration and types of ions and osmolytes present (Draper, 2004;Lambert and Draper, 2007;Mullen et al., 2012). It is known that many abiotic stresses, such as high salinity, heat, and drought, impact cellular concentrations of ions and osmolytes in plants (Seki et al., 2007;Mullen et al., 2012;Chan et al., 2013), which might impact RNA secondary structure and regulate splicing and other aspects of post-transcriptional gene regulation. Hence, structural elements in pri-miRNAs might function as sensors of abiotic stresses, such as temperature, drought, and salinity (Reddy et al., 2013).
In the case of intronless MIR166a and MIR167 genes, an increase in the level of both mature miRNAs and their precursors was detected upon heat stress. These results strongly suggest the transcriptional regulation of the expression of these miRNAs; however, the effect of pri-miRNA transcript stability cannot be excluded.
The present study has shown that barley miRNA precursors as well as mature miRNA molecules were heat inducible, which might be required for barley thermotolerance through down-regulation of their target genes. The involvement of plant miRNAs in response to high temperature was previously shown in plants such as Arabidopsis Guan et al., 2013), Populus tomentosa (Chen et al., 2012), or Brassica rapa . Interestingly, the highthroughput sequencing of the wheat miRNAome allowed identifiction of a set of eight miRNAs that were up-regulated in heat stress conditions, including tae-miR160, 166, and 167 (Xin et al., 2010). This finding is consistent with the present results; however, miR160 is strongly down-regulated during heat stress in Populus (Chen et al., 2012). These results suggest the existence of a common set of miRNAs involved in response to heat in plants, though the response of miRNAs might differ significantly depending on the species studied. The transcriptome analysis of wheat revealed a large number of heat stress-responsive genes encoding transcription factors and proteins involved in phytohormone biosynthesis/signalling, calcium and sugar signalling pathways, RNA metabolism, ribosomal proteins, primary and secondary metabolism, as well as proteins related to other stresses (Qin et al., 2008). However, only a few of the identified heat stress-responsive genes were tested and confirmed as being regulated by miRNAs in heat stress conditions (Xin et al., 2010).
When the level of a miRNA is up-regulated, usually its target gene is down-regulated post-transcriptionally. In the present study, the conserved and the novel targets of barley miRNAs that function in diverse pathways and processes were identified. Barley as a crop plant is sensitive to high temperature, especially during the grain-filling stage (Hakala et al., 2012). Frequently identified features of plants exposed to high temperatures are a significant decrease in height and biomass reduction. In addition, heat stress which occurs during the early stages of barley development is known to reduce tillering, thereby decreasing crop yield (Savin et al., 1997;Wallwork et al., 1998;Alam et al., 2007). Considering barley growth and development, heat stress was reported to decrease the days to visible awns, days to heading, and days to ripe harvest (Wallwork et al., 1998). Heat stress not only affects the aerial parts of the plants but also influences the growth and morphology of the roots, causing the arrest of root elongation and branching (Marschner, 1995).
The module of the putative network of heat-responsive barley miRNAs and their respective target genes involved in the heat stress response is presented in Fig. 5. PHAVOLUTA, REVOLUTA, and HOX9 are the HD-Zip transcription factors whose expression is regulated by miR166a, an evolutionarily conserved plant miRNA (Barik et al., 2014). These transcription factors have been shown to be required for the establishment of the apical meristem, proper pattern formation in lateral organs, and leaf development Nagasaki et al., 2007). In plants such as maize and Arabidopsis, the adaxial/abaxial leaf polarity is established by fine-tuned interplay between an abaxial gradient of miR166, adaxial presence of a ta-siRNA (trans-acting small interfering RNA), and their targets (Juarez et al., 2004;Nogueira et al., 2007;Husbands et al., 2009). Thus, the changes in the level of miR166 in heat stress might influence the leaf morphology. Shoot morphology and regeneration in vitro are also affected by changes in the miR160a level (Qiao et al., 2012). It was determined in this study that an increase in barley miR160a during heat stress considerably down-regulates the expression level of the auxin response transcription factors ARF17 and ARF13. ARF17 is a regulator of GH3-like early auxin response genes, and its function was reported in emerging leaf symmetry, premature inflorescence development, and flower development (Mallory et al., 2005). Auxin is an important hormone involved in many plant processes including cell elongation and division, and thus growth (Overvoorde et al., 2010;Gallavotti, 2013;Sauer et al., 2013;de Wit et al., 2014). As the heat stress effects on plants often include growth inhibition, the auxin signalling pathway might be an important component of plant thermotolerance.
It was also confirmed that another auxin-responsive transcription factor gene, ARF8, a target gene of miR167h, was down-regulated by heat in barley. In Arabidopsis, ARF8 regulates floral organ development and vascular patterning of petals (Tabata et al., 2010). Moreover, ARF8 together with AFR6 and ARF17 regulates adventitious root formation in Arabidopsis (Gutierrez et al., 2009(Gutierrez et al., , 2012. Furthermore, the involvement of miR160 and miR167 families targeting ARF transcription factors in lateral and adventitious root formation, as well as root cap development, has been shown in rice (Meng et al., 2009(Meng et al., , 2010. The three ARF genes display overlapping expression domains, interact genetically, and regulate each other's expression at both the transcriptional and post-transcriptional level by modulating miR160 and miR167 availability (Gutierrez et al., 2009). In addition, the fine-tuning of adventitious root initiation requires the GH3 genes that modulate jasmonic acid homeostasis (Gutierrez et al., 2012). As root growth inhibition is one of the major consequences of the effect of high temperature on plants, it is of great importance to understand this cross-talk between plant hormones, transcription factors, and miRNAs.
Interestingly, a novel target gene of barley miR167h, Nek5-like kinase, was also identified in the present study. In Arabidopsis, the NEK protein subfamily, including NEK4, NEK5, and NEK6, plays an important role in epidermal cell expansion and microtubule functions (Motose et al., 2011). It has been reported by Nakajima and Suzuki (2013) that heat stress induces misfolded protein aggregation in plant cells in a microtubule-dependent manner. In addition, a link between miRNA activities and microtubules was also observed in plants and animals, suggesting that miRNA-based translational repression or miRNA activities in general may rely on cytoskeleton dynamics in diverse species (Ma et al., 2013).
In conclusion, the present study revealed that four barley miRNAs-miR160a, 166a, 167h, and 5175a-respond to heat stress. The expression level of barley pri-miRNAs as well as the splicing of miRNA-bearing introns was induced in plants subjected to heat stress, suggesting the transcriptional as well post-transcriptional regulation of miRNA expression. The accumulation of miRNAs triggered the changes in the expression level of conserved and novel target genes characterized in this study. The findings showed that barley miRNAs, together with target genes whose expression they regulate, might function in a complex regulatory network developed by barley to cope with stressful conditions. The findings will contribute to explaining the mechanisms of plant thermotolerance.

Supplementary data
Supplementary data are available at JXB online. Figure S1. Schematic representation of barley MIR genes and the hairpin structures of barley/plant orthologues of the pre-miRNA precursors. Figure S2. The expression profiles at the level of pri-miR-NAs and mature miRNAs during barley development.   Figure S4. Expression of the heat stress marker gene HSP17 in barley plants subjected to heat stress. Table S1. Primers and hybridization probes used in the experiments.