CRISPR / Cas9-Induced Mutagenesis of Semi-Rolled Leaf1,2 Confers Curled Leaf Phenotype and Drought Tolerance by Inﬂuencing Protein Expression Patterns and ROS Scavenging in Rice ( Oryza sativa L.)

: Rice leaf morphology is an essential agronomic trait to develop drought-tolerant genotypes for adequate and stable crop production in drought-prone areas. Here, rolled leaf mutant plants were acquired by CRISPR / Cas9-based mutagenesis of Semi-rolled leaf1,2 ( SRL1 and SRL2 ) genes, and isobaric tags for relative and absolute quantiﬁcation (iTRAQ) based proteomic analysis was performed to analyze the subsequent proteomic regulation events. Homozygous mutants exhibit decreased chlorophyll content, transpiration rate, stomatal conductance, vascular bundles (VB), stomatal number, and agronomic traits with increased panicle number and bulliform cells (BCs). Under drought stress, mutant plants displayed lower malondialdehyde (MDA) content while higher survival rate, abscisic acid (ABA) content, superoxide dismutase (SOD), catalase (CAT) activities, and grain ﬁlling percentage compare with their wild type (WT). Proteomic results revealed that 270 proteins were signiﬁcantly downregulated, and 107 proteins were upregulated in the mutant line compared with WT. Proteins related to lateral organ boundaries’ (LOB) domain (LBD) were downregulated, whereas abiotic stress-responsive proteins were upregulated in the CRISPR mutant. LBD proteins (Q5KQR7, Q6K713, Q7XGL4, Q8LQH4), probable indole-3-acetic acid-amido synthetase (Q60EJ6), putative auxin transporter-like protein 4 (Q53JG7), Monoculm 1 (Q84MM9) and AP2 (Apetala2) domain-containing protein (Q10A97) were found to be hub-proteins. The hybrids developed from mutant restorers showed a semi-rolled leaf phenotype with increased panicle number, grain number per panicle, and yield per plant. Our ﬁndings reveal the intrinsic value of genome editing and expand the knowledge about the network of proteins for leaf rolling and drought avoidance in rice. Our study provides some novel evidence supporting the regulatory module in leaf rolling and drought tolerance in rice. yield. These results might pave the way for understanding the gene action and inheritance pattern of CRISPR / Cas9 generated mutants.


Introduction
The leaf is a primary photosynthetic organ in plants having a decisive role in plant growth, development, and survival. Leaf rolling is the key character that determines plant architecture that significantly affects crop yield and the accumulation of dry matter in the later period [1]. Rolled leaves are determined erected leaf phenotype, carbon fixation and gas exchange for photosynthesis, reduced transpiration and canopy interception of solar radiation. Studies have shown that rolled leaf genotypes reduce water loss rate and enhance drought tolerance while semi-rolled leaf plants exhibit more rapid water loss and lower water content in leaves and reduced drought tolerance [2][3][4][5][6]. Drought is a prolonged dry period in the natural climate cycle that impairs rice production and challenging tasks due to its unpredictable and complex nature. Under drought stress, plants undergo a range of biochemical and physiological pathways and excessive accumulations of ROS (reactive oxygen species). Antioxidant enzymes, such as peroxidase (POD), catalase (CAT), and superoxide dismutase (SOD), play a vital role in dealing with ROS accumulation. Abscisic acid (ABA) plays the main role in regulating stomatal closure, while malondialdehyde (MDA) production indicates cell membrane damage from water-deficit stress. Leaf rolling decreases stomatal conductance and reduced water loss under severe water-deficit conditions. Therefore, the development of leaf rolling genotypes may provide a source to enhance rice yield under drought conditions. In China, breeding for rolled leaf genotypes is preferred for high-yielding rice varieties, including super hybrid rice [2,7,8]. Therefore, generation and identification of the genetic components that modulate leaf rolling may provide tools to enhance rice yield.
Many genes in rice control rolled leaf phenotype by regulating the number, size and the arrangement of bulliform cells (BCs), including OsZHD1 [9], ACL1 and ACL2 [10], REL1 [11], REL2 [12], OsCSLD4/NRL1 [13][14][15], RL14 [16], and OsMYB103L [17]. Semi-rolled leaf 1 (SRL1) is allelic to curled leaf and dwarf 1 (cld1), which encodes a putative glycosylphosphatidylinositol-anchored protein and regulates leaf rolling on the adaxial side by increasing the number of BCs. The SRL1 mutant with loss-of-function showed lesser lignin contents and cellulose in the epidermis of BCs [6,18]. SRL2 rice mutants showed curved and narrow leaves mainly expressed in the leaf sheaths, vascular bundles (VB) of leaf blades, and roots and encode a protein of unknown biochemical function. SRL2 is also responsible for the sclerenchyma cells (SCs) or mesophyll cells (MCs), which play a critical role in leaf rolling by severely impairing cuticle development and organization of the epidermis [19].
CRISPR/Cas9 (clustered regularly interspaced short palindromic repeat CRISPR-associated endonuclease 9) is an efficient tool for targeted mutagenesis, which has been successfully applied in rice [20,21]. The CRISPR/Cas9 system has also been employed for the simultaneous editing of multiple genes [22], as well as target gene repression and activation [23].
To detect the changes in proteins induced by the mutations, including internal, external, or both stimuli, proteomic analysis proved to be a powerful tool. The isobaric tags for relative and absolute quantification (iTRAQ) based proteomic analysis is a newly developed technique that can quantitatively analyze protein abundance [24]. It has been extensively used to recognize different biological processes in rice and other plant species [25].
Many prediction tools and experimental approaches have been developed, but the vastly complex mechanisms underlying the cell's reply to genetic mutations are still to be comprehended. Several studies have been performed related to rice rolled leaf phenotypes, but according to our knowledge, there is currently no iTRAQ-based proteomic study available to understand the rolled leaf prototype induced by CRISPR/Cas9. This study is an application of proteomics to assess the outcome of CRISPR-induced events in rolled leaf rice mutants. We generated CRISPR mutants with several morphological imperfections, and our results revealed that the SRL1 and SRL2 loss-of-function mutations led to a rolled leaf phenotype that showed drought tolerance. Our study provides some novel evidence supporting the regulatory module in leaf rolling and drought tolerance in rice.

Plant Material and Growth Conditions
In this study, three restorer rice lines (GXU16, GXU20, and GXU28), were used for genetic transformations. Plants were grown in standard greenhouse conditions (16 h light at 30 • C/8 h dark at 22 • C), and in the experimental field of Guangxi University (22 • N, 108 • E), Guangxi Province, China.

CRISPR/Cas9 Plasmid Construction and Generation of Mutant Plants
Two target sites for each gene were selected by the online tool CRISPR-GE (http://skl.scau.edu. cn/home/) [21]. The first target for SRL1 was designed in the second exon (669-688 bp), and the second target was designed in the third exon (1409-1428 bp) ( Figure S1A,B). The first target of SRL2 was designed in the first exon (1494-1516 bp) and second target in the second exon (1844-1863 bp) ( Figure S1C,D). Targets were selected with high efficiency and low off-target scores (Table S2). The structure of sgRNA's ( Figure S2) was constructed by an online tool, CRISPR-P (version 2.0) (http://crispr.hzau.edu.cn/cgi-bin/CRISPR2/CRISPR). The sgRNA expression cassette was constructed through overlapping PCR. The cloning of the sgRNA expression cassette into the pYLCRISPR/Cas9 binary vector was performed by restriction-ligation reactions. Finally, the ligation product was transformed into Escherichia coli (DH5α) by the heat shock method, and the bacterial solution was applied to LB (Luria-broth) agar plates containing 50 ug/mL kanamycin for 12 h. Simple PCR was performed using specific primers (PB-L and PB-R) and the bacterial solution as a template [22]. Rice calli transformation was completed following the method by Hiei et al. (1994) [26], and calli were transferred to regeneration media to obtain green plants.

Genotyping of T 0 Generation and Off-Target Effect Detection
The genomic DNA of T 0 generation was extracted, and site-specific primers (SRL1T1-F/R, SRL1T2F/R and SRL2T1F/R, SRL2T2F/R) (Table S1) were used to amplify target sites. The product confirming the target fragment size was sent to the Beijing Genomics Institute (BGI) for sequencing. The mutations were analyzed by the online analysis tool, DSDecode (http://skl.scau.edu.cn/ dsdecode/) [27]. The presence/absence of T-DNA was evaluated using Cas9 gene-specific primers (Cas9-F/Cas9-R) ( Table S1). The potential off-target sites were selected from an online tool, CRISPR-GE. Four off-target sites for each target were analyzed by direct Sanger sequencing with site-specific primers (POT1-POT16) (Table S3).

Morphophysiological and Microscopic Analysis
T 2 plants were grown under normal irrigation in the paddy field until maturity to record the data for plant height, panicle number, panicle length, flag leaf length, flag leaf width, grain number per panicle, seed setting rate, 1000-grain weight, grain length, and grain width. In normal irrigation, the standing water was maintained from transplanting to 15 days before maturity. The transpiration rate and stomatal conductance were measured from sixty days old plants during sunny days (09:00 am, 11:00 am, and 15:00 pm) by using the Li-Cor 6400 Portable Photosynthesis System. Leaf samples (0.1 g) were taken, and chlorophyll was extracted with 80% acetone and measured according to the previous method, as reported by Lichtenthaler (1987) [28]. The leaf samples were taken at the tillering stage, and cells were observed and counted under an Olympus CX40 microscope and photographed with mshot-mc50.

Measurement of ABA Contents, MDA Contents, and Antioxidant Enzyme Assays Under Drought Stress
The 15 days old seedlings were grown in a growth chamber (16 h light at 28 • C and 8 h dark at 20 • C with 80% humidity), in nutrient solution containing 10% polyethylene glycol (PEG)-6000. The nutrient solution was refreshed after every seven days, and after five weeks, the plants were treated with 10% (w/v) PEG6000. The ABA level was determined using radioimmunoassay (RIA) [29]. The 0.2 g seedlings were ground with liquid nitrogen and homogenized in 1 ml ddH 2 O and then shaken overnight at 4 • C. The supernatant was used for ABA assay after centrifugation. The MDA (malondialdehyde) content was determined according to a previously established method [29]. The results were expressed as mg g −1 FW (fresh weight) of the seedlings. The seedlings were crushed into a powder in a mortar with the pestle with liquid nitrogen. The crude enzyme of the powder was extracted in 50 mM chilled sodium phosphate buffer (pH 7.8) and 1% (w/v) polyvinylpolypyrrolidone (PVPP) at 4 • C. The homogenate was centrifuged at 12,000 g at 4 • C for 15 min. The SOD activity was determined by the nitro blue tetrazolium (NBT) photochemical reduction by 50% according to the previously established method [30]. CAT activity was assayed by measuring the rate of decomposition of H 2 O 2 at 240 nm, and units were expressed as min −1 g −1 FW, following the previous method [31]. POD activity was assayed according to the previously established method [32], at 470 nm, and units were expressed as min −1 g −1 FW. Wild type (WT) and mutant lines were also grown in controlled conditions to induce drought stress. The drought stress was imposed from 50 days (panicle initiation) to 80 days after sowing (heading). After the termination of drought stress, water was applied evenly to the planting area and grain filling percentage was calculated.

Protein Extraction, Digestion, and Labeling
Total protein was extracted from the leaf samples of WT (GXU16) and its mutant line (GXU16- 19).
Leaf samples were ground fully in liquid nitrogen and homogenized with borax/polyvinyl-polypyrrolidone/phenol (BPP) in a ratio of 1:3 and vortexed for 10 min at 4 • C. Tris-saturated phenol (v:v =1:1) was added and vortexed at 4 • C for 10 min and five volumes of 100 mM ammonium acetate in methanol at 4 • C was used for collection and precipitation of phenolic phase. The sample was centrifuged (20 min, 6000 rpm, 4 • C), and pellets were suspended in 10 mL of methanol soon after collection and repeated in 10 mL of acetone, 10 mL of methanol, and 1 mL of acetone again, respectively. The sample was centrifuged at 12,000 rpm for 20 min and pellets were collected, air-dried, and suspended in 150 L of radioimmunoprecipitation assay (RIPA) lysis buffer containing 1% (w/v) sodium deoxycholate, 0.1% (v/v) TritonX-100, 0.1% (w/v) SDS, 50 mM Tris-HCl (pH 8.0), and 150 mM NaCl. The protein concentration was measured by the BCA method (bicinchoninic acid). Trypsin (Promega, Madison, WI, USA) at 37 • C at a ratio of 1:50 (enzyme/substrate) was used for protein digestion overnight, and iTRAQ labeling was completed following the manufacturer's protocol (Applied Biosystems, Sciex, Foster City, CA, USA). The samples were labeled as (16_GXU)-116 and (GXU16_19)-118. The further experimental procedure was carried out according to the previously established method [33].

Protein-Protein Interaction (PPI) and Data Analysis
Proteins were identified by the Mascot search engine (Matrix Science, London, UK; version 2.3.02), and the Osativa_204 (40869 entries) database was selected. Protein Pilot Software 4.0 (AB SCIEX, Redwood City, CA, USA) was used for data analysis with the 2-fold cutoff for quantitative changes to asses up-and downregulated proteins. The DEPs (differentially expressed proteins) were then imported to the GO (http://geneontology.org/) and KEGG database (http://www.genome.jp/ kegg/pathway.html) for metabolic pathway and biological processes examination. The STRING database (version 10.0) (https://string-db.org/) was searched for protein-protein interaction (PPI) for all upregulated/downregulated proteins. This database contains interactions from earlier published interaction studies, which include physical (direct) and functional (indirect) associations; they are based on computational prediction, from the knowledge transfer among the organisms, and from the aggregated interactions from other databases. PPIs belonging to the uploaded dataset were selected to minimize false positive/negative interactions based on the confidence score, and the lines represent associations among different proteins. The experimentally validated interactions with the highest combined scores were selected as significant and DEPs. The nodes with degree, closeness, and betweenness scores higher than the mean, as calculated by the Cytoscape (version 3.7.2) plugin, Cytohubba, were considered hub nodes. The network of the top ten highly expressed proteins was constructed based on the degree of protein expression. The data analysis for all morphological characters was performed using the SPSS 16.0 Statistical Software Program.

Expression Analysis of Target Genes and Validation of Proteomic Data
Total RNA was isolated from the fresh leaves by using the TaKaRa MiniBEST Plant RNA Extraction Kit according to the manufacturer's instructions. The rice actin gene was used as a control. The expression level of SRL1 and SRL2 was confirmed by using specific primers. A total of five proteins were randomly selected for dynamic transcriptional expression profiles. Primers used in this study are listed in Table S1. Real-time quantitative PCR analysis on a Bio-Rad CFX96 Real-time System using the qRT-PCR kit (Takara PrimeScript RT-PCR Kit-Perfect Real Time), 20 µL reaction system (SYBRT Green 10µL, RNase/DNase-free Water 3.8 µL, 0.6 forward and reverse primers) was undertaken. The amplification reaction cycle was adjusted at 95 • C for 30 s, 95 • C for 5 s, 60 • C for 30 s, 65 • C for 5 s, 40 cycles. The 2 −∆∆Ct method was used to calculate relative expression levels.

Development of Hybrids and Assessment of Main Agronomic Characters
The WTs (GXU16, GXU20, and GXU28) and rolled leaf mutants (GXU16-19, GXU20-8, and GXU28-12) were used as a restorer (R) and crossed with flat-leaf cytoplasmic male sterile (CMS) lines (36A and 52A) to develop F 1 hybrids. All the F 1 hybrids were planted in the field conditions, and the hybrid combinations and parents were used for evaluating the main agronomic characters.

Assembly of Targets in pYLCRISPR/Cas9Pubi-H
The amplification of the sgRNA expression cassette for all targets was confirmed by overlapping PCR ( Figure 1A). The PCR product was mixed and purified by TaKaRa MiniBEST Purification Kit Ver.4.0. The CRISPR/Cas9 binary vector was successfully constructed, and all target sequences were confirmed in the vector ( Figure 1B) by using the specific primers (SP-L1 and SP-R) (Table S1).

Frequency of Mutations and Off-target Analysis
After the CRISPR/Cas9 constructed transformation into rice embryogenic calli, the positive mutant lines were selected by using specific hygromycin phosphotransferase primers (HPTF/HPTR) (Table S1). Nucleotide sequencing analysis of T 0 plants showed that mutant plants showed different mutation rates for each target. For GXU16, a total of 22 plants were obtained in T 0 generation with a mean mutation rate of 22.7% homozygous and 45.4% heterozygous. A total of 19 positive mutant plants was obtained for GXU20, including 19.3% homozygous and 37.5% heterozygous mutations. For GXU28, a total of 25 plants was obtained, with 23.9% homozygous and 47.7% heterozygous mutations ( Figure S3). The DNA of 30 mutant plants was amplified for four most likely sites for each target with the highest-ranking off-target potential. Results showed that there was no off-target mutation were detected (Table S4).

Screening of T-DNA-Free Plants
The results showed that 55 mutant plants were amplified to the Cas9 vector sequence, while 37 plants were not amplified and thus named as T-DNA-free plants ( Figure S4). We selected three homozygous T-DNA free plants (GXU16-19, GXU20-8, and GXU28-12) for further investigations. GXU16-19 showed homozygous mutations with 5 pb, 9 bp, 1 bp, and 4 bp deletions on first, second, third, and fourth target sites, respectively. Homozygous mutant GXU20-8 showed 10 bp, 8 bp, 6 bp, and 6 bp deletions on first, second, third, and fourth target sites, respectively. The mutant line GXU28-12 resulted in homozygous mutations with 11 bp, 7 bp, 5 bp, and 3 bp deletions on first, second, third, and fourth target sites, respectively ( Figure 2A). Sequencing results confirmed that homozygous mutants showed stable mutations in T 1 and subsequent T 2 generations ( Figure 2B).

Phenotyping Under Normal and Drought Conditions
The mutant lines showed a rolled leaf phenotype until maturity, while the WT leaves were flat during the whole phase of development ( Figure 3A,B). The leaf rolling index of the SRL1 mutant was 87.5%, while WT showed a flat-leaf phenotype ( Figure 3C). Mutant plants showed a lower transpiration rate than WT at the three different time points during sunny days ( Figure 3D). The stomatal conductance was significantly lower in mutant plants ( Figure 3E). Under normal conditions, the plant height of the mutant lines GXU16, GXU20, and GXU28 was about 32.0 cm, 34.9 cm, and 37.8 cm, lower than that of WT at mature stage, respectively (Table 1); the panicles number were increased from 7.7 to 11.3; the panicle length was reduced by 42%; the flag leaf length and width was decreased by 26.5% and 27.3%, respectively; the grain number per panicle and seed setting rate was decreased by 21.0% and 23.9%, respectively; the 1000 grain weight was reduced by 14.2%; the grain width and grain length of the mutants were reduced by about 16% and 15.78%, respectively, compared with WT (Table 1). Stomatal conductance was higher in WT as compared to mutant plants ( Figure 3E). The grain filling percentage under normal conditions was 83% and 80% in WT and mutant plants, respectively, while under severe drought, the percentage of filled grain in WT plants was nearly 2.5%, whereas that in mutant plants was 27.5% ( Figure 3F). LT (leaf type); WT (wild type); PH (plant height); PN (panicle numbers); PL (panicle length); FLL (flag leaf length); FLW (flag leaf width); GNPP (grain number per panicle); SSR (seed setting rate); GW (1000-grain weight); GL (grain length); GW (grain width); GYPP (grain yield per plant). * represent a significant difference, and ns represents a non-significant difference, Student's t-test, p ≤ 0.01. The data is the mean of five independent samples from T 2 generation.

Microscopic Analysis
There was a significant difference between WT and mutant leaf cells ( Figure 4A,B). Microscopic analysis showed that the stomata number in mutant plants was less than WT ( Figure 4C,D).

Mutant Plants Enhance ABA Content and Improve Antioxidant Enzyme Activities Under Drought Stress
In comparison to WT, CRISPR mutants showed a higher ABA level while MDA content was increased in WT after drought stress ( Figure 5A,B). The SOD, POD, and CAT activities were higher in the mutant plants, especially SOD and CAT ( Figure 5C-E). A drought tolerance assay showed that the leaves of WT became wilted after 6 days while the leaves of mutant plants were not wilted even after 15 days of drought treatment. After 15 days of drought treatment, water was again applied to the plants to recover normal growth. The survival rate of mutant plants was increased by 45% as compared to WT ( Figure 5F). WT plants showed normal phenotype under controlled conditions, with more physiological damage than mutant plants after drought stress and re-watering ( Figure 5G-I). represents wild type and mutants, respectively. Data are presented as means ± sd (n = 5). * and ns indicates a significant and non-significant difference, respectively, Student's t-test, p ≤ 0.01.

Protein Identification and Quantitation of Rolled Leaf Mutant and Wild Type
A total of 441,569 spectra were produced for both samples. The raw library size was 535,647 and 806,761 and normalized library size was 667,640.122 and 674,665.262 for WT (GXU16) and CRISPR mutant (GXU16-19), respectively. The box-and-whisker plot of proteomic count data distribution showed a difference between the medians of both samples ( Figure 6A). The distribution of protein count data histogram showed that the frequency of proteins was highest at 5 to 8 of log 2 (counts + 1) ( Figure 6B). By arranging the data, we identified 6538 proteins, and after adjusting the filter cutoff 10, we obtained 6452 proteins. The comparison showed 377 DEPs, 270 downregulated, and 107 upregulated ( Figure 6C). Some DEPs were found to be highly associated with leaf rolling and abiotic stress tolerance and considered as candidate proteins ( Table 3). The detailed information of all identified proteins is listed in (Additional File 2). The top twenty proteins were also selected with higher adjusted p-value cutoff to select proteins with the highest mean expression. A2YP90 (uncharacterized protein), B8BNS8 (ribosome assembly factor mrt4), and A2YI21 (PsbP domain-containing protein) showed higher expression levels in the mutant line with lower expression levels in WT. Q1KL27 (peptidyl-prolyl cis-trans isomerase), and B8BFK2 (uncharacterized protein) showed higher expression levels in WT and lower expression levels in mutant line ( Figure 6D).

Gene Ontology (GO) Annotation of Up-and Downregulated Proteins
GO annotation for upregulated proteins related to molecular function (MF) showed that most of the DEPs participated in binding, catalytic, transporter, and antioxidant activities. For cellular component (CC), the upregulated DEPs were largely enriched in the cell part, membrane part, organelle part, organelle, protein-containing complex, and extracellular region. For the biological process (BP), the upregulated DEPs were mainly enriched in the cellular process, metabolic process, biological regulation, and response to a stimulus ( Figure 8A). GO terms for MF showed that most of the downregulated DEPs were involved in binding and catalytic activity. CC terms showed that the downregulated DEPs were mainly enriched in the cell part, organelle, membrane part, protein-containing complex, and organelle part. The GO terms for BP showed that most of the downregulated DEPs were involved in the cellular and metabolic processes ( Figure 8B).

Kyoto Encyclopedia of Genes and Genomes (KEGG) Pathway Enrichment Analysis of Up-and Downregulated Proteins
Enrichment analysis was performed with p-value ≤ 0.05 as significantly enriched in DEPs. KEGG analysis of upregulated DEPs revealed metabolic pathways, biosynthesis of secondary metabolites, and carbon metabolism, while other pathways showed a smaller number of proteins involved ( Figure 9A). The downregulated DEPS were significantly enriched in metabolic pathways, ribosome, biosynthesis of secondary metabolites, glutathione metabolism, carbon metabolism, fatty acid metabolism, biosynthesis of amino acids, and carbon fixation, while other pathways showed few proteins involved ( Figure 9B).

Expression Analysis of Target Genes Quantitative Real-Time-PCR Validation of DEPs
Results showed that the expression of SRL1 and SRL2 was significantly suppressed in mutant lines ( Figure 11A). To verify the proteomic results, five DEPs associated genes were randomly selected, including proteasome subunit alpha type-2 (PAB1), the mediator of RNA polymerase II transcription subunit 17 (MED17), ethylene response factor (Snorkel 2), pathogenesis-related protein (Pr1b), and pathogenesis-related protein (Pr1a). The RT-qPCR results were consistent with the proteomic analyses ( Figure 11B).

Performance of Semi-Rolled Leaf Hybrids Produced by Wild Type and Mutant Rolled Leaf Restorers
The performance of the hybrids was investigated and was compared with those hybrids produced by the WT restorers crossed with the same CMS lines. The hybrids (36A/GXU16, 52A/GXU16, 36A/GXU20, 52A/GXU20, 36A/GXU28, and 52A/GXU28) produced from WT restorer lines showed flat-leaf type, but the hybrids (36A/GXU16-19, 52A/GXU16-19, 36A/GXU20-8, 52A/GXU20-8, 36A/GXU28-12, and 52A/GXU28-12) developed from mutant restorer lines showed the semi-rolled leaf phenotype ( Figure S5). The results revealed that the other agronomic characteristics, such as plant height, panicle length, seed set rate, grain weight, grain length, and grain width, were the same with the hybrids produced by WT restorers (Table 4). However, two major yield components, the panicle number and grain number per panicle, were significantly increased from 7.8 to 10.5 and 134 to 170, respectively, in the hybrids produced by mutant restorers. The yield per plant was also significantly increased from 29.7g to 44.5g, in hybrids produced by mutant restorers. These results showed that the hybrids produced by mutant restorers resulted in significant yield heterosis (Table 4). LT (leaf type); SR (semi rolled); WT (wild type); PH (plant height); PN (panicle numbers); PL (panicle length); GNPP (grain number per panicle); SSR (seed setting rate); GW (1000-grain weight); GL (grain length); GW (grain width); GYPP (grain yield per plant). Data are the mean of five replicates. * and NS represents the significant and non-significant differences, respectively. Student's t-test, p ≤ 0.01.

Discussion
CRISPR/Cas9 based mutations were successfully generated, and homozygous, heterozygous, and WT plants were obtained in the T 0 generation. These results suggest the efficiency of the CRISPR/Cas9 system in rice and the generation of homozygous mutants in the T 0 generation without any off-target events. Previous studies also suggested that Cas9 seldom tempts off-target mutations and concluded that the off-target rate in rice is very low [20,21]. In our study, we screened the T-DNA free plants at a frequency of 40.2%, and expression of SRL1 and SRL2 was significantly depressed in the mutant plants. Under normal conditions, mutant lines displayed an 87% leaf rolling index, decreased agronomic character, transpiration rate, stomatal conductance, chlorophyll content, VB, number of stomata, while panicles number, and BCs were significantly increased. Rolled leaf mutants possess low net photosynthesis rate and chlorophyll content, resulting in less yield compare to WT at low plant intensity [3]. The mutant plants conferred drought tolerance and showed higher ABA content, anti-oxidants activities, and survival rate at the seedling stage. Under drought stress, the mutant plants showed higher grain filling percentage than WTs at maturity. Previous studies also showed that leaf rolling and increased BCs have a positive role in stress tolerance by adjusting transpiration [34,35]. ABA is significantly induced by drought stress and plays a role in triggering the enhanced ROS generation [36]. In previous studies, mutations in the rice OsCHR4 gene induce narrow and rolled leaves with increased cuticular wax resulting in reduced water loss rate and enhanced drought tolerance [5]. The knockdown of rice microRNA166, using the short tandem target mimic (STTM) system, resulted in rolled-leaf phenotype, morphological changes that confer drought resistance. STTM166 plants had smaller BCs, reduced stomatal conductance, and showed decreased transpiration rates [37].
To ensure the fidelity of biological conclusions, it is important to confirm the success of any gene editing at the protein level. Recently, researchers have successfully edited OsGA20ox2 in rice and studied the proteome-wide changes that were induced [38]. Here, the iTRAQ-based proteomic approach was applied both to confirm protein level changes induced by gene editing and to measure the proteome-wide changes that occur in response to the system perturbation in a single experiment. A total of 377 differentially abundant proteins were identified, containing 107 upregulated proteins and 270 downregulated proteins. The KEGG analysis showed that up-and downregulated proteins were mostly enriched in the biosynthesis of secondary metabolites, metabolic pathways, and carbon metabolism. Metabolic responses play a significant role in drought stress within a variety of plant species, and different metabolic adjustments were observed under water-deficit stress [39][40][41]. Secondary metabolites, such as phenolic acid and flavonoids, described as markers of abiotic and biotic stress tolerance, are ubiquitous in the plant kingdom [42]. Carbon metabolism helps the plants decrease photosynthetic carbon assimilation, which results in enhanced ROS production by disrupting cellular homeostasis [43].
GO analysis displayed that most of the DEPs were involved in binding, catalytic, transporter, and antioxidant activities. Proteins having binding and catalytic activities, such as DNA-binding proteins and tRNA-binding proteins, may be involved in the regulation of intra-and extra-cellular signals associated with sugar biosynthesis, development, and growth. Antioxidant components enhancement is a strategy to prevent oxidative damage under drought stress [36].
The proteins related to stress, including POD, SOD, ABA-ripening protein 5, HSP domain-containing protein, hydrolase_4 domain-containing protein, ParB domain-containing protein, NAD(P)-bd domain-containing protein, were upregulated. Antioxidants, POD and SOD are vital ROS-scavenging components in plants, and their increased expression resulted in drought tolerance in rice [36]. ABA-induced proteins are regulated by abiotic stresses, which are hydrophilically charged, low-molecular-weight, plant-specific proteins [44]. Heat shock proteins (HSPs) are ubiquitous proteins that function as molecular chaperons encoded by a multigene family and play important roles in abiotic stress tolerance [45]. NAD(P) domain proteins and respiratory burst oxidase homolog proteins (rbohs) are homologs that play a significant regulatory role in promoting ABA-induced stomatal closing, Ca 2+ increase, and ROS production in plants [46].
Lateral organ boundaries (LOB) domain protein-like (LBD) (Q8L3S3, Q7XGL4, Q852M3, and Q8L4M5), glycosyltransferase, phosphoglycerate kinase, peptidylprolyl isomerase, and acyl-coenzyme A oxidase, were significantly accumulated in the rolled leaf mutant. The LBD gene family is a plant-specific transcription factor family with a coiled structure, having a potential role in lateral organ development, such as flowers, leaves, and roots. Forty-three members of the LBD gene family in rice have been identified. LOB domain proteins are negative regulators of cell proliferation in the adaxial side of leaves. [47]. Glycosyl hydrolases, acyl-coenzyme A oxidase, peptidylprolyl isomerase, and Phosphoglycerate kinase are common enzymes in plants having a crucial role in the degradation of biomass such as starch, cellulose, and hemicellulose [6].
In protein network analysis and hub-protein analysis, we found some proteins related to abiotic stress and leaf development. The proteins, including WUSCHEL-related homeobox 4, putative auxin transporter-like protein 4, LOB domain proteins, probable indole-3-acetic acid-amido synthetase GH3.4, protein Monoculm 1, AP2/EREBP transcription factor baby boom1 and AP2 domain-containing protein, were top co-expressed ( Figure 10). Wuschel (WUS) plays a role in shoot apical meristem (SAM) maintenance, and its mutants display aborted SAM maintenance during embryonic and later developmental stages [48]. Monoculm 1 (MOC1) is the first cloned gene which is responsible for tiller formation in rice [49]. AP2/EREBPs (Apetala2/ethylene-responsive element-binding proteins) are involved in regulatory networks of response to pathogen attack, hormones, and environmental signals involving ERFs (ethylene-responsive factors) and DREBs (dehydration responsive element binding proteins) [50]. Hub-proteins reveal multiple dysregulated pathways. Highly connected hub nodes are central to a network's architecture and function. Unsupervised hierarchical clustering analysis based on the hub protein expression profiles showed that the identified top hub proteins serve as a molecular signature to differentiate CRISPR/cas9 mutants and WTs. The proteomic results were validated by performing RT-qPCR for five selected genes (PAB1, MED17, Snorkel 2, Pr1b, and Pr1a). As expected, the results were found to be correlated with proteomic data ( Figure 11A). In the top twenty proteins, the ribosome assembly factor mrt4 and PsbP domain-containing protein were highly expressed in mutant while peptidyl-prolyl cis-trans isomerase showed higher expression level in WT. The psbP proteins have been reported to be differentially expressed under abiotic stress, while ribosomal proteins have been differentially regulated under both biotic and abiotic stresses [51,52]. Peptidyl prolyl cis-trans isomerase (PPIase) plays a role under drought stress by different regulatory pathways in the presence of both cyclophilins (Cyps) and FK506-binding proteins (FKBPs) [53].
The hybrids produced from mutant restorers showed semi-rolled leaf phenotype, increased panicle number, grain number per panicle, and grain yield per. This indicates that the semi-rolled leaf hybrids have more reasonable plant types and the use of one roll leaf parent crossing with another flat-leaf parent to produce a semi-rolled leaf hybrid is helpful to produce an ideal plant type. The hybrids showed an intermediated leaf phenotype to the WTs and mutant plants, which may be due to incomplete dominance. These incomplete dominance loci should play an important role in leaf rolling, which implies that dominance is an important contributor to heterosis in hybrid breeding. The generation of numerous F 1 crossed will provide efficient ways to study the mechanism underlying this phenomenon. This study opens an era of research that deserves more attention in the future, including the generation of more suitable mutant populations via CRISPR/Cas9 and the development of more powerful genomic and statistical analysis frameworks.
Our results provide a unique molecular mechanism underlying rice leaf rolling. According to our knowledge, many studies have been carried out on leaf rolling, but there was no study found based on the CRISPR/Cas9-generated rolled leaf mutants and iTRAQ based proteomic analysis of directed mutations. The employment of CRISPR/Cas9 and proteomics together will open a new way to understand the effect of single or multiple gene mutations on the whole proteome. Many identified proteins in our study were different from previous studies because the rolled leaf mutants were generated by nucleases-based mutagenesis at different frequencies of mutations. The hybrids produced in this study are usually more valuable to commercial growers and have worthwhile advantages to plant breeders to produce new and better cultivars by employing modern technologies. Further studies are warranted to explore the molecular regulatory mechanism underlying rolled leaf mutants and to devise strategies for improving the rice phenotype and drought tolerance.

Conclusions
In summary, the present study added some novel insights into the functional role of SRL1 and SRL2. CRISPR mutants influence the expression of many proteins that directly affect the different BP, CC, and MF pathways. The up-and downregulation of stress-related proteins showed that rolled leaf mutants play a role in abiotic stress tolerance. Further study with these mutants will help to gain a better understanding of regulatory pathways involved in leaf rolling and abiotic stress tolerance. The downregulation of LBD proteins in rolled leaf mutants will bring a sharp focus in future molecular breeding of rice. Identified uncharacterized and known hub-genes will provide the correlative support to study the leaf rolling and abiotic stress tolerance in rice. Undoubtedly, homozygous mutants with heritable mutations lacking off-target effects and T-DNA will provide the basic breeding material for the development of new rice varieties for water deficit areas. We provided the molecular evidence that rolled leaf CRISPR mutants affected the expression of different proteins and showed that drought tolerance could be improved by developing rolled leaf genotypes. Moreover, the semi-rolled leaf hybrids produced from mutant restorers provide the basic idea to fast track the breeding system and to improve the plant architecture and yield. These results might pave the way for understanding the gene action and inheritance pattern of CRISPR/Cas9 generated mutants.

Supplementary Materials:
The following are available online at http://www.mdpi.com/2073-4395/9/11/728/s1, Table S1: List of primers used in this study; Table S2: Efficiency score and positions all targets; Table S3: Primers used for off-target analysis; Table S4: Detection of mutations on the putative off-target sites; Figure S1: (A) Sequence and position of both targets of and target-specific primers in SRL1. (B) The structure of the SRL1 gene and an illustration of both target sites in exon2 and exon3. (C) Sequence and position of both targets of and target-specific primers in SRL2. (D) The structure of SRL2 and illustration of both target sites in exon1 and exon1. Yellow highlighted sequences represent target sites; PAM (protospacer adjacent motif) is green highlighted, and target-specific primers are purple highlighted. The black box represents the exon regions, and the lines between them indicate introns. T1 and T2 represent target1 and target2, respectively; Figure S2: Schematic representation of secondary structures of (A) sgRNA1; (B) sgRNA2; (C) sgRNA3 and (D) sgRNA4 used in the experiment. Figure S3: (A) Mutation rates of all targets in T 0 generation; Figure S4: The gel electrophoresis results of T-DNA free plants; Figure