• <tr id="yyy80"></tr>
  • <sup id="yyy80"></sup>
  • <tfoot id="yyy80"><noscript id="yyy80"></noscript></tfoot>
  • 99热精品在线国产_美女午夜性视频免费_国产精品国产高清国产av_av欧美777_自拍偷自拍亚洲精品老妇_亚洲熟女精品中文字幕_www日本黄色视频网_国产精品野战在线观看 ?

    Identification of microRNAs regulating grain filling of rice inferior spikelets in response to moderate soil drying post-anthesis

    2022-08-16 09:25:24ZhenningTengYinkeChenYouqingYunYqiongPengYkeYiHuihuiYuZhenxieYiJinchngYngYnPengMeijunDunJinhuZhngNenghuiYe
    The Crop Journal 2022年4期

    Zhenning Teng, Yinke Chen, Youqing Yun, Yqiong Peng, Yke Yi, Huihui Yu, Zhenxie Yi,Jinchng Yng, Yn Peng, Meijun Dun,b,, Jinhu Zhng, Nenghui Ye,b,e,

    a College of Agriculture, Hunan Agricultural University, Changsha 410128, Hunan, China

    b Hunan Provincial Key Laboratory of Rice Stress Biology, Hunan Agricultural University, Changsha 410128, Hunan, China

    c Changjun High School of Changsha, Changsha 410128, Hunan, China

    d Key Laboratory of Crop Genetics and Physiology of Jiangsu Province, Yangzhou University, Yangzhou 225009, Jiangsu, China

    e Department of Biology, Hong Kong Baptist University, Kowloon 999077, Hong Kong, China

    f School of Life Sciences and State Key Laboratory of Agrobiotechnology, the Chinese University of Hong Kong, Shatin 999077, Hong Kong, China

    Keywords:Moderate soil drying Grain filling Inferior spikelet MicroRNA Rice

    A B S T R A C T The grain filling of inferior spikelets is much less complete than that of superior spikelets in rice cultivars with large panicles and numerous spikelets and is promoted by moderate soil drying(MD)post-anthesis.A growing body of evidence has shown that microRNAs function in regulating grain development.However, little is known about the mechanism of microRNA control of grain filling of inferior spikelets in response to MD. In this study, grain filling of inferior spikelets was promoted by MD treatment in Nipponbare.Small-RNA profiling at the most active grain-filling stage was conducted in inferior spikelets under control (CK) and MD treatment. Of 521 known and 128 novel miRNAs, 38 known and 9 novel miRNAs were differentially expressed between the CK and MD treatments. Target genes of differentially expressed miRNAs were involved in multiple developmental and signaling pathways associated with catalytic activity, carbohydrate metabolism, and other functions. Both miR1861 and miR397 were upregulated by MD,leading to a decrease in OsSBDCP1 and OsLAC,two negative regulators of SSIIIa activity and BR signaling,respectively.In contrast,miR1432 abundance was reduced by MD,resulting in upregulation of OsACOT and thus an elevated content of both ABA and IAA.These results suggest that both starch synthesis and phytohormone biosynthesis are regulated by differentially expressed miRNAs in inferior spikelets in response to MD treatment. Our results suggest the molecular mechanisms by which miRNAs regulate grain filling in inferior spikelets of rice under moderate soil drying, providing potential application in agriculture to increase rice yields by genetic approaches.

    1. Introduction

    Rice (Oryza sativaL.) is a staple food for half of the global population and comprises 23%of the calories consumed worldwide[1].China is the largest rice producer and consumer in the world and feeds 20% of the global population with only approximately 5% of the planet’s water resources and 7% of its arable land [2]. Rapid population growth and economic development demand increased food production [3]. Hybrid rice can yield 25% more than conventional cultivars and increase China’s food production capability[4,5].However,hybrid rice cultivars often fail to achieve their high yield potential, owing to poor grain filling of later-flowering inferior spikelets, in particular in the newly bred ‘‘super rice” with numerous spikelets on its panicles [6].

    Poor grain filling in inferior spikelets is common in cereal crops and is sensitive to environmental changes and artificial treatments,including air temperature and CO2concentration in the atmosphere [7,8], post-anthesis moderate soil drying (MD) [9,10],chemical applications such as spermine and spermidine [11], and plant hormones and their inhibitors [12-14]. Apical primary branches,which flower earlier and generate larger and heavier kernels, are ‘‘superior” spikelets. In contrast, spikelets located on the proximal secondary branches that are pollinated later and generate smaller or even no kernels, are defined as ‘‘inferior” spikelets(Fig. S1) [11]. Approximately 90% of a rice kernel is made up of starch, making starch biosynthesis the most important biological process in kernels during grain filling. Manipulations of this pathway by modern genetic techniques will alter grain filling. Positional cloning [15] revealed that the ADP-glucose pyrophosphorylase (AGP) large subunit 2 (OsAPGL2) gene was defective in thegif2(grain incomplete filling 2) mutant, resulting in a slower filling rate and a marked decrease in starch content in comparison with kernels from wild type (WT) plants. However,overexpression of the ADP-glucose pyrophosphorylase gene (glgC)in CS8 transgenic rice did not increase grain weight, owing to the upregulation of the starch binding domain-containing protein 1(OsSBDCP1) gene, encoding a repressor of starch synthase IIIa(SSIIIa) enzyme by binding with SSIIIa protein, suggesting one or more rate-limiting steps downstream of AGPase in the starch synthesis pathway [16]. Indeed, mutation ofOsSSIIIachanged starch granule formation and led to white-core floury endosperm in rice grain [17], implying that regulation of enzyme activity of starch synthesis is also critical for grain filling. MD treatment has been shown [10,18] to promote enzyme activities in sucrose-to-starch conversion in inferior spikelets, mediated by an elevated abscisic acid(ABA)content.Indole-3-acetic acid(IAA)and brassinosteroids(BRs) were also reported [19-22] to act in regulating grain filling by regulating endosperm cell division during the early stage of grain filling. These investigations suggest that plant hormones and starch biosynthesis-associated enzymes are both necessary for grain filling.

    MicroRNAs (miRNAs), a type of endogenous noncoding small RNA produced from stem-loop-structured precursors, function in plantlet morphogenesis [23,24], stress responses [25], nutrient absorption [26,27] and grain or fruit development [20,28-30]. In rice, dynamic expression of various miRNAs during grain filling suggests the involvement of miRNA in regulating this fundamental process [31,32]. Small-RNA sequencing [33] showed that many miRNAs were differentially expressed in inferior versus superior spikelets, including miR167, which in turn regulated the accumulation of IAA and resulted in poor grain filling in inferior spikelets.Overexpression of miR397 increased grain size and promoted panicle branching by suppressing the expression of the laccase-like protein (OsLAC) gene, which encodes a laccase protein and is involved in sensitivity to brassinosteroids (BRs), resulting in an increase in overall grain yield [20]. In contrast, suppression of rice miR1432 increased grain weight by upregulating the expression of the acyl-CoA thioesterase (OsACOT) gene. Contents of both auxin and abscisic acid were increased in an miR1432-defective mutant and in OXmACOT plant whose miR1432 target site was mutated[34]. osa-miR1848 directs the cleavage of obtusifoliol 14αdemethylase (OsCYP51G3) mRNA to regulate the biosynthesis of phytosterol and BR. Increased osa-miR1848 and decreasedOsCYP51G3expression reduced BR content, resulting in a BRdeficiency phenotype with respect to leaf angle, grain size, and seed quality [35]. These studies indicated that different miRNAs were key regulators for plant hormones in rice kernels and thus acted in regulating grain development.

    Post-anthesis MD has been shown to promote grain filling of rice. miRNAs also function in controlling seed development. However, the function of miRNAs that are mobilized by MD treatment in regulating grain filling of inferior spikelets in rice has not been reported.The objective of the present study was to identify differentially expressed miRNAs in inferior spikelets under MD treatment by small RNA-seq at the most active stage of grain filling.Our results will be the first to identify the miRNA-mediated grain filling of inferior spikelets under MD treatment, which might be a target for engineering yield improvement in crops.

    2. Materials and methods

    2.1. Plant materials and growth conditions

    Experiments were conducted in the greenhouse at Hunan Agricultural University during the rice growing season(June-October).The WT rice (Oryza sativaL. subspeciesjaponica) cultivar used in this study was Nipponbare (Nip). Seeds were germinated in darkness at 28°C for 2-3 days until the roots measured 0.5 cm.Germinated seedlings were sown in the paddy field. Thirty-day-old seedings were then transplanted into plastic pots. Seedings were raised at a hill spacing of 0.15 m × 0.15 m with two seedlings per hill. Each treatment was applied to 20 pots and 12 hills per pot. Fertilizer application followed normal agricultural practice,as described by Zhang et al. [10].

    2.2. Treatments and sampling

    Nip rice was grown in plastic pots filled with paddy soil and sheltered under a plastic shelter. All plants were well-watered until 6 days after anthesis (DAA). The application of water to the pots was then adjusted to produce two treatments as described by Wang et al. [36]. The well-watered control (CK) plants were watered daily to maintain a water depth of 2-5 cm. For the moderate drying(MD)treatment,pots were kept at a soil water potential of -25 kPa using tensiometers manufactured by Institute of Soil Science, Chinese Academy of Sciences, Nanjing, China. Two tensiometers were installed at 15 cm soil depth in each pot.Readings of the tensiometers were recorded four times daily from 8:00 AM to 4:00 PM.When soil water potential dropped to a mean value of-25 kPa,0.5 L of tap water was added to maintain soil moisture.

    Three hundred panicles that headed on the same day in each treatment were tagged to give an accurate record of the flowering date. Twenty tagged panicles from each treatment were sampled every-three days from anthesis to maturity. The sampled panicles were divided into three groups as three replicates. Superior and inferior kernels were separated from the panicles following Chen et al. [11,37].The sampled kernels were dried at 80 °C to constant weight for 72 h,weighed,and stored for measuring non-structural carbohydrates(NSC)and quantifying grain filling rate.A further 30 tagged panicles were sampled at 6, 9, and 15 DAA, and superior and inferior kernels were separated from the panicles, immersed in liquid nitrogen,and stored in a-80°C freezer for RNA extraction and measurement of ABA, IAA, and enzyme activities. Twelve plants were sampled from each treatment at maturity for yieldtrait analysis. Numbers of filled and empty kernels were recorded and the grain weight (mg) of superior and inferior spikelets and grain yield per plant (g) were calculated.

    2.3. Assays of grain weight, carbohydrate, ABA, IAA and enzyme activity

    The grain weight of superior and inferior kernels and the content of soluble sugar and starch of the inferior kernels under CK and MD treatments were measured.A total of 75 superior and inferior kernels each were used for measurement of grain dry weight.The grain-filling process was fitted with the Richards [38] growth equation following Zhu et al.[39].The samples used for measuring the content of soluble sugar and starch were ground into powder to pass a 100-mesh sieve. To a 15 mL centrifuge tube, 50 mg of ground sample was added with 5 mL of 80%(v/v)ethanol and kept in a water bath at 80°C for 30 min.The tube was then centrifuged at 4500 r min-1for 30 min after cooling.The supernatant was collected and the extraction was repeated three times. The sugar extract was then diluted to 25 mL with distilled water, and the soluble sugar and sucrose contents were determined following Yang et al. [39]. The residues left in the centrifuge tubes after extracting sugar were dried at 80 °C for starch extraction with HClO4following Yang et al. [40].

    ABA and IAA contents were assayed by Wuhan MetWare Biotechnology Co., Ltd. (Wuhan, China). Inferior kernels were frozen in liquid nitrogen, ground into powder, and extracted with methanol/water/formic acid(15:4:1,v/v/v).The combined extracts were evaporated to dryness under a nitrogen gas stream,reconstituted in 80% (v/v) methanol,and filtered (polytetrafluoroethylene,0.22 μm; Anpel, Shanghai, China). The sample extracts were separated on an LC-ESI-MS/MS system (HPLC, Shim-pack UFLC Shimadzu CBM30A system, www.shimadzu.com.cn; MS, Applied Biosystems 6500 Triple Quadrupole, www.appliedbiosystems.-com.cn).

    Enzyme activities of StSase with three biological replications were measured with StSase reagent kits (Beijing Solarbio Science& Technology Co., Ltd.). The activity of StSase was measured after grinding 0.1 g grain samples in 1 mL of extraction solution and centrifuging for 10 min at 10,000 r min-1at 4°C.Supernatant(200 μL)was mixed with working solution according to the instructions and the absorbance at 340 nm was recorded. Activities of enzymes were expressed as U g-1min-1.

    2.4. Small-RNA sequencing, identification of miRNAs, and target gene prediction

    Inferior kernels sampled at 9 DAA under CK and MD treatments were used for small-RNA sequencing with three biological replicates. Total RNA was extracted with TRIzol reagent according to the manufacturer’s instructions(Invitrogen,CA,USA)and genomic DNA was removed with DNase I RNase-free (TaKaRa, Dalian,China). RNA concentration was measured using ND-2000 (Nano-Drop Technologies,MA,USA)and RNA integrity was assessed with a 2100 Bioanalyzer (Agilent Technologies, Santa Clara, CA, USA).Only high-quality RNA samples (OD260/280= 1.8-2.2,OD260/230≥2.0, RIN ≥7, 28S:18S ≥1.0, >1 μg) were used to construct sequencing libraries.

    Total RNA(3 μg per sample)was used as input material for the small-RNA library. Sequencing libraries were generated with a NEBNext Multiplex Small RNA Library Prep Set for Illumina (NEB,USA). Activated 5′and 3′adaptors were ligated to the total RNA(Illumina, San Diego, CA, USA). Adaptor-ligated RNA was transcribed into first-strand cDNA using reverse transcriptase and random primers. PCR was performed for 11-12 cycles using primers complementary to the two adaptors for library enrichment. Fragments of appropriate size were isolated on a 6% Novex TBE PAGE gel. Each library was quantified by TBS380. Bridge PCR amplification was performed using the cBot Cluster Generation System.Amplified cDNA fragments were sequenced with an Illumina HiSeqxTen/NovaSeq6000 instrument. The raw data were first quality-controlled with the fastx toolkit software (http://hannonlab.cshl.edu/fastx_toolkit/) to obtain clean small-RNA reads by removal of low-quality bases (with Sanger base quality <20). All identical sequences of sizes ranging from 18 to 32 nt were then BLAST searched against the Rfam database (http://rfam.sanger.ac.uk/) to remove non-miRNA sequences (ribosomal RNA (rRNA),transfer RNA (tRNA), small nucleolar RNA (snoRNA), etc. The remaining reads were used to predict known miRNAs via BLAST search of miRbase (http://www.mirbase.org/), and novel miRNAs were predicted with Mireap (http://sourceforge.net/projects/mireap) software. The expression level of each miRNA was calculated by the transcripts per million reads (TPM) method. Significantly differentially expressed (DE) miRNAs were extracted withPadjust <0.05 and |log2FC| ≥1 by DEseq2.

    psRobot (http://omicslab.genetics.ac.cn/psRobot/index.php),TargetFinder (https://github.com/carringtonlab/TargetFinder), and miRTarBase (http://mirtarbase.mbc.nctu.edu.tw/php/index.php)were used to predict potential miRNA targets. Only targets that were found by all three programs were accepted as target genes of miRNAs. The predicted target genes were aligned with BLAST(http://blast.ncbi.nlm.nih.gov/) and annotated with GO (http://www.geneontology.org/) and KEGG (http://www.genome.jp/kegg/)database.

    2.5. Quantitative real-time PCR and stem-loop qRT-PCR

    For qRT-PCR, total RNA extract from kernels was isolated and reverse transcribed to first-strand cDNA with a HiScript firststrand cDNA synthesis kit (Vazyme, Piscataway, NJ, USA). Transcript levels of selected genes were measured by qRT-PCR using StepOne Plus(Applied Biosystems,MA,USA)with Hieff qPCR SYBR Green Master Mix (Yeasen Biotech, Shanghai, China).Actin1was used as an internal standard.

    qRT-PCR analysis of miRNAs was based on stem-loop reverse transcription PCR. RNA samples were reverse-transcribed with a miRNA 1st Strand cDNA Synthesis Kit (Vazyme). Expression profiles of miRNAs were examined with a miRNA Universal SYBR qPCR Master Mix (Vazyme) using StepOne Plus. U6 snRNA was used as the reference gene for miRNA validation in rice. For each sample,the mean value of relative expression level from three qRT-PCRs was used to calculate transcript abundance. The primers used for qRT-PCR are listed in Table S1.

    2.6. Statistical analysis

    SPSS 23 software (IBM Corp., Armonk, NY, USA) was used to perform a significant difference test. Post-hoc comparisons were tested by Duncan’s test.

    2.7. Accession numbers

    All small-RNA raw sequence reads have been deposited in NCBI under accession number PRJNA736574. Sequence data can be found in the GenBank data libraries or RGAP (Rice Genome Annotation Project) database under the following accession numbers:OsSBDCP1, LOC_Os01g63810;OsSSIIIa, LOC_Os08g09230;OsACOT,LOC_Os04g35590;OsGH3-4, LOC_Os05g42150;OsGH3-8, LOC_Os07g40290;OsLAC, LOC_Os05g38420; andActin1, LOC_Os03 g50885.

    3. Results

    3.1.Physiological characterization of inferior spikelets under moderate soil drying

    The dynamic grain filling process in spikelets is based on the order within a single panicle. The superior kernels showed no significant difference in weight between the CK and moderate soil drying (MD) treatments (Fig. S2a; Table S2). However, exposure to MD treatment increased inferior grain weight by 5.07%(Fig. 1a), and grain yield by 14.97% relative to the CK treatment(Fig.1b;Table S2).The increased seed-setting rate also contributed to increased yield(Fig.S2b).As shown in Fig.1c and d,grain filling and grain-filling rate of inferior spikelets were increased by MD treatment, with the grain-filling period shortened by 2.01 days and the grain-filling rate increased by 22.61% compared with CK at 9 DAA (Fig. 1c, d). However, there was no significant difference between CK and MD treatments in superior kernels (Fig. S2c, d).During grain filling of inferior spikelets,the NSC content and starch content were significantly higher under the MD treatment(Fig.2a,b),while the soluble sugar content was lower under the MD treatment than under the CK treatment (Fig. 2c), suggesting a higher efficiency of grain filling induced by MD treatment. Consistently,the contents of NSCs, starch, and soluble sugars in superior spikelets showed no difference between the CK and MD conditions(Fig. S3).

    Fig.1. Moderate soil drying at 6 DAA promoted grain filling and increased grain weight of inferior kernels.(a)Dry weight of inferior kernels exposed to control condition(CK)and moderate soil drying (MD). (b) Grain yield per plant exposed to CK and MD. (c, d) Dynamic grain weight and grain filling rate of inferior spikelets under CK and MD.Vertical bars indicate mean±standard error(n=12 plants in panels a,b;n=3 in panels c,d).Arrows in panels(a)and(b)means improvement of grain weight and grain yield.Post-hoc comparisons were tested by Duncan’s test, * and **, significant at the P <0.05 and P <0.01 probability levels, respectively.

    3.2. Analysis of small-RNA libraries

    small-RNA analysis of samples at the most active grain filling stage (9 DAA) in inferior kernels showed that, respectively,11,335,002 and 11,841,917 reads with error rates of 0.02%, Q30 values of 95.88%, and GC contents of 48.26% were obtained for CK and MD sRNA libraries. Removal of ambiguous reads left 10,862,401 (95.83%) and 11,403,353 (96.30%) high-quality reads.Of these, 5,703,492 (62.2%) and 5,546,422 (56.4%) were mapped to the rice genome version 7 of the MSU Rice Genome Annotation Project (RGAP) (Table S3). We categorized all mapped reads according to their properties, including the length and alignment position on the genome, using miRbase and Rfam among the sequences of CK and MD treatments. Of these miRNAs, 171,443(3.01%) and 225,992 (4.07%) were similar in sequence to known miRNAs. Other sequences, including rRNAs (CK: 957,111, 16.78%;MD: 585,962, 10.56%), tRNAs (109,996, 1.93%; 59,765, 1.08%),snoRNAs (24,546, 0.43%; 21736, 0.39%), snRNAs (5060, 0.09%;2336, 0.04%), exons (296,122, 5.19%; 374,262, 6.75%), introns(311,453, 5.46%; 324,504, 5.85%) and unknown sequences(3,806,794, 66.74%; 3,934,036, 70.93%) were also identified(Table S3). Analyses of the unique sRNA sequence revealed that 21-nt (CK, 9.85%; MD, 11.26%) and 24-nt (CK, 40.35%; MD,41.68%) small RNAs were predominant (Fig. S3).

    3.3. Differentially expressed miRNAs and their targets in inferior kernels in response to MD treatment

    Respectively 649 (521 known and 128 novel) and 635 (514 known and 121 novel) miRNAs were obtained from the CK and MD treatment libraries (Tables S4, S5; Fig. S5). A total of 47 (38 known and 9 novel) miRNAs differentially between the CK and MD treatments in inferior kernels were obtained (Fig. 3).

    Thirty-eight known miRNAs were differentially expressed,with 21 miRNAs downregulated (Fig. 3a) and 17 upregulated (Fig. 3b).Some differentially expressed miRNAs between inferior spikelets and superior spikelets were downregulated or upregulated by MD treatment. miR156, miR159, miR167, miR397 and miR1432,the most abundant miRNAs that were more highly expressed in inferior than in superior spikelets [33], were suppressed by MD treatment (Fig. 3a). In contrast, expression levels of miR1861b,e,f,g,i,k,l,m,miR1863c,miR1866,and miR1871,which displayed lower expression levels in inferior than in superior spikelets,were inhibited by MD treatment, especially miR1866, which was highly expressed in inferior spikelets(Fig. 3b). Nine novel miRNAs differentially expressed between CK and MD treatments in inferior kernels were obtained (Fig. 3c). The miRNAs differentially expressed in response to MD treatment indicate the potential role of miRNAs in regulating MD-promoted grain filling of inferior spikelets.

    Fig. 2. Moderate soil drying alterede NSC concentration in inferior kernels. NSC content (a), starch accumulation (b), and soluble sugar content (c) in inferior spikelets subjected to MD treatment during grain filling. Vertical bars indicate standard errors (n =3).Post-hoc comparisons were tested by Duncan’s test. *and **, significant at the P <0.05 and P <0.01 probability levels, respectively.

    Fig.3. Identification of differentially expressed miRNAs and most enriched KEGG pathways of miRNA targets in inferior kernels under CK and MD treatments. Heat map of differentially expressed known miRNAs down-regulated (a) or up-regulated (b) by MD treatment. (c) Heat map of differently expressed novel miRNAs of inferior spikelets from Nip at 9 DAA under CK and MD treatment.The figure was plotted using the normalized values of the expression of each gene in different samples.‘‘0”(blue color)refers to the low value and the‘‘4”(red color)represents the high value.(d)KEGG annotation analysis of differentially expressed miRNA targets of inferior spikelets of Nip at 9 DAA under CK and MD treatment.

    Plant miRNAs recognize their target mRNAs through perfect or near perfect base pairing, and function by cleaving targets by translational repression [41]. For this reason, identification of miRNA target genes is a prerequisite for the prediction of miRNA function [42]. To understand the functions of identified miRNAs,4786 putative targets for all identified miRNAs from both treatments were predicted by psRobot, targetfinder, and RNAhybrid,including 251 targets predicted by the differentially expressed miRNAs (Table S6), including auxin-responsive factor (ARF) genes,MADS-box family gene (MADS) genes, NAC domain-containing protein genes, MYB family transcription factor genes, and other transcription factors or genes.Gene ontology(GO)and Kyoto Encyclopedia of Genes and Genomes (KEGG) annotation and enrichment analysis (Fig. 3d, S6 and S7) were performed to identify the potential functions of these miRNA target genes. Of these, several miRNA target genes were enriched in functions of ‘‘binding”, ‘‘catalytic activity”, ‘‘cell part”, ‘‘cellular process”, ‘‘membrane part”,‘‘metabolic process”, ‘‘organelle” and ‘‘biological regulation”(Fig. S6). KEGG analysis of the miRNA target genes revealed that most of the target genes were enriched in functions of‘‘metabolism”, ‘‘genetic information processing” and ‘‘cellular processes” (Fig. 3d). ‘‘Transport and catabolism” and ‘‘carbohydrate metabolism” were also represented in the most enriched KEGG pathways, suggesting that the differentially expressed miRNAs induced by MD treatment might regulate grain filling by regulating carbohydrate transport and metabolism via their target genes.

    Fig. 4. Moderate soil drying increased starch biosynthesis via the miR1861-OsSBDCP1-OsSSIIIa signal pathway of inferior kernels. (a-d) qRT-PCR analysis of expression patterns of osa-miR1861,OsSBDCP1,and OsSSIIIa in 9-DAA inferior spikelets.(e)Enzyme activities of StSase in 9-DAA inferior spikelets.(f)Rice osa-miR1861 cleavage site in OsSBDCP1 mRNA. Vertical bars indicate mean ± standard error (n = 3). Post-hoc comparisons were tested by Duncan’s test. * and **, significant at the P <0.05 and P <0.01 probability levels, respectively.

    3.4. Upregulation of miR1861 led to inhibition of the OsSBDCP1 gene,in turn releasing SSIIIa for starch accumulation

    MiR1861 was induced by MD treatments in inferior kernels in small-RNA-seq results (Fig. 3b), a finding further confirmed by qRT-PCR. Both miR1861e,k,m and miR1861b,f,i,l were induced by MD treatment (Fig. 4a, b).OsSBDCP1(LOC_Os01g63810), encoding a starch binding domain-containing protein, was predicted to be the target gene of miR1861 (Table S6), which contained cleavage sites for both miR1861e,k,m and miR1861b,f,i,l (Fig. 4f). As expected, the expression of theOsSBDCP1gene was reduced in MD-treated inferior spikelets(Fig.4c),displaying a reverse expression pattern from that of miR1861 (Fig. 4a, b). Recently, it was reported [16] thatOsSBDCP1acted as a repressor of SSIIIa activity by both directly binding to OsSSIIIa and reducing the abundance of SSIIIa protein. In the present study, both expression ofOsSSIIIagene and SSIIIa activity were significantly increased by MD treatment in inferior spikelets(Fig.4d,e),suggesting that MD treatment might suppress the expression ofOsSBDCP1by upregulating the expression of miR1861,thus increasing the abundance and activity of the SSIIIa enzyme and eventually the accumulation of starch in inferior spikelets.

    3.5. The miRNA-phytohormone pathway was also involved in MDpromoted grain filling of inferior spikelets

    MiR1432-5p was one of the most significantly suppressed miRNAs in inferior spikelet by MD treatment(Fig.3a),as confirmed by qRT-PCR (Fig. 5a). miR1432 was reported [34] to negatively regulate grain filling of rice by cleaving its target gene,OsACOT(LOC_Os04g35590), a regulator of phytohormone biosynthesis.BecauseOsACOTwas also predicted to be a target gene of miR1432 (Table S6; Fig. 5f), we examined the expression of theOsACOTgene. It was significantly induced in MD-treated inferior spikelets (Fig. 5b), resulting in elevated levels of both ABA and IAA in inferior kernels under MD treatment (Fig. 5c, d). As auxin response genes,OsGH3-4andOsGH3-8showed increased expression in MD-treated inferior spikelets (Fig. 5e), in agreement with the finding that the content of IAA was increased by MD treatment.Both ABA and IAA are necessary for grain filling[34],and MD treatment might increase the expression ofOsACOTby reducing the abundance of miR1432, which would increase the contents of ABA and IAA and thus promote grain filling.

    3.6. Moderate soil drying increased inferior grain weight via the miR397-OsLAC signaling pathway

    As reported previously [20], miR397 increases grain yield by downregulating its target,OsLAC. miR397 was one of the most abundant miRNAs in inferior spikelets,and small-RNA-seq analysis showed that miR397 was induced by MD treatment in inferior spikelets (Fig. 6a), a finding further confirmed by qRT-PCR. The expression of both miR397a and miR397b was significantly increased by MD (Fig. 6b, c). As shown in Fig. 6e, both miR397a and miR397b can recognize and cleave theOsLACgene which plays a negative role in grain filling. The finding that the expression ofOsLACwas significantly reduced in inferior spikelets under MD treatment (Fig. 6d) suggested that MD treatment might suppress the expression ofOsLACvia upregulation of miR397, eventually leading to increased grain filling in inferior spikelets.

    4. Discussion

    4.1. Moderate soil drying after anthesis promotes grain filling of inferior spikelets in rice

    Rice is one of the most important staple foods in the world.The poor grain filling efficiency of inferior spikelets is becoming a severe problem in rice production, especially in some super rice varieties with huge panicles[6,37].In this study we demonstrated that MD treatment increased the grain filling rate and grain weight of inferior spikelets in Nip (Fig. 1). Consistently, the contents of non-structural carbohydrates and starch in inferior spikelets were also increased by MD, while the soluble sugar content was much lower in inferior spikelets under MD than under CK (Fig. 2), suggesting that MD treatment promoted the conversion of soluble sugar to starch. In agreement,moderate soil drying irrigation during the mid-to-late grain filling stages can increase the grain filling rate and grain weight of inferior spikelets [9,43].

    Fig. 5. Moderate soil drying increased ABA and IAA concentration via the miR1432-OsACOT signal pathway to promote inferior grain filling. (a-b) qRT-PCR analysis of expression patterns of osa-miR1432-5p and OsACOT in 9 DAA inferior spikelets. (c-d) ABA and IAA content in inferior spikelet of rice cultivar Nip under MD treatment. (e)qRT-PCR analysis of IAA auxin-responsive genes in inferior spikelets under MD treatment at 9 DAA. (f) Rice osa-miR1432-5p cleavage site in OsACOT mRNA. Vertical bars indicate mean ± standard error (n = 3). Post-hoc comparisons were tested by Duncan’s test. * and **, significant at the P <0.05 and P <0.01 probability levels, respectively.

    Fig.6. MiR397-OsLAC,a signal pathway for increainge inferior grain weight under moderate soil drying.(a)Small-RNA sequencing result of osa-miR397a and osa-miR397b.(b-d)qRT-PCR analysis of expression patterns of osa-miR397a,osa-miR397b,and OsLAC in 9 DAA inferior spikelet.(e)Rice osa-miR397 cleavage site in OsLAC mRNA.Vertical bars indicate mean±standard error(n=3).Post-hoc comparisons were tested by Duncan’s test.*and**,significant at the P <0.05 and P <0.01 probability levels,respectively.

    4.2. Differentially expressed miRNAs in response to MD treatment contribute to increased grain filling in inferior spikelets

    MiRNAs are epigenetic regulators of plant growth and development, as well as response to stress [44]. The regulation of seed development by miRNAs represents a research focus of miRNA function study [45]. A growing body of evidence [33,34,46,47]has demonstrated that miRNAs also play crucial roles in controlling grain filling.Mutation of miR156,which was highly expressed during grain filling [33] increased 1000-kernel weight in comparison with WT seeds [47]. Expression of miR156 was higher in inferior than in superior spikelets, especially in the mid-to-late stage of grain filling [48]. Both miR159 and miR167 were reported [33,49]to show higher expression in inferior than in superior spikelets and to regulate seed development via auxin and BRs, respectively. In the present study, 38 known miRNAs and 9 novel miRNAs were differentially expressed between CK and MD treatments in inferior spikelets(Fig.3).Among them,miR156,miR159,and miR167 were reduced by MD treatment in inferior kernels (Fig. 3a), suggesting that these miRNAs interfere with grain filling such that their repression by MD treatment would lead to increased grain filling of inferior spikelets.

    In contrast,the expression of miR1861 and miR1866 was moderately high and was induced by MD treatment in inferior spikelets(Fig.3b).It has been reported[50]that both miR1861 and miR1866 are regulated by salt stress. Overexpression of miR1861 increased salt tolerance in rice [51]. Both miR1861 and miR1866 were also reported [28,46,48] to be preferentially expressed in rice seeds and involved in regulating grain filling. Target gene analysis revealed two candidate targets of miR1866,LOC_Os02g13870(aquaporin protein gene,putative)andLOC_Os06g44060(phospholipase D gene,OsPLD)(Table S6),which were also reported[52-54]to be involved in regulating yield formation and other developmental processes. These reports have demonstrated that several miRNAs, together with their target genes, function in controlling grain filling. The finding that some of these seed developmentrelated miRNAs were induced or suppressed by MD treatment suggests that miRNAs differentially expressed in response to MD might contribute to improved grain filling in inferior spikelets.

    4.3.MD-regulated miRNAs participate in starch biosynthesis in inferior spikelets

    Starch in endosperm makes up approximately 90% of the final dry weight of brown rice [55], making its biosynthesis the most important process during grain filling. Biosynthesis of starch in kernels is catalyzed by numerous enzymes, of which starch synthase (StS) has been shown [16] to be a key regulator of starch accumulation.SSIIIa is the major form expressed in rice endosperm[56], and its abundance and activity are regulated by SBDCP1 protein, a starch binding domain-containing protein which is a negative regulator in starch biosynthesis during grain filling.Expression ofOsSBDCP1is suppressed by artificial miRNA,resulting in restored protein levels and enzyme activity of SSIIIa [16]. These results suggest the possibility of miRNA-regulated starch biosynthesis via the SBDCP1-SSIIIa complex during grain filling.OsSBDCP1was predicted to be the target gene for miR1861, which was induced by MD in inferior spikelets(Table S6;Fig.4b).The expression of miR1861 is much higher in superior than in inferior kernels,suggesting that miR1861 might be involved in regulating grainfilling differences between superior kernels and inferior kernels[33,48]. The inverse expression pattern of miR1861 andOsSBDCP1between CK and MD treatment in our results (Fig. 4a-c) further supports the notion thatOsSBDCP1is regulated by miR1861 during grain filling.Upregulation of miR1861 eventually led to an increase inOsSSIIIagene expression and StS enzyme activity (Fig. 4d, e).Thus,MD promotion of starch biosynthesis was mediated partially by the miR1861-OsSBDCP1-OsSSIIIapathway,resulting in increased grain filling in inferior spikelets.

    4.4. MiRNA-regulated plant hormones in response to MD treatment also contributed to increased grain filling in inferior spikelets

    Numerous studies have shown an interaction between epigenetic modifications and plant hormones in regulating plant growth and development,including modulation of plant hormone biosynthesis, transport, and signal transduction by miRNAs [57]. Plant hormones have been shown [18] to be crucial in regulating grain filling, and the accumulation of ABA, IAA, and BRs in rice grain was found[19,22,58]to be important for endosperm development.However, the relationship between miRNAs and plant hormones during grain filling has yet to be established. Recently [34],miR1432 was shown to negatively regulate grain filling via its target gene,OsACOT,which encodes an acyl-CoA thioesterase.Further study revealed that the miR1432-OsACOTmodule regulated fatty acid metabolism and phytohormone biosynthesis. Both ABA and IAA were increased in the miR1432 mutant andOsACOToverexpression lines. The expression of miR1432 is regulated by drought stress [59]. miR1432 displayed higher abundance in inferior than in superior spikelets [33], suggesting a relationship of miR1432 with drought stress and grain filling. In the present study, moderate soil drying reduced the expression of miR1432 in inferior kernels(Figs.3a,5a),a response accompanied by upregulation ofOsACOT(Fig. 5b). The contents of both ABA and IAA were also elevated in MD-treated inferior kernels (Fig. 5c, d). Our results suggest that MD treatment might upregulate the expression ofOsACOTby suppressing miR1432,resulting in elevated ABA and IAA and thus promoting grain filling in inferior spikelets.

    Fig.7. Signal pathways of miRNA involved in inferior grain filling under moderate soil drying.Three miRNA-target pathways were involved in MD-promoted grain filling.The first one is miR1861-OsSBDCP1-OsSSIIIa pathway,in which miR1861 was up-regulated by MD treatment,resulting in the downregulation of OsSBDCP1 and thus the activation of SSIIIa enzyme. In the miR397-OsLAC module, MD treatment also increased the expression of miR397, which cleaved the OsLAC gene and thereby increased BR signal. In contrast to miR1861 and miR397,miR1432 was reduced by MD treatment,resulting in up-regulation of OsACOT,leading to increased ABA and IAA contents.Increased starch synthesis, as well as elevated plant hormone biosynthesis and signal transduction by MD-regulated miRNAs, increase grain filling in inferior spikelets. Red color indicates increased content or expression and blue color indicates reduced expression of genes.Arrows represent positive effects;lines ending in cross bars represent negative effects.Dotted line indicates an unknown pathway.

    Besides miR1432, expression of miR397 was also changed by MD in inferior spikelets. Expression of miR397 was promoted in inferior kernels under MD treatment(Fig.6a,b).MiR397 is preferentially expressed in rice kernels [46] and superior spikelets have higher levels of miR397 than inferior spikelets [33], implying that miR397 is involved in regulating grain filling differences between superior and inferior spikelets. Furthermore, upregulation of miR397 by overexpression or in theOsAGO17-overexpressing line increased grain size and weight by reducing the expression ofOsLACgenes and cleavage targets for miR397 [20,60].OsLACnegatively regulated BR signaling and suppression ofOsLACincreased grain size and weight in rice [20,61,62]. Consistently, expression ofOsLACwas reduced by MD treatment in inferior spikelets(Fig.6c),suggesting the involvement of the miR397-OsLACmodule in mediating MD-promoted grain filling. Expression ofOsLACwas also reduced by MD treatment in inferior spikelets (Fig. S8). Thus,MD treatment might activate the miR397-OsLACmodule to promote grain filling via the BR signaling pathway in inferior kernels.

    5. Conclusions

    Post-anthesis moderate soil drying differentially affected the expression of miRNAs many of which have been reported to function in regulating grain size and grain weight. Target genes of differentially expressed miRNAs were involved in the carbohydrate metabolism pathway. The three most important miRNA-target pathways were involved in MD-promoted grain filling. The first is the miR1861-OsSBDCP1-OsSSIIIapathway, in which miR1861 was upregulated by MD treatment, resulting in downregulation ofOsSBDCP1and thus the activation of SSIIIa enzyme. In the miR397-OsLACmodule, MD treatment also increased the expression of miR397, which cleaved theOsLACgene, thereby increasing the BR signal. In contrast to miR1861 and miR397, miR1432 was reduced by MD treatment,resulting in upregulation of theOsACOTgene,leading to increased ABA and IAA contents(Fig.7).Thus,promotion of starch synthesis, as well as elevated plant hormone biosynthesis and signal transduction by MD-regulated miRNAs,increase grain filling in inferior spikelets.

    CRediT authorship contribution statement

    Zhenning Teng:Methodology,Investigation,Writing-original draft, Writing - review & editing.Yinke Chen:Investigation,Methodology, Formal analysis.Youqing Yuan:Investigation, Formal analysis.Yaqiong Peng:Investigation.Yake Yi:Investigation.

    Huihui Yu:Investigation.Zhenxie Yi:.Jianchang Yang:Writing- review & editing.Yan Peng:Methodology.Meijuan Duan:Conceptualization, Funding acquisition.Jianhua Zhang:Conceptualization, Funding acquisition, Writing - review & editing.Nenghui Ye:Project administration, Visualization, Funding acquisition,Writing - review & editing.

    Declaration of competing interest

    The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.

    Acknowledgments

    This work was supported by the National Natural Science Foundation of China (32171927, 31971924, U21A201755), the Natural Science Foundation of Hunan Province (2021JJ30349), the Key Research and Development Program of Hunan Province(2018NK1010), Science and Technology Plan of Changsha City(kq2004034), Scientific Research Project of Education Department of Hunan Province (19A245), and the Hong Kong Research Grant Council (AoE/M-05/12, AoE/M-403/16, GRF12103219, 12103220,14177617).

    Appendix A. Supplementary data

    Supplementary data for this article can be found online at https://doi.org/10.1016/j.cj.2021.11.004.

    av播播在线观看一区| 春色校园在线视频观看| 中国三级夫妇交换| 亚洲精品视频女| 精品熟女少妇av免费看| 亚洲无线观看免费| 成人综合一区亚洲| 新久久久久国产一级毛片| 蜜桃国产av成人99| 国语对白做爰xxxⅹ性视频网站| 高清不卡的av网站| 男女边吃奶边做爰视频| 七月丁香在线播放| 欧美日韩在线观看h| 天堂中文最新版在线下载| 女人久久www免费人成看片| 国产精品蜜桃在线观看| 一区二区三区乱码不卡18| 国产精品一国产av| 午夜福利网站1000一区二区三区| 午夜福利,免费看| 另类亚洲欧美激情| 中文字幕人妻熟人妻熟丝袜美| 特大巨黑吊av在线直播| 亚洲精品第二区| 久久久久网色| 一区二区三区免费毛片| 黄片播放在线免费| 久久人人爽人人爽人人片va| 亚洲国产精品专区欧美| 观看av在线不卡| 22中文网久久字幕| 午夜福利在线观看免费完整高清在| 免费高清在线观看日韩| 亚洲精品美女久久av网站| 欧美激情 高清一区二区三区| 两个人免费观看高清视频| 在线观看人妻少妇| 高清毛片免费看| 亚洲精品国产av蜜桃| 亚洲精品中文字幕在线视频| 少妇熟女欧美另类| 下体分泌物呈黄色| 老司机影院成人| 一区二区三区免费毛片| 亚洲,一卡二卡三卡| 欧美日韩视频高清一区二区三区二| 大片电影免费在线观看免费| 亚洲精品色激情综合| 国产伦精品一区二区三区视频9| 亚洲天堂av无毛| 欧美日韩视频精品一区| av在线观看视频网站免费| 纵有疾风起免费观看全集完整版| 寂寞人妻少妇视频99o| 人体艺术视频欧美日本| 人人妻人人澡人人爽人人夜夜| 亚洲久久久国产精品| 麻豆成人av视频| 18+在线观看网站| 99久久精品一区二区三区| 成人综合一区亚洲| 久久狼人影院| 欧美一级a爱片免费观看看| 亚洲第一av免费看| 成人综合一区亚洲| 精品视频人人做人人爽| 日本欧美视频一区| 精品亚洲成国产av| 欧美精品高潮呻吟av久久| 亚洲av成人精品一区久久| a 毛片基地| 男男h啪啪无遮挡| av国产精品久久久久影院| 国产黄片视频在线免费观看| 国产日韩一区二区三区精品不卡 | 成人影院久久| 中文字幕最新亚洲高清| 中文字幕最新亚洲高清| 在线观看国产h片| 激情五月婷婷亚洲| 欧美国产精品一级二级三级| 9色porny在线观看| 精品午夜福利在线看| 亚洲av.av天堂| 3wmmmm亚洲av在线观看| 久久久久精品久久久久真实原创| 国产欧美日韩综合在线一区二区| 99热国产这里只有精品6| 日本午夜av视频| 69精品国产乱码久久久| 日韩精品有码人妻一区| 国产一区二区在线观看av| 亚洲综合色惰| 各种免费的搞黄视频| 久久久国产一区二区| 日韩av免费高清视频| 国产男人的电影天堂91| 欧美精品国产亚洲| 亚洲av福利一区| 国产精品免费大片| 色哟哟·www| 有码 亚洲区| 国产国拍精品亚洲av在线观看| 久久久久视频综合| 欧美精品国产亚洲| 精品久久国产蜜桃| 国国产精品蜜臀av免费| 制服人妻中文乱码| 天堂8中文在线网| 天堂8中文在线网| 国产黄片视频在线免费观看| 亚洲欧美日韩卡通动漫| 男男h啪啪无遮挡| 一本大道久久a久久精品| 老司机亚洲免费影院| 最近2019中文字幕mv第一页| 99热6这里只有精品| 国产亚洲一区二区精品| 国产精品蜜桃在线观看| 久久精品国产鲁丝片午夜精品| 久久久久久久久久久丰满| 日本91视频免费播放| 国产免费又黄又爽又色| 久热这里只有精品99| 午夜免费男女啪啪视频观看| 午夜福利,免费看| 国产成人freesex在线| 看免费成人av毛片| 久久女婷五月综合色啪小说| av.在线天堂| 亚洲精品久久久久久婷婷小说| 3wmmmm亚洲av在线观看| 亚洲av.av天堂| 国产男女内射视频| 天堂8中文在线网| 18禁在线无遮挡免费观看视频| 一区二区三区乱码不卡18| 国内精品宾馆在线| 免费黄色在线免费观看| 中文字幕亚洲精品专区| 极品少妇高潮喷水抽搐| 日韩视频在线欧美| kizo精华| 九草在线视频观看| 欧美三级亚洲精品| 亚洲精品乱码久久久久久按摩| av在线观看视频网站免费| 日本午夜av视频| 亚洲一级一片aⅴ在线观看| 久久午夜综合久久蜜桃| 美女福利国产在线| 高清欧美精品videossex| 久久久久精品性色| 91精品国产九色| 老女人水多毛片| 99久久精品一区二区三区| 色5月婷婷丁香| 国产白丝娇喘喷水9色精品| av.在线天堂| 免费高清在线观看视频在线观看| 国产熟女午夜一区二区三区 | 热re99久久精品国产66热6| 午夜精品国产一区二区电影| 精品亚洲成a人片在线观看| 男人添女人高潮全过程视频| 亚洲欧美成人精品一区二区| 看免费成人av毛片| 两个人的视频大全免费| 哪个播放器可以免费观看大片| 九九久久精品国产亚洲av麻豆| 亚洲av福利一区| 99九九线精品视频在线观看视频| 91久久精品国产一区二区成人| 亚洲综合精品二区| 亚洲人与动物交配视频| 亚洲av不卡在线观看| 国产亚洲av片在线观看秒播厂| 亚洲第一av免费看| 亚洲性久久影院| 亚洲第一区二区三区不卡| 久久久午夜欧美精品| 少妇丰满av| 有码 亚洲区| 在线观看免费高清a一片| 大香蕉97超碰在线| 日韩一区二区三区影片| av线在线观看网站| 黄色毛片三级朝国网站| 久久国内精品自在自线图片| 青春草国产在线视频| 精品亚洲乱码少妇综合久久| 青春草视频在线免费观看| 久久亚洲国产成人精品v| 免费高清在线观看日韩| 亚洲精品日韩在线中文字幕| 国产一区二区在线观看av| 日韩强制内射视频| 最近手机中文字幕大全| 青春草国产在线视频| 久久久a久久爽久久v久久| 午夜福利,免费看| 99视频精品全部免费 在线| 亚洲精品视频女| 久久久久久久久大av| 精品一区在线观看国产| 久久午夜福利片| 一边摸一边做爽爽视频免费| 乱人伦中国视频| 久久久欧美国产精品| 国产免费福利视频在线观看| 欧美bdsm另类| 多毛熟女@视频| 91精品国产国语对白视频| 国产精品99久久99久久久不卡 | 亚洲精品第二区| 欧美成人午夜免费资源| 赤兔流量卡办理| 国产免费现黄频在线看| 一区二区三区精品91| 美女国产视频在线观看| 国产精品女同一区二区软件| 大香蕉97超碰在线| 亚洲精品国产av蜜桃| 80岁老熟妇乱子伦牲交| 亚洲成色77777| 成人手机av| xxx大片免费视频| av国产精品久久久久影院| 欧美xxxx性猛交bbbb| 国产深夜福利视频在线观看| 91精品一卡2卡3卡4卡| 全区人妻精品视频| 欧美日韩精品成人综合77777| 又粗又硬又长又爽又黄的视频| av在线播放精品| 久久久久国产精品人妻一区二区| 日日啪夜夜爽| 亚洲无线观看免费| 高清黄色对白视频在线免费看| 少妇丰满av| 久久影院123| 一区二区三区乱码不卡18| 如何舔出高潮| 在线观看免费视频网站a站| 美女脱内裤让男人舔精品视频| 国产极品天堂在线| 中文精品一卡2卡3卡4更新| 久久精品久久精品一区二区三区| 欧美xxxx性猛交bbbb| 精品亚洲成国产av| 成人免费观看视频高清| 超色免费av| 国产国拍精品亚洲av在线观看| 卡戴珊不雅视频在线播放| 亚洲精品成人av观看孕妇| 一个人免费看片子| 成人免费观看视频高清| 汤姆久久久久久久影院中文字幕| 777米奇影视久久| 卡戴珊不雅视频在线播放| 国产乱来视频区| 亚洲av成人精品一区久久| 国产成人精品在线电影| 91aial.com中文字幕在线观看| 色哟哟·www| 2021少妇久久久久久久久久久| 99热这里只有精品一区| 亚洲久久久国产精品| av网站免费在线观看视频| 黑人高潮一二区| 国产成人a∨麻豆精品| 免费看av在线观看网站| 另类亚洲欧美激情| 三级国产精品欧美在线观看| 亚洲综合色惰| 日日摸夜夜添夜夜爱| 国产精品国产三级国产专区5o| 两个人免费观看高清视频| 日本-黄色视频高清免费观看| 少妇人妻 视频| 国产高清有码在线观看视频| 一边摸一边做爽爽视频免费| 精品一区二区三卡| 中文字幕人妻丝袜制服| 亚洲av成人精品一二三区| 青青草视频在线视频观看| 久久精品夜色国产| 午夜免费观看性视频| 中文字幕人妻丝袜制服| 亚洲伊人久久精品综合| 一区二区三区四区激情视频| 秋霞在线观看毛片| 国产av一区二区精品久久| 纵有疾风起免费观看全集完整版| 亚洲一级一片aⅴ在线观看| 亚洲成人手机| 特大巨黑吊av在线直播| 欧美日本中文国产一区发布| 91aial.com中文字幕在线观看| 欧美少妇被猛烈插入视频| 日韩亚洲欧美综合| 欧美精品高潮呻吟av久久| 国产熟女欧美一区二区| 久久久久久伊人网av| 在线观看免费高清a一片| www.av在线官网国产| 亚洲一区二区三区欧美精品| 大片电影免费在线观看免费| 午夜老司机福利剧场| 亚洲精品乱码久久久v下载方式| www.av在线官网国产| 欧美日韩在线观看h| 啦啦啦中文免费视频观看日本| 成人二区视频| 成年人免费黄色播放视频| 亚洲熟女精品中文字幕| 视频中文字幕在线观看| 亚洲精品中文字幕在线视频| 日韩成人伦理影院| 只有这里有精品99| 欧美性感艳星| 建设人人有责人人尽责人人享有的| 欧美变态另类bdsm刘玥| 亚洲精品乱码久久久久久按摩| 激情五月婷婷亚洲| 亚洲精品第二区| 尾随美女入室| 亚洲美女黄色视频免费看| 日韩成人av中文字幕在线观看| 九色亚洲精品在线播放| 美女中出高潮动态图| 午夜福利网站1000一区二区三区| 国产国拍精品亚洲av在线观看| 美女国产高潮福利片在线看| 日韩视频在线欧美| a级毛片在线看网站| 久久ye,这里只有精品| 久久99热6这里只有精品| kizo精华| 黄片无遮挡物在线观看| 精品少妇久久久久久888优播| 亚洲国产成人一精品久久久| 国产在线免费精品| www.色视频.com| 国模一区二区三区四区视频| 国产日韩欧美亚洲二区| 在线看a的网站| 精品久久国产蜜桃| 一级黄片播放器| 午夜福利,免费看| 哪个播放器可以免费观看大片| 亚洲欧美日韩卡通动漫| 精品99又大又爽又粗少妇毛片| 成人免费观看视频高清| 我的女老师完整版在线观看| 国产精品久久久久成人av| 国产免费又黄又爽又色| 美女大奶头黄色视频| 国产亚洲最大av| 超碰97精品在线观看| 国产成人精品福利久久| 国产精品秋霞免费鲁丝片| 成人影院久久| 成人无遮挡网站| 极品少妇高潮喷水抽搐| 全区人妻精品视频| 黑人高潮一二区| 免费大片18禁| 国产无遮挡羞羞视频在线观看| 国产高清三级在线| 伦理电影大哥的女人| 日韩制服骚丝袜av| 人妻人人澡人人爽人人| 国产淫语在线视频| 视频在线观看一区二区三区| 黄片播放在线免费| 午夜精品国产一区二区电影| 亚洲成人av在线免费| 色吧在线观看| 亚洲精品日本国产第一区| 国产成人午夜福利电影在线观看| 一区在线观看完整版| 精品少妇内射三级| 欧美精品高潮呻吟av久久| 亚洲精品aⅴ在线观看| 男女国产视频网站| 欧美丝袜亚洲另类| 免费高清在线观看日韩| 高清av免费在线| 国产日韩欧美视频二区| 欧美+日韩+精品| 日本午夜av视频| 美女内射精品一级片tv| 国产一区有黄有色的免费视频| 在线看a的网站| 久久午夜综合久久蜜桃| 亚洲国产成人一精品久久久| 精品国产露脸久久av麻豆| 九色成人免费人妻av| 热re99久久精品国产66热6| 久久久久久久久久久免费av| 我的女老师完整版在线观看| 国产成人aa在线观看| 人人妻人人爽人人添夜夜欢视频| 日韩不卡一区二区三区视频在线| 国产精品熟女久久久久浪| 久久亚洲国产成人精品v| 国产一区二区在线观看日韩| 人妻人人澡人人爽人人| 日韩欧美一区视频在线观看| 欧美精品一区二区免费开放| 久久韩国三级中文字幕| 久久久国产欧美日韩av| 三级国产精品欧美在线观看| 黄片播放在线免费| 中文字幕人妻熟人妻熟丝袜美| 高清黄色对白视频在线免费看| 中国国产av一级| 老熟女久久久| 9色porny在线观看| 蜜桃久久精品国产亚洲av| 高清午夜精品一区二区三区| 日本猛色少妇xxxxx猛交久久| 久久久精品区二区三区| 女性生殖器流出的白浆| 大片电影免费在线观看免费| 建设人人有责人人尽责人人享有的| 人妻制服诱惑在线中文字幕| 免费不卡的大黄色大毛片视频在线观看| 九色亚洲精品在线播放| 黄色配什么色好看| 色婷婷av一区二区三区视频| 国产精品一区www在线观看| 亚洲久久久国产精品| 97超视频在线观看视频| 少妇 在线观看| 少妇高潮的动态图| xxx大片免费视频| 国产无遮挡羞羞视频在线观看| 色网站视频免费| 黑人猛操日本美女一级片| 欧美 日韩 精品 国产| 国产精品熟女久久久久浪| kizo精华| 日本vs欧美在线观看视频| 熟女人妻精品中文字幕| 久久鲁丝午夜福利片| 亚洲综合色网址| 婷婷成人精品国产| av.在线天堂| 美女xxoo啪啪120秒动态图| 精品亚洲成a人片在线观看| 欧美日韩成人在线一区二区| 大片免费播放器 马上看| 少妇的逼水好多| 国产一区二区在线观看日韩| av不卡在线播放| 亚洲国产精品999| 亚洲av.av天堂| 国产无遮挡羞羞视频在线观看| 国产有黄有色有爽视频| 成人18禁高潮啪啪吃奶动态图 | 国产探花极品一区二区| 日韩大片免费观看网站| 热re99久久精品国产66热6| 亚洲国产精品成人久久小说| 老女人水多毛片| 久久99蜜桃精品久久| 婷婷色综合www| 人妻一区二区av| 国产精品人妻久久久影院| 国产亚洲一区二区精品| 美女脱内裤让男人舔精品视频| 欧美性感艳星| 亚洲国产毛片av蜜桃av| 超色免费av| 久久国内精品自在自线图片| 国产在线视频一区二区| 国产成人精品婷婷| 99久久人妻综合| 18在线观看网站| 精品卡一卡二卡四卡免费| 黄色欧美视频在线观看| 欧美日韩一区二区视频在线观看视频在线| 成人影院久久| 精品亚洲成a人片在线观看| 制服诱惑二区| 日韩亚洲欧美综合| 黑丝袜美女国产一区| 国产成人av激情在线播放 | av有码第一页| 国产视频首页在线观看| 亚洲欧美色中文字幕在线| 日本色播在线视频| 久久精品夜色国产| 99精国产麻豆久久婷婷| 天美传媒精品一区二区| 久久狼人影院| 中国三级夫妇交换| 国产精品久久久久久精品古装| 伊人亚洲综合成人网| 亚洲高清免费不卡视频| 亚洲性久久影院| 免费av中文字幕在线| 乱码一卡2卡4卡精品| 五月伊人婷婷丁香| 国产成人aa在线观看| 国产有黄有色有爽视频| 免费观看无遮挡的男女| 18+在线观看网站| 国产国拍精品亚洲av在线观看| 熟妇人妻不卡中文字幕| 国产精品成人在线| 久久久久精品性色| 精品国产一区二区三区久久久樱花| 国产黄片视频在线免费观看| 久久婷婷青草| 精品久久蜜臀av无| 久久久久久人妻| av播播在线观看一区| 97精品久久久久久久久久精品| 精品国产露脸久久av麻豆| 啦啦啦视频在线资源免费观看| 18+在线观看网站| 亚洲综合色网址| 老女人水多毛片| 亚洲美女视频黄频| 嘟嘟电影网在线观看| 久久久亚洲精品成人影院| 日韩亚洲欧美综合| 欧美日韩综合久久久久久| 国产黄片视频在线免费观看| 日日摸夜夜添夜夜添av毛片| 国产精品欧美亚洲77777| 在线观看免费视频网站a站| 国产精品国产三级国产av玫瑰| 中文字幕制服av| 免费观看a级毛片全部| 热re99久久精品国产66热6| 国产精品国产三级专区第一集| 一本一本综合久久| 精品国产一区二区三区久久久樱花| av在线播放精品| 国产探花极品一区二区| 久久人人爽av亚洲精品天堂| 少妇的逼水好多| 国产精品三级大全| 大又大粗又爽又黄少妇毛片口| 色吧在线观看| 黑人巨大精品欧美一区二区蜜桃 | 国产黄频视频在线观看| 99热全是精品| 这个男人来自地球电影免费观看 | 亚洲欧洲日产国产| 涩涩av久久男人的天堂| 国产国语露脸激情在线看| 人妻系列 视频| 亚洲精品久久久久久婷婷小说| 狂野欧美白嫩少妇大欣赏| 永久免费av网站大全| 国产黄片视频在线免费观看| 人妻一区二区av| 18禁在线播放成人免费| 久久久欧美国产精品| 91久久精品国产一区二区三区| 日韩视频在线欧美| 日韩大片免费观看网站| 新久久久久国产一级毛片| 国产成人91sexporn| 精品亚洲乱码少妇综合久久| 国产欧美日韩综合在线一区二区| 久久精品久久久久久久性| 超色免费av| 久久久精品区二区三区| 国产精品成人在线| 欧美丝袜亚洲另类| av播播在线观看一区| 欧美丝袜亚洲另类| 国产午夜精品久久久久久一区二区三区| 日韩一区二区三区影片| 2022亚洲国产成人精品| 日韩一区二区三区影片| 超色免费av| 国产亚洲一区二区精品| 街头女战士在线观看网站| 妹子高潮喷水视频| 91aial.com中文字幕在线观看| 超碰97精品在线观看| 国产午夜精品久久久久久一区二区三区| 在线观看www视频免费| 日本av免费视频播放| 国产亚洲一区二区精品| 免费人妻精品一区二区三区视频| 晚上一个人看的免费电影| 国产 一区精品| 特大巨黑吊av在线直播| 少妇被粗大猛烈的视频| 99热这里只有精品一区| 国产伦精品一区二区三区视频9| 国产深夜福利视频在线观看| 亚洲国产av影院在线观看| 亚洲美女搞黄在线观看| 国产伦精品一区二区三区视频9| 免费黄网站久久成人精品| 少妇丰满av| 亚洲不卡免费看| 国产一区二区在线观看av| 精品人妻一区二区三区麻豆| av免费在线看不卡| 亚洲国产毛片av蜜桃av| 777米奇影视久久| 99久久精品一区二区三区| 亚洲欧洲精品一区二区精品久久久 | 久久99一区二区三区|