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

    Transcriptomes of early developing tassels under drought stress reveal differential expression of genes related to drought tolerance in maize

    2018-06-06 09:12:53WANGNanLlLiangGAOWenweiWUYongboYONGHongjunWENGJianfengLlMingshunZHANGDeguiHAOZhuanfangLlXinhai
    Journal of Integrative Agriculture 2018年6期

    WANG Nan, Ll Liang , GAO Wen-wei, WU Yong-bo YONG Hong-jun WENG Jian-feng Ll Mingshun ZHANG De-gui HAO Zhuan-fang , Ll Xin-hai

    1 Institute of Crop Sciences, Chinese Academy of Agricultural Sciences, Beijing 100081, P.R.China

    2 College of Agriculture, Xinjiang Agricultural University, Urumqi 830000, P.R.China

    1. lntroduction

    Drought is one of the major environmental stresses affecting maize production, and inadequate rainfall, due to climate change, will negatively impact agricultural productivity. The improvement of maize drought tolerance has long been important for enhancing the stability of food production,although conventional breeding most often takes advantage of favorable genetic recombination among germplasm to improve yield. The comprehensive elucidation of molecular regulatory mechanisms by which maize responds to drought stress is fundamentally important for maize genetics and breeding (Xionget al. 2002; Nakashimaet al. 2014).Abiotic stresses including drought, high temperature, and salinity all influence the expression of a large number of genes that induce or regulate other genes through complex transcriptional networks (Todakaet al. 2015).Transcriptome analyses of plant tissues and organs under drought-stressed conditions will reveal relevant genes and pathways and thereby helping unravel the mechanisms of drought tolerance.

    Before, some drought-inducible genes have been identified by analyzing the transcriptomes of maize reproductive organs using microarray technology (Zhuanget al. 2007). These genes encode cell stress-protective proteins and regulatory proteins, as well as signal transduction components (Shinozakiet al. 2003). Among techniques available for transcriptome analysis, microarrays require the design and synthesis of probes for known genes.This method has some technical limitations, including false positives generated by hybridization of multiple gene products to the arrays, which has limited their usefulness(Okoniewski and Miller 2006; Shenet al. 2014). Recently,along with the rapid and economic sequencing technology,sequencing method based on RNA level and its research results have systematically altered our view of the extent and complexity of eukaryotic transcriptomes (Wanget al.2009). A total of 619 maize drought response genes were identified using RNA-Seq, and 20 genes regulating flowering time were further analyzed, such asPRR37andCONSTANS(Leeet al. 2017). Drought-induced abortion of maize embryos suggested that signaling events such as increased abscisic acid (ABA) levels and decreased glucose levels interfered with ABA and sugar signaling pathways, which might have activated programmed cell death or senescence and arrested the cell cycle in the developing maize ovary(Kakumanuet al. 2012).

    Drought or heat severely affects normal pollination and decreases yield in maize. The period from meiosis to tetrad break-up in anthers and anthesis are key sensitivity periods to drought for maize. Drought stress during the first stage inhibits development of microspores or pollen, then causing male sterility (Saini 1997). Zhuanget al. (2007) identified 1 513 transcripts that were differentially expressed during meiosis in immature tassels and ears under drought-stressed and well-watered conditions, and 33 up-regulated transcripts were identified in both immature ears and tassels after drought stress using oligo microarray analysis. However, the metabolic pathways involved and the reasons that drought stress results in male sterility are still not clearly understood. In the present study, we used RNA-Seq to explore drought-induced gene expression in early developing tassels from 10 inbred lines,which selected from five heterotic groups of China (Liuet al.2015) with different drought resistance (Haoet al. 2011), under both well-watered and drought-stressed conditions. The goal of the present study was: i) to identify transcripts that were differentially expressed under drought-stress condition, ii) to explore the important gene networks and metabolic pathways induced by drought in early developing tassels, and iii) to propose transcriptional and other mechanisms of droughtinduced male sterility in maize.

    2. Materials and methods

    2.1. Plant materials and watering treatments

    The 10 inbred maize (Zea maysL.) lines were planted in 10-L pots containing a 2:2:1:1 (v/v) mix of field soil:peat:fowl manure:perlite with 1 g of fritted trace elements (Table 1).More than 20 parallel pots of each inbred line were planted for getting the same size of developing tassels after peeling off in the summer of 2014. Plants were irrigated by hand every other day to maintain soil water content close to capacity after the 3-leaf stage, and pesticides were applied every 2 weeks. When the plants reached to 9- or 10-leaf (V9 or V10)stage, some of them were subjected to drought-stressed condition, and some were still under well-watered condition(control). After treatment, when plants reached to floralorgan differentiation stage (V11 to V13 stages, varied among inbreds), early developing young tassels were sampled for RNA-Seq. At the same time, soil moisture content (%) was recorded to monitor the level of drought stress applied to each maize line. Steady-state leaf fluorescence was also recorded on leaves to reflect the corresponding response of plants under drought-stressed or well-watered conditions(Mini-PAM, Effeltrich, Germany). When the soil moisture content in drought-stressed treatments fell to about 25%,the sheath leaves of maize were stripped from plants, then the early developing tassels were sampled from all inbred lines in the experiment.

    Table 1 The information of 10 inbred lines experimented in this study1)

    2.2. lsolation of total RNA from tassel samples and RNA-Seq

    A total of 20 tassel samples from each of 10 inbred lines under drought-stressed and well-watered conditions were homogenized in liquid nitrogen. Total RNA was extracted from each sample using TRIzol?Reagent (Invitrogen,USA), according to the manufacturer’s protocol. As all the plants would be identified corresponding drought-tolerant responses under physiological dehydration, we took three inbreds at least as replications to select significantly differentially expressed genes for further analysis.

    Total RNA samples were sequenced at Beijing Genomics Institute, Shenzhen, China. The poly-A+RNA was sheared into short fragments of approximately 200 bp, and cleaved RNA fragments were copied into the first-strand of cDNA using random-hexamer primers and reverse transcriptase(Qiagen, Duesseldorf, Germany). The second-strand of cDNA was synthesized in reactions containing DNA polymerase I, RNaseH, buffer, and dNTPs (TransGen Biotech, Beijing, China). The cDNA fragments were purified using AMPure XP beads (Beckman Coulter, USA). Adaptorligated cDNAs were enriched by PCR in reactions to create the final cDNA library. Libraries were sequenced on the Illumina HiSeq 2000 Platform (Illumina, USA) according to the manufacturer’s recommendations.

    2.3. Data processing and analysis

    Raw reads were subjected to quality control, and the lowquality reads (low-quality bases in raw reads of >30%)and reads containing vector or adapter sequences were removed before data analysis. Clean reads acquired were aligned to the maize B73 reference genome (AGPv2) using SOAPaligner/SOAP2 (www.maizesequence.org) (Liet al.2009), only allowing mismatches of less than five bases between alignments. Digital gene expression counts were normalized by log2transformation and compared using the reads per kilobase transcriptome per million mapped reads method (RPKM) (Mortazaviet al. 2008). The specific algorithm developed by BGI for calculating gene expression was showed previously in Fuet al. (2013). Comparative analysis of differentially expressed genes (DEGs) from 10 inbred lines was performed using Microsoft Office 2013(Microsoft Corporation, US). Here, an false discovery rate(FDR)≤0.001 (Benjamini and Yekutieli 2001), an absolute value of the Log2Ratio≥1, andP≤0.001 were used as the thresholds by which to assess the significance of differences in gene expression.

    All DEGs were mapped to Gene Ontology (GO) terms under the top-level terms or domains at the GO database(http://www.geneontology.org), and transcript numbers were calculated for each term. Then a hypergeometric test with a Bonferroni correction for multiple comparisons was used to identify significantly enriched GO terms in the input list of DEGs using GO::Term Finder with a threshold ofP≤0.05 (Abdi 2007).Kyoto encyclopedia of genes and genomes (KEGG) was used to identify metabolic or signal transduction pathways that were significantly enriched between DEGs sampled from droughtstressed or well-watered plants relative to the transcriptomic background of early developing tassels (Kanehisaet al. 2008).

    2.4. Verification of transcript expression using qRT-PCR

    Seeds of the inbred line Qi 319 were germinated on moist filter paper in an incubator (28°C/22°C for 14 h/10 h). After four days, the germinated seedlings were then planted in quartz sand in a greenhouse. The seedling endosperms were peeled and transferred to Hoagland’s nutrient solution at 2-leaf (V2) stage. At 3-leaf (V3) stage, the seedlings were treated with 20% polyethylene glycol 6 000 (PEG-6000).

    Total RNA was extracted using Trizol (Transgen Biotech,Beijing, China) from leaves, stems, and roots of seedling after 0 (control), 1, 3, 6, and 12 h of PEG treatment. qRTPCR was performed using SuperReal PreMix Plus (SYBR Green) Kit (Tiangen) and products were quantified using a iQ5/MyiQ2 Real-time PCR Detection System (Bio-Rad,California, USA). Expression of the glyceraldehyde-3-phosphate dehydrogenase (GAPDH) gene from maize was used as an internal control for gene expression. Relative gene expression was calculated using the 2–ΔΔCTmethod(Livak and Schmittgen 2001).

    3. Results

    3.1. Phenotypic and photosynthetic responses to drought stress

    After five days without watering at floral organ differentiation stage of young tassel in the drought-stressed treatment, soil moisture contents in the pots decreased to an average of 24.08% with the highest of 26.80% in Chang 7-2 and the lowest of 22.00%in X178. In contrast, soil moisture contents in well-watered treatment remained constantly above 62.00%, with the highest of 64.30% in Qi 319 and the lowest of 62.00% in Dan 340, indicating that the control was in well-watered condition (Table 2).Because plant development is dynamic, the steadystate chlorophyll fluorescence of the uppermost leaf of each line was also monitored to determine the degree of photosynthetic inhibition due to drought stress. Average steady-state leaf fluorescence decreased by 0.127 under drought stress compared with that of control. As predicted, the tassel lengths sampled from all plants under drought-stressed condition were smaller than those under well-water condition (control). On average, the tassel length decreased from 4.11 to 3.01 cm after drought stress(Table 2 and Fig. 1). These results indicated that tassel growth was inhibited and gene expression was changed under drought stress.

    3.2. Maize transcriptome revealed gene expression changes under drough-stress condition

    A total of 940 million high-quality reads were generated from the transcriptomes of 20 immature tassels, which derived from 10 maize inbred lines under control and drought-stressed conditions.About 76.88% of the reads were mapped to the B73 reference genome (AGPv2) and 78.42%were mapped to the annotated maize genes,which revealed that there were still 23.12% unique sequences (50 487 novel transcripts) in these lines when compared with that of B73. Among these annotated genes, 76.16% genes covered >50% of the entire gene (Fig. 2). An average of 78.42% of genes could be uniquely mapped and 29 582 expressed genes were identified. After filtering, a total of 19 028 significant DEGs (FDR≤0.001, Log2Ratio≥1,andP≤0.001) were identified from each line as enriched in some metabolic pathways under droughtstressed conditions. Of these, 8 518 differentially expressed genes were up-regulated and 10 510 genes were down-regulated (Table 3). Interestingly,the more drought-tolerant genotypes were, the more up-regulated genes were. For example, the downregulated genes were 8.05 and 1.94 times than the up-regulated genes in drought-tolerant inbreds Qi 319 and Tie 7922, respectively, but 0.26 and 0.92 times in drought-susceptible inbreds CA339 and Ji 81162, respectively.

    Fig. 1 The sampled young tassels of 10 inbred lines. The left and right tassels of each line was under drought-stressed and wellwatered conditions (control), respectively.

    Fig. 2 An average of gene coverage percent by reads.

    The significant pathways identified in at least three inbred lines were further summarized in Table 4. RNA transport and mRNA surveillance were the most two pathways enriching in genes with binding functions, and were highly expressed in plant tassels under drought-strssed condition, except for CA339 with 23.62% of the enriched transcripts. Droughttolerant line Tie 7922 had the largest proportion (35.77%) of genes encoding proteins involved in RNA transport in droughtstressed tassels. The following enriched metabolic pathways were biosynthesis of secondary metabolites, amino acids,and macromolecular metabolites such as carbohydrates,lipids, glycans, terpenoids. Among these, pathways for secondary metabolites showed the greatest change in DEG profiles under drought stress, and genes encoding enzymes involved in the biosynthesis of these compounds accounted for greater than 11.34% of all DEGs in six inbreds, including Qi 319, X178, Ye 478, B73, Dan 340, and Ji 81162. Additionally,genes encoding enzymes involved in ‘phenylpropanoid and stilbenoid biosynthesis’ and ‘diarylheptanoid and gingerol biosynthesis’ in pathways for secondary metabolites biosynthesis were also up-regulated during drought response in early developing tassels (Table 4 and Fig. 3). Moreover,primary metabolic pathway categories affected by drought stress included genes expression for carbohydrate and lipid metabolism. Greater than 1.36% of all DEGs in six inbred lines (X178, Yu 12, Chang 7-2, Ye 478, B73, and Dan 340)corresponded to genes encoding enzymes were involved in galactose metabolism, and greater than 0.37% of all DEGs in five inbred lines (Qi 319, X178, Ye 478, B73, and Ji 81162)corresponded to genes encoding enzymes involved in fatty acid elongation (Table 4 and Fig. 3).

    3.3. GO term enrichment of the drought-stressed maize transcriptome

    The frequencies of DEGs under the top three most significantly enriched GO terms are shown in Fig. 4. In the‘molecular function’ category, most of the DEGs belonged to the ‘binding’ and ‘catalytic activity’ categories. Other terms such as ‘a(chǎn)ntioxidant activity’, ‘electron carrier activity’ and‘enzyme regulator activity’ represented an average of less than 10% of the DEGs affected by drought. For the ‘cellular component’ domain, the GO terms ‘cell’, ‘cell part’ and‘organelle’ represented over 60% of DEGs, indicating that drought stress had a great impact on cellular function. And we found that three membrane-related GO terms, were also enriched in our study, as Osakabeet al. (2014) showed that membrane transport played an important role in maintaining cellular homeostasis under drought stress. In the ‘biological process’ domain, 22 GO terms were enriched. Among these, except for ‘metabolic process’ and ‘cellular process’,comprising more than half of the DEGs, were affected by drought stress, some other biological processes, such as‘single-organism process’, and ‘response to stimulus’, were also affected by drought stress to some degree.

    Table 3 Number of differentially expressed genes in 10 maize inbred lines under drought-stressed condition

    Table 4 Drought-induced pathways and related differentially expressed gene with annotation in 10 inbred lines

    Fig. 3 The major response metabolisms under drought-stressed in young tassel derived from RNA-Seq. The three or more than three significantly common pathways in 10 inbred lines under drought stress were collected in black rectangles. These branch pathways was a portion of backbone metabolisms in blue boxes.

    Fig. 4 Main enrichment analysis of Gene Ontology (GO) annotation of differentially expressed genes (DEGs) in cellular component,molecular function and biological process. At least three tassels for each of inbred line with common GO annotations under drought stress were collected, the horizontal axis represents the percentages of DEGs involved in one GO term. The error bars show standard deviations for different lines.

    3.4. ldentification of DEGs under drought stress in developing maize tassels

    The largest proportion of DEGs from drought-stressed tassels was assigned to the KEGG pathways (Table 4).Of them, the exon-junction complex (EJC)-related genes were an important part of the ‘RNA transcript’ and ‘mRNA surveillance’ pathways that influence translation (Tangeet al. 2004). In this study, almost all annotated DEGs in the ‘mRNA surveillance pathway’ were enriched for genes encoding two components of the EJC, RNPS1 (RNA-binding protein with serine-rich domain 1) and SRm160 (Serine/arginine repetitive matrix protein 1). A total of 87 DEGs related to RNPS1 and SRm160 were identified in droughtstressed tassels in five or more inbred lines (Appendix A).This pathway annotation indicated that in drought-stressed young tassels, the EJC could regulate splicing, MAPK signaling, and transport by affecting the expression of genes encoding RNPS1 and SRm160.

    Nineteen DEGs representing genes encoding putative transcription factors were also noted in drought-stressed young tassels, including zinc finger proteins (4), bHLH(4), NAC (3), WRKY (2), bZIP (1), HB (1), SPL (1), LLY(1), Trihelix (1), and GRF (1). Genes encoding three transcription factors (hb8, GRMZM2G135447; c3h39,GRMZM5G801627; and dof26, GRMZM2G061292) showed increased expression in drought-stressed young tassels in at least four inbred lines. The expression ofdof26respondses to salt stress in maize seedlings (Chen and Cao 2015).However, the transcript abundance of genes encoding three members of the abscisic acid-, stress-, and ripening-induced(ASR)gene family (ZmASR1, GRMZM2G136910; ZmASR2,GRMZM5G854138; and ZmASR5, GRMZM2G052100),decreased in early developing tassels from most of inbred lines under drought stress. Transcripts of these genes accumulated in response to environmental stresses and might have effects in protecting cellular components and regulating gene expression (Perez-Diazet al. 2014). The exact regulation mechanism of these DEGs should be further investigated.

    3.5. Verification and expression of regulative DEGs in leaf, stem and root under drought stress

    For verification of transcription factors expressed in drought-stressed tassels, three of them,hb8,c3h39anddof26, were further tested with the inbred Qi 319 (Fig. 5).After treating with PEG-6000, all the three genes showed down-expressed in roots after 12 h treatment. Interestingly,c3h39anddof26showed nearly the same expression patterns in all treated tissues, which expression increased gradually, peaked at 3 h PEG treatment, then declined in stem. But the expression ofhb8was a little earlier than other genes, with the highest level at 1 h PEG treatment in all treated tissues, then it rapidly declined in roots and stems and remained higher expression level in leaves.

    4. Discussion

    4.1. The transcriptome of maize early developing tassels changes significantly under drought stress

    The floral organ differentiation stage is important for tassel and anther development. Except RNA regulation at the transcriptional level played important role in drought-stressed response in maize, a large number of DEGs encoding proteins involved in several metabolic pathways related to drought stress response were identified in early developing maize tassels under drought stress. Although gene expression profiles were highly variable among inbred lines, common and important regulation genes and their underlying pathways were identified in at least three inbreds. As reported in other species, the major metabolic pathways that were identified in tassels under drought stress might play important roles in the adaptation of drought, especially in the late period of drought stress. In drought-stressed tobacco root,transcripts of genes encoding enzymes involved in fatty acid,carbohydrate, and lipid metabolism-related pathways, as well as oxidoreductases, were greatly enriched, and a total of 17 pathways were significantly affected by drought stress (Yinet al. 2015). In diploidPaulowina australis, ‘photosynthesis’,‘biosynthesis of secondary metabolites’, and ‘phenylpropanoid biosynthesis’ pathways were significantly enriched in response to drought (Donget al. 2014). Transcriptome sequencing of moso bamboo (Phyllostachys edulis) at different flowering stages showed that transcripts related to‘genetic information processing’ were the most abundant,including those in ‘RNA transcript and mRNA surveillance’pathway (Gaoet al. 2014). In developing rice embryos, DEGs in pathways for carbohydrate, lipid, and energy metabolism were enriched, including those encoding enzymes involved in starch and sucrose metabolism, glycolysis, and fatty acid metabolism (Xuet al. 2012).

    Anthocyanin biosynthesis was related to steadystate leaf fluorescence of maize. In this study, two genes (GRMZM2G106443 and GRMZM5G896260),encoding glycosyltransferase and involved in anthocyanin biosynthesis, were down-regulated in Yu 12, which might keep the steady-state leaf fluorescence of Yu 12 after drought stress. GRMZM5G896260 was also down-regulated in Tie 7922 with relatively slightly decreased steady-state leaf fluorescence. These results suggested that the differential expression of genes encoding proteins involved in the regulation of RNA behavior, the biosynthesis of carbohydrates and lipids, and energy metabolism could likely be important for plant development and drought adaption.

    Fig. 5 Changes in transcript levels of three transcription factors under drought stress. Error bars represented standard deviation for the three biologic replicates of qRT-PCR experiments.

    4.2. Responses of genes encoding enzymes involved in carbohydrate metabolism to drought stress

    A large number of DEGs related to galactose metabolism were observed in drought-stressed tassels of six inbred lines (B73, Chang 7-2, Dan 340, X178, Ye 478, and Yu 12)(Appendix B and Fig. 6). Among these DEGs, the transcript abundance of gene (GRMZM2G052336), which encoded UDP-glucose 4-epimerase, decreased in drought-stressed tassels in four lines (B73, Dan 340, Ye 478, and Yu 12).Genes encoding invertase and alpha-glucosidase (bins A4 and A5 of Fig. 6, Appendix B), could break sucrose down into glucose and fructose (Liet al. 2011) and were down-regulated in drought-stressed tassels. Invertase activity is very important in early kernel development for converting sucrose originating from source leaves into hexoses to support endosperm and embryo cell division(Bateet al. 2004). Andersenet al. (2002) suggested that low acid invertase activity could decrease starch biosynthesis and lead to depletion of starch reserves.Enzyme 6-phosphofructokinase, involved in the regulation of glycolysis in eukaryotes, can phosphorylate fructose 6-phosphate to fructose-1,6-bisphosphate with ATP (Nariya and Inouye 2002). In our study, a large proportion of genes that encode 6-phosphofructokinase showed decreased transcript abundance, which could affect glucose and energy metabolism in drought-stressed tassels (bin A6 of Fig. 6,Appendix B). DEGs encoding 6-phosphofructokinase and invertase, includingMin1,Ivr2, and GRMZM2G095725(similar to cell wall invertase 2) were significantly downregulated in drought-stressed developing tassels. These results indicated that invertase activity might result in depleted glucose levels and ATP supplies, which might be lethal to tassels or anthers under drought.

    Fig. 6 Effects of drought stress on galactose metabolism. Drought stress responses associated with galactose metabolism in young tassel were shown, with a sucrose and glycolysis metabolism procedure in gray and green box, respectively. A1, hexose-1P-uridylyltransferase; A2, raffinose synthase; A3, α-D-galactosidase; A4, invertase; A5, glucosidase; A6, 6-phosphofructokinase.

    4.3. Responses of genes encoding enzymes involved in lipid metabolism to drought stress

    Lipid biosynthesis in maize embryo uses hexoses produced by invertases and reversible sucrose synthesis (Liet al.2011), which suggests that changes in expression of genes encoding invertases might also affect lipid metabolism. In our study, decreases in the expression of a large number of genes related to lipid metabolism were observed in early developing tassels under drought stress, including those in pathways for ‘fatty acid elongation’ and ‘cutin, suberin,and wax biosynthesis’ (Fig. 7, Table 5). The expression of genes encoding 3-ketoacyl-CoA synthase, the first step in fatty acid elongation that is important in long-chain fatty acid biosynthesis, decreased significantly in the early developing tassels of inbred lines B73, Ji 81162, Qi 319, X178, and Ye 478 under drought stress (bin A of Fig. 7). Verylong-chain fatty acids are precursors for cuticular wax biosynthesis activated by fatty acid synthesis and elongation in the endoplasmic reticulum by enzymes including three ketoacyl-CoA synthase (Lokeshet al. 2013). The fatty acid elongase function of 3-ketoacyl-CoA synthase determines fatty acid chain length and controls synthesis of subsequent elongated products including ketones and wax esters(Lokeshet al. 2013). Mutants such asCER1orWAX2in wax biosynthesis pathway genes generate conditional male sterile or glossy stem phenotypes (Aartset al. 1995;Wanget al. 2015).

    Fig. 7 The major genes to response drought in lipid metabolism. Drought stress responses associated with fatty acid elongation(A) and cutin, suberin and wax biosynthesis (B1, B2, C1, C2, and C3) in lipid metabolism were shown, with a table showing gene ID, annotation functions and fold change values of each line. The values in red and green indicated fold increase and decrease in expression of different lines in drought-induced tassels, respectively. The boxes in grayness indicated that the values of mapped reads in well-watered (control) and drought-stressed were less than 30. The N in green indicated that the values of expression in drought-stressed were zero.

    Table 5 Key fatty acid elongation genes and wax biosynthesis response to drought stress in young tassel

    Two DEGs, representing members of the cytochrome P-450 omega hydroxylase superfamily (GRMZM2G323830 and GRMZM2G022996) that encoding fatty acid hydroxylases(Matsunoet al. 2009), decreased in transcript abundance in early developing tassels under drought stress. These enzymes would be involved in ω-hydroxy fatty acid and longchain ketone biosynthesis (bins B1 and C2 of Fig. 7). The transcript abundances of two other genes in this superfamily(GRMZM2G062946 and GRMZM2G400129) increased slightly in line X178 (read values <30 under control and drought-stressed conditions). In addition, the expression of genes (GRMZM2G086210 and GRMZM2G141214)encoding glucose-methanol-choline (GMC) oxidoreductase,which oxidizes aliphatic alcohols to aldehydes during ω-dicarboxylic acid biosynthesis (Etxebesteet al. 2012),decreased significantly in early developing tassels under drought stress. Another fatty acid hydroxylase gene(GRMZM2G066578) also showed decreased expression in the drought-stressed tassels. However, the DEG encoding fatty aldehyde decarbonylase (GRMZM2G029912)showed increased transcript abundance in line X178.Interestingly, the expression of genes (GRMZM2G036217,GRMZM2G120938, andGRMZM2G480516), which encode fatty acyl reductases decreased significantly in early developing tassels under drought stress, which could have affected long-chain wax ester biosynthesis in tassels (bin C3 of Fig. 7). Fatty acyl reductases are essential for pollen development and could be related to male sterility under drought stress (Aartset al. 1997).Overall, these data suggest that the down-regulation under drought stress of the expression of these genes encoding enzymes involved in fatty acid elongation and suberin and wax biosynthetic pathways could have affected pollen development and could have led to male sterility in drought-stressed tassels.

    5. Conclusion

    Transcriptome of early developing tassel associated with drought tolerance in maize showed that most of the active metabolic pathways occurred at RNA level, especially the regulation of exon-junction complexes, thus we deduced that transcription regulation might play key roles in abiotic tolerance. Moreover, the carbohydrate and lipid metabolisms related to osmotic adjustment under water-stressed condition were involved in, with most encoding enzymes down-regulated. At the same time, a lot of important genes were identified, which would help to elucidate transcription regulation mechanisms related to drought tolerance.

    Acknowledgements

    This research was jointly funded by the National Natural Science Foundation of China (31661143010) and the Pilot Project of Breeding of the Seven Major Crops, China(2016YFD0101803).

    Appendicesassociated with this paper can be available on http://www.ChinaAgriSci.com/V2/En/appendix.htm

    Aarts M G M, Hodge R, Kalantidis K, Florack D, Wilson Z A,Mulligan B J, Stiekema W J, Scott R, Pereira A. 1997. TheArabidopsisMALE STERILITY 2 protein shares similarity with reductases in elongation/condensation complexes.The Plant Journal,12, 615–623.

    Aarts M G M, Keijzer C J, Stiekema W J, Pereira A. 1995.Molecular characterization of theCER1gene ofArabidopsisinvolved in epicuticular wax biosynthesis and pollen fertility.The Plant Cell,7, 2115–2127.

    Abdi H. 2007. The bonferonni and ?idák corrections for multiple comparisons. [2017-04-12]. http://www.utdallas.edu/~herve/Abdi-Bonferroni2007-pretty.pdf

    Andersen M N, Asch F, Wu Y, Jensen C R, Naested H,Mogensen V O, Koch K E. 2002. Soluble invertase expression is an early target of drought stress during the critical, abortion-sensitive phase of young ovary development in maize.Plant Physiology,130, 591–604.

    Bate N J, Niu X P, Wang Y W, Reimann K S, Helentjaris T G. 2004. An invertase inhibitor from maize localizes to the embryo surrounding region during early kernel development.Plant Physiology,134, 246–254.

    Benjamini Y, Yekutieli D. 2001. The control of the false discovery rate in multiple testing under dependency.The Annals of Statistics,29, 1165–1188.

    Chen Y Z, Cao J. 2015. Comparative analysis of dof transcription factor family in maize.Plant Molecular Biology Reporter,33, 1245–1258.

    Dong Y P, Fan G Q, Zhao Z L, Deng M J. 2014. Transcriptome expression profiling in response to drought stress inPaulownia Australis.International Journal of Molecular Sciences,15, 4583–4607.

    Etxebeste O, Herrero-Garcia E, Cortese MS, Garzia A,Oiartzabal-Arano E, de los Rios V, Ugalde U, Espeso E A. 2012.GmcAis a putative glucose-methanol-choline oxidoreductase required for the induction of asexual development inAspergillus nidulans.PLoS ONE,7, e40292.

    Fu J J, Cheng Y B, Linghu J J, Yang X H, Kang L, Zhang Z X, Zhang J, He C, Du X M, Peng Z Y, Wang B, Zhai L H,Dai C M, Xu J B, Wang W D, Li X R, Zheng J, Chen L,Luo L H, Liu J J,et al. 2013. RNA sequencing reveals the complex regulatory network in the maize kernel.Nature Communications,4, 2832.

    Gao J, Zhang Y, Zhang C L, Qi F Y, Li X P, Mu S H, Peng Z H. 2014. Characterization of the floral transcriptome of Moso Bamboo (Phyllostachys edulis) at different flowering developmental stages by transcriptome sequencing and RNA-Seq analysis.PLoS ONE,9, e98910.

    Hao Z F, Li X H, Liu X L, Xie C X, Li M S, Zhang D G, Zhang S H. 2010. Meta-analysis of constitutive and adaptive QTL for drought tolerance in maize.Euphytica,174, 165–177.

    Hao Z F, Li X H, Su Z J, Xie C X, Li M S, Liang X L, Weng J F, Zhang D G, Li L, Zhang S H. 2011. A proposed selection criterion for drought resistance across multiple environments in maize.Breeding Science,61, 101–108.

    Kakumanu A, Ambavaram M M R, Klumas C, Krishnan A,Batlang U, Myers E, Grene R, Pereira A. 2012. Effects of drought on gene expression in maize reproductive and leaf meristem tissue revealed by RNA-Seq.Plant Physiology,160, 846–867.

    Kanehisa M, Araki M, Goto S, Hattori M, Hirakawa M, Itoh M, Katayama T, Kawashima S, Okuda S, Tokimatsu T,Yamanishi Y. 2008. KEGG for linking genomes to life and the environment.Nucleic Acids Research,36, 480–484.

    Li E, Wang S, Liu Y, Chen J G, Douglas C J. 2011. OVATE FAMILY PROTEIN4 (OFP4) interaction with KNAT7 regulates secondary cell wall formation inArabidopsis thaliana.The Plant Journal,67, 328–341.

    Li R Q, Yu C, Li Y R, Lam T W, Yiu S M, Kristiansen K, Wang J. 2009. SOAP2: An improved ultrafast tool for short read alignment.Bioinformatics,25, 1966–1967.

    Liu C L, Hao Z F, Zhang D G, Xie C X, Li M S, Zhang X C,Yong H J, Zhang S H, Weng J F, Li X H. 2015. Genetic properties of 240 maize inbred lines and identity-by-descent segments revealed by high-density SNP markers.Molecular Breeding,35, 1–12.

    Livak K J, Schmittgen T D. 2001. Analysis of relative gene expression data using real-time quantitative PCR and the 2(T) (-Delta Delta C) method.Methods,25, 402–408.

    Lokesh, kiranmai, Pandurangaiah M, Babu S, Kumar N,Sudhakar C. 2013. Role of plant fatty acid elongase (3 keto acyl-CoA Synthase) gene in cuticular wax biosynthesis.Journal of Agriculture and Allied Sciences,2, 35–42.

    Matsuno M, Compagnon V, Schoch G A, Schmitt M, Debayle D, Bassard J E, Pollet B, Hehn A, Heintz D, Ullmann P, Lapierre C, Bernier F, Ehlting J, Werck-Reichhart D.2009. Evolution of a novel phenolic pathway for pollen development.Science,325, 1688–1692.

    Mortazavi A, Williams B A, Mccue K, Schaeffer L, Wold B. 2008.Mapping and quantifying mammalian transcriptomes by RNA-Seq.Nature Methods,5, 621–628.

    Nakashima K, Yamaguchi-Shinozaki K, Shinozaki K. 2014. The transcriptional regulatory network in the drought response and its crosstalk in abiotic stress responses including drought, cold, and heat.Frontiers in Plant Science,5,621–627.

    Nariya H, Inouye S. 2002. Activation of 6-phosphofructokinaseviaphosphorylation by Pkn4, a protein Ser/Thr kinase ofMyxococcus xanthus.Molecular Microbiology,46,1353–1366.

    Okoniewski M J, Miller C J. 2006. Hybridization interactions between probesets in short oligo microarrays lead to spurious correlations.BMC Bioinformatics,7, 276–279.

    Osakabe Y, Yamaguchi-Shinozaki K, Shinozaki K, Tran L S. 2014. ABA control of plant macroelement membrane transport systems in response to water deficit and high salinity.The New Phytologist,202, 35–49.

    Perez-Diaz J, Wu T M, Perez-Diaz R, Ruiz-Lara S, Hong C Y,Casaretto J A. 2014. Organ- and stress-specific expression of theASRgenes in rice.Plant Cell Reports,33, 61–73.

    Saini H S. 1997. Effects of water stress on male gametophyte development in plants.Sex Plant Reprot,10, 67–73.

    Shen C X, Li D, He R H, Fang Z, Xia Y M, Gao J, Shen H, Cao M L. 2014. Comparative transcriptome analysis of RNA-seq data for cold-tolerant and cold-sensitive rice genotypes under cold stress.Journal of Plant Biology,57,337–348.

    Shinozaki K, Yamaguchi-Shinozaki K, Seki M. 2003. Regulatory network of gene expression in the drought and cold stress responses.Current Opinion in Plant Biology,6, 410–417.

    Tange T, Nott A, Moore M J. 2004. The ever-increasing complexities of the exon junction complex.Current Opinion in Cell Biology,16, 279–284.

    Todaka D, Shinozaki K, Yamaguchi-Shinozaki K. 2015. Recent advances in the dissection of drought-stress regulatory networks and strategies for development of drought-tolerant transgenic rice plants.Frontiers in Plant Science,6, 84.

    Wang W J, Liu X W, Gai X S, Ren J J, Liu X F, Cai Y L, Wang Q, Ren H Z. 2015.Cucumis sativusL.WAX2plays a pivotal role in wax biosynthesis, influencing pollen fertility and plant biotic and abiotic stress responses.Plant and Cell Physiology,56, 1339–1354.

    Wang Z, Gerstein M, Snyder M. 2009. RNA-seq: A revolutionary tool for transcriptomics.Nature Reviews Genetics,10,57–63.

    Xiong L, Schumaker K S, Zhu J K. 2002. Cell signaling during cold, drought, and salt stress.The Plant Cell,14, 165–183.

    Xu H, Gao Y, Wang J B. 2012. Transcriptomic analysis of rice(Oryza sativa) developing embryos using the RNA-Seq technique.PLoS ONE,7, e30646.

    Yin F Q, Qin C, Gao J, Liu M, Luo X R, Zhang W Y, Liu H J, Liao X H, Shen Y U, Mao L K, Zhang Z M, Lin H J, Lubberstedt T,Pan G T. 2015. Genome-wide identification and analysis of drought-responsive genes and microRNAs in tobacco.Inter national Journal of Molecular Sciences,16, 5714–5740.

    Zhuang Y, Ren G J, Yue G D, Li Z X, Qu X, Hou G H, Zhu Y, Zhang J R. 2007. Effects of water-deficit stress on the transcriptomes of developing immature ear and tassel in maize.Plant Cell Reports,26, 2137–2147.

    国产精品一二三区在线看| 小蜜桃在线观看免费完整版高清| 国内精品一区二区在线观看| 日本色播在线视频| 综合色av麻豆| 97热精品久久久久久| 亚洲精品456在线播放app| 日日摸夜夜添夜夜添小说| 网址你懂的国产日韩在线| 日韩人妻高清精品专区| 九色成人免费人妻av| 久久久久九九精品影院| 国产精品一区二区免费欧美| 亚洲第一区二区三区不卡| 国产蜜桃级精品一区二区三区| 在线观看美女被高潮喷水网站| 国产真实乱freesex| 亚洲在线自拍视频| 色哟哟哟哟哟哟| 国产精品一及| 夜夜爽天天搞| 午夜久久久久精精品| 老熟妇乱子伦视频在线观看| 搡老岳熟女国产| or卡值多少钱| 在线看三级毛片| 亚洲av中文字字幕乱码综合| 22中文网久久字幕| 日韩强制内射视频| 色视频www国产| 亚洲成av人片在线播放无| 亚洲激情五月婷婷啪啪| 精品久久国产蜜桃| 99热全是精品| 不卡视频在线观看欧美| 亚洲欧美成人综合另类久久久 | 日韩一本色道免费dvd| 在线免费观看的www视频| 嫩草影院入口| 久久精品国产亚洲网站| 欧美成人免费av一区二区三区| 国产精品一区二区免费欧美| 91在线观看av| 亚洲性久久影院| 亚洲专区国产一区二区| 在线观看av片永久免费下载| 不卡一级毛片| 精品人妻偷拍中文字幕| 日本五十路高清| 精品一区二区免费观看| 国产精品日韩av在线免费观看| 日本a在线网址| 午夜福利视频1000在线观看| 少妇熟女欧美另类| 色尼玛亚洲综合影院| 啦啦啦观看免费观看视频高清| 亚洲性久久影院| 成年女人看的毛片在线观看| 最好的美女福利视频网| 久久久久久伊人网av| 卡戴珊不雅视频在线播放| 啦啦啦啦在线视频资源| 国产色婷婷99| 欧美激情国产日韩精品一区| 婷婷精品国产亚洲av| 一个人看的www免费观看视频| 亚洲一区二区三区色噜噜| 少妇人妻一区二区三区视频| 天堂√8在线中文| 卡戴珊不雅视频在线播放| 最近2019中文字幕mv第一页| 成熟少妇高潮喷水视频| 亚洲丝袜综合中文字幕| 不卡视频在线观看欧美| 一进一出抽搐gif免费好疼| 久久国产乱子免费精品| 久久人人爽人人爽人人片va| 最近的中文字幕免费完整| 人人妻人人看人人澡| 波多野结衣高清无吗| 久久久精品94久久精品| 久久久国产成人精品二区| 欧美激情国产日韩精品一区| 3wmmmm亚洲av在线观看| 亚洲va在线va天堂va国产| 国产综合懂色| 22中文网久久字幕| 亚洲性夜色夜夜综合| 成年免费大片在线观看| 国产精品99久久久久久久久| 九九久久精品国产亚洲av麻豆| 97热精品久久久久久| 十八禁国产超污无遮挡网站| 久久热精品热| 在线观看美女被高潮喷水网站| 如何舔出高潮| 黄色一级大片看看| 免费不卡的大黄色大毛片视频在线观看 | 男人狂女人下面高潮的视频| .国产精品久久| 香蕉av资源在线| 最近手机中文字幕大全| 精品一区二区三区av网在线观看| 免费看a级黄色片| 中文字幕熟女人妻在线| 一卡2卡三卡四卡精品乱码亚洲| 亚洲av二区三区四区| 狠狠狠狠99中文字幕| 女的被弄到高潮叫床怎么办| 亚洲电影在线观看av| 床上黄色一级片| 国产成人aa在线观看| 人人妻人人看人人澡| 亚洲欧美日韩无卡精品| videossex国产| 亚洲在线自拍视频| 久久精品人妻少妇| 精品免费久久久久久久清纯| 乱码一卡2卡4卡精品| 亚洲电影在线观看av| 女生性感内裤真人,穿戴方法视频| 国产成人aa在线观看| 国产午夜精品论理片| 五月伊人婷婷丁香| 高清毛片免费观看视频网站| 中文字幕av在线有码专区| 天天躁夜夜躁狠狠久久av| 18禁在线播放成人免费| 久久久久久久久久久丰满| 高清午夜精品一区二区三区 | 久久午夜亚洲精品久久| 久久久久久久久中文| 亚洲乱码一区二区免费版| 69人妻影院| 在线播放国产精品三级| 嫩草影视91久久| 超碰av人人做人人爽久久| 欧美成人免费av一区二区三区| 精品不卡国产一区二区三区| 国产一区二区三区在线臀色熟女| 国内久久婷婷六月综合欲色啪| 人妻夜夜爽99麻豆av| 久久久精品94久久精品| 亚洲国产精品久久男人天堂| 久久精品人妻少妇| videossex国产| 精品国内亚洲2022精品成人| 97人妻精品一区二区三区麻豆| 久久韩国三级中文字幕| 看片在线看免费视频| 国产 一区精品| 国产精品久久视频播放| 午夜视频国产福利| 国产成人91sexporn| 一本一本综合久久| 少妇熟女aⅴ在线视频| 婷婷六月久久综合丁香| 亚洲人与动物交配视频| 人妻制服诱惑在线中文字幕| 有码 亚洲区| 麻豆av噜噜一区二区三区| 国产黄片美女视频| 色综合站精品国产| 亚洲国产精品久久男人天堂| 婷婷精品国产亚洲av在线| 亚洲七黄色美女视频| 丝袜喷水一区| 久久人人精品亚洲av| 18禁在线无遮挡免费观看视频 | av在线天堂中文字幕| 亚洲国产精品久久男人天堂| 亚洲激情五月婷婷啪啪| 日韩欧美精品免费久久| 联通29元200g的流量卡| 久久精品久久久久久噜噜老黄 | 一个人观看的视频www高清免费观看| 日本一本二区三区精品| 午夜福利18| 日韩高清综合在线| 欧美日韩国产亚洲二区| av免费在线看不卡| 亚洲一区高清亚洲精品| 永久网站在线| 两个人视频免费观看高清| 国产综合懂色| 精品久久久久久久久亚洲| 国内精品宾馆在线| 成年版毛片免费区| 欧美一区二区国产精品久久精品| 国内久久婷婷六月综合欲色啪| 成人av一区二区三区在线看| 国产精品久久久久久av不卡| 搡老熟女国产l中国老女人| 免费av毛片视频| 熟女人妻精品中文字幕| 中文亚洲av片在线观看爽| 久久精品综合一区二区三区| 99热这里只有精品一区| 白带黄色成豆腐渣| 少妇丰满av| 综合色丁香网| 麻豆av噜噜一区二区三区| 麻豆一二三区av精品| 午夜福利高清视频| 变态另类丝袜制服| 国产精品美女特级片免费视频播放器| 久久久久精品国产欧美久久久| 欧美日韩国产亚洲二区| 国产高清不卡午夜福利| 欧美又色又爽又黄视频| АⅤ资源中文在线天堂| 午夜精品在线福利| 久久久精品大字幕| 欧美三级亚洲精品| 亚洲自偷自拍三级| 日韩高清综合在线| 天堂影院成人在线观看| 亚洲精品一卡2卡三卡4卡5卡| 色噜噜av男人的天堂激情| 欧美激情在线99| 欧美国产日韩亚洲一区| .国产精品久久| 淫秽高清视频在线观看| 99热这里只有是精品50| 欧美日韩在线观看h| 给我免费播放毛片高清在线观看| 国语自产精品视频在线第100页| 有码 亚洲区| 女的被弄到高潮叫床怎么办| 波多野结衣高清无吗| 少妇猛男粗大的猛烈进出视频 | 国产男人的电影天堂91| 亚洲欧美精品综合久久99| av.在线天堂| 深爱激情五月婷婷| 欧美潮喷喷水| 91麻豆精品激情在线观看国产| 欧美3d第一页| 不卡一级毛片| 免费搜索国产男女视频| 久久久久久久久久成人| 九九爱精品视频在线观看| 大型黄色视频在线免费观看| 国产精华一区二区三区| 看十八女毛片水多多多| 久久精品综合一区二区三区| av黄色大香蕉| 国产不卡一卡二| 午夜激情欧美在线| 国产成人精品久久久久久| 女人十人毛片免费观看3o分钟| 国产大屁股一区二区在线视频| 无遮挡黄片免费观看| 精品一区二区三区视频在线| 国产精品久久久久久久电影| 亚洲高清免费不卡视频| 黄色配什么色好看| 我要看日韩黄色一级片| 啦啦啦韩国在线观看视频| 婷婷色综合大香蕉| 少妇人妻精品综合一区二区 | 日本五十路高清| 国内精品宾馆在线| 日本成人三级电影网站| 俄罗斯特黄特色一大片| 国内久久婷婷六月综合欲色啪| 全区人妻精品视频| 国产一区二区在线观看日韩| 超碰av人人做人人爽久久| 免费看av在线观看网站| 黄色配什么色好看| 久久久久久久久大av| 精品欧美国产一区二区三| 真人做人爱边吃奶动态| 国产精品无大码| av免费在线看不卡| 亚洲精品国产成人久久av| 美女高潮的动态| 三级经典国产精品| 深夜精品福利| 熟女人妻精品中文字幕| 国产精品人妻久久久久久| av国产免费在线观看| 国产精品久久视频播放| 大香蕉久久网| 两个人视频免费观看高清| 免费av观看视频| 亚洲精品色激情综合| 国产91av在线免费观看| 国产 一区精品| 在线免费观看不下载黄p国产| 国产精品福利在线免费观看| 男插女下体视频免费在线播放| 色综合色国产| 亚洲国产精品成人综合色| 国产精品人妻久久久久久| 亚洲av成人av| 大又大粗又爽又黄少妇毛片口| 免费av毛片视频| 亚洲av不卡在线观看| 此物有八面人人有两片| 非洲黑人性xxxx精品又粗又长| 人妻夜夜爽99麻豆av| 国产成人影院久久av| 白带黄色成豆腐渣| 久久韩国三级中文字幕| 亚洲国产精品合色在线| 国产精品精品国产色婷婷| 精品国内亚洲2022精品成人| 国产国拍精品亚洲av在线观看| 在线播放国产精品三级| 国产成人91sexporn| 校园春色视频在线观看| 18禁黄网站禁片免费观看直播| 一进一出抽搐动态| 黄色一级大片看看| 日本与韩国留学比较| 九九在线视频观看精品| 伦理电影大哥的女人| 久久99热这里只有精品18| 亚洲精品久久国产高清桃花| 亚洲在线自拍视频| 亚洲真实伦在线观看| 亚洲aⅴ乱码一区二区在线播放| 国产综合懂色| 亚洲精品日韩在线中文字幕 | 有码 亚洲区| 晚上一个人看的免费电影| 亚洲电影在线观看av| 听说在线观看完整版免费高清| 国产91av在线免费观看| 麻豆久久精品国产亚洲av| 最近中文字幕高清免费大全6| 免费av观看视频| 成年女人永久免费观看视频| 成人精品一区二区免费| 久久久精品欧美日韩精品| 成人欧美大片| 女生性感内裤真人,穿戴方法视频| 国产色婷婷99| 热99在线观看视频| 给我免费播放毛片高清在线观看| 深爱激情五月婷婷| 欧美另类亚洲清纯唯美| 亚洲激情五月婷婷啪啪| 插阴视频在线观看视频| 两性午夜刺激爽爽歪歪视频在线观看| 亚洲欧美成人综合另类久久久 | 国产免费一级a男人的天堂| 丰满的人妻完整版| 精品免费久久久久久久清纯| 国产精品永久免费网站| 免费观看在线日韩| 免费看美女性在线毛片视频| 午夜免费激情av| 黄色欧美视频在线观看| 欧美日本亚洲视频在线播放| 性色avwww在线观看| 日韩成人伦理影院| av免费在线看不卡| 夜夜夜夜夜久久久久| 在线天堂最新版资源| 久久久久久久午夜电影| 高清毛片免费看| 国产一级毛片七仙女欲春2| 久久久久国内视频| 在线观看免费视频日本深夜| 日本欧美国产在线视频| 亚洲中文字幕一区二区三区有码在线看| 美女cb高潮喷水在线观看| 国产精品久久久久久久电影| 国产欧美日韩一区二区精品| 国产高清视频在线播放一区| 三级毛片av免费| 亚洲精品日韩在线中文字幕 | a级一级毛片免费在线观看| 国国产精品蜜臀av免费| 在线播放无遮挡| 欧美性感艳星| 久久鲁丝午夜福利片| 久久6这里有精品| 国产精品99久久久久久久久| 丝袜美腿在线中文| 波野结衣二区三区在线| 美女cb高潮喷水在线观看| 日本-黄色视频高清免费观看| 成人无遮挡网站| 高清午夜精品一区二区三区 | 嫩草影院新地址| 搞女人的毛片| 国产成人a∨麻豆精品| 极品教师在线视频| 菩萨蛮人人尽说江南好唐韦庄 | 插阴视频在线观看视频| av在线蜜桃| 久99久视频精品免费| 97碰自拍视频| 免费看日本二区| 91av网一区二区| 精品99又大又爽又粗少妇毛片| 男女做爰动态图高潮gif福利片| 久久久久久久久久久丰满| 我的女老师完整版在线观看| 啦啦啦韩国在线观看视频| 97在线视频观看| 国产精品人妻久久久久久| 久久久精品94久久精品| 欧美激情国产日韩精品一区| 嫩草影院入口| 寂寞人妻少妇视频99o| 一本一本综合久久| 亚洲精品日韩在线中文字幕 | 亚洲欧美日韩东京热| 少妇猛男粗大的猛烈进出视频 | 免费不卡的大黄色大毛片视频在线观看 | 91午夜精品亚洲一区二区三区| 狂野欧美白嫩少妇大欣赏| 日韩精品青青久久久久久| 99热这里只有是精品在线观看| 少妇人妻一区二区三区视频| 国产精品永久免费网站| 免费高清视频大片| 久久亚洲国产成人精品v| 国产精品女同一区二区软件| 国产白丝娇喘喷水9色精品| 日韩高清综合在线| 国内少妇人妻偷人精品xxx网站| 免费在线观看成人毛片| 美女免费视频网站| 一本精品99久久精品77| 精品久久国产蜜桃| 熟妇人妻久久中文字幕3abv| 国产成人a∨麻豆精品| 熟女电影av网| 51国产日韩欧美| 日本五十路高清| 晚上一个人看的免费电影| 婷婷六月久久综合丁香| 亚洲三级黄色毛片| 久久久精品欧美日韩精品| 国产精品乱码一区二三区的特点| 日韩高清综合在线| 中文字幕av成人在线电影| 啦啦啦啦在线视频资源| 精品一区二区免费观看| 欧美潮喷喷水| 搡老妇女老女人老熟妇| av天堂中文字幕网| 久99久视频精品免费| 亚洲乱码一区二区免费版| 国产精品久久电影中文字幕| 亚洲自偷自拍三级| 久久久久久久久中文| 91在线精品国自产拍蜜月| 精品不卡国产一区二区三区| 亚洲成人久久爱视频| 人人妻人人澡欧美一区二区| 亚洲图色成人| 有码 亚洲区| 欧美激情国产日韩精品一区| 亚洲精品久久国产高清桃花| 观看免费一级毛片| 一进一出好大好爽视频| 亚洲国产精品久久男人天堂| 亚洲人成网站在线观看播放| 人妻夜夜爽99麻豆av| 午夜福利视频1000在线观看| 婷婷精品国产亚洲av| 校园人妻丝袜中文字幕| 久久久欧美国产精品| 九色成人免费人妻av| 国产精品人妻久久久久久| 亚洲欧美日韩卡通动漫| 国内精品美女久久久久久| 日本撒尿小便嘘嘘汇集6| 午夜a级毛片| 国产成人a区在线观看| 中文字幕熟女人妻在线| 黄色日韩在线| 可以在线观看的亚洲视频| 一级a爱片免费观看的视频| 日本在线视频免费播放| 国产成年人精品一区二区| 国产在线男女| 69av精品久久久久久| 国产成人精品久久久久久| 在现免费观看毛片| 欧美三级亚洲精品| 久久久久九九精品影院| 国产亚洲精品久久久久久毛片| 99久久无色码亚洲精品果冻| 少妇熟女aⅴ在线视频| 噜噜噜噜噜久久久久久91| 日韩欧美国产在线观看| 日本 av在线| 九九爱精品视频在线观看| 午夜a级毛片| 免费黄网站久久成人精品| 日韩精品有码人妻一区| 亚洲欧美精品综合久久99| 欧美日韩精品成人综合77777| 国产单亲对白刺激| 2021天堂中文幕一二区在线观| 最近的中文字幕免费完整| 欧美+亚洲+日韩+国产| 国产精品美女特级片免费视频播放器| 欧美一区二区精品小视频在线| av女优亚洲男人天堂| 1000部很黄的大片| 久久久国产成人免费| 国产精品一区二区性色av| 国产片特级美女逼逼视频| 国产高清视频在线播放一区| 国产精品一区www在线观看| 亚洲国产高清在线一区二区三| 免费在线观看成人毛片| 一进一出好大好爽视频| 在线播放无遮挡| 午夜福利在线观看吧| 亚洲精品亚洲一区二区| 亚洲美女搞黄在线观看 | 久久久久久久亚洲中文字幕| 欧美一区二区亚洲| 亚洲国产高清在线一区二区三| 国产视频内射| 在线免费观看的www视频| 尾随美女入室| 国产亚洲av嫩草精品影院| 亚洲自拍偷在线| 又黄又爽又刺激的免费视频.| 国产国拍精品亚洲av在线观看| 可以在线观看毛片的网站| 国产精品一区二区三区四区免费观看 | 亚洲18禁久久av| 丝袜喷水一区| 国产成人福利小说| 久久久久久久久中文| 十八禁网站免费在线| 日本一本二区三区精品| 国产成年人精品一区二区| 精品国内亚洲2022精品成人| a级毛片免费高清观看在线播放| 久久久久国产精品人妻aⅴ院| 在线免费十八禁| 色噜噜av男人的天堂激情| 国产欧美日韩精品一区二区| 欧美成人a在线观看| 男人舔奶头视频| 国产乱人偷精品视频| 丰满乱子伦码专区| 国产黄片美女视频| 狂野欧美白嫩少妇大欣赏| 色综合亚洲欧美另类图片| 中文资源天堂在线| 亚洲精品国产av成人精品 | 色尼玛亚洲综合影院| 国产一区二区三区在线臀色熟女| 高清午夜精品一区二区三区 | 久久鲁丝午夜福利片| 在线免费十八禁| 搡女人真爽免费视频火全软件 | 免费高清视频大片| 亚洲av五月六月丁香网| 国产中年淑女户外野战色| 国产蜜桃级精品一区二区三区| 91在线观看av| 给我免费播放毛片高清在线观看| 亚洲精品成人久久久久久| 国产av麻豆久久久久久久| 亚洲婷婷狠狠爱综合网| 日日啪夜夜撸| 99热这里只有是精品50| 成人午夜高清在线视频| 又粗又爽又猛毛片免费看| 亚洲天堂国产精品一区在线| 少妇裸体淫交视频免费看高清| 我的老师免费观看完整版| 内地一区二区视频在线| 国产精品精品国产色婷婷| 亚洲中文字幕一区二区三区有码在线看| 女的被弄到高潮叫床怎么办| 久久久久久久午夜电影| 国产亚洲精品久久久久久毛片| 久久草成人影院| 国产一区二区亚洲精品在线观看| av福利片在线观看| 性欧美人与动物交配| 有码 亚洲区| 成人亚洲欧美一区二区av| 性欧美人与动物交配| 国产成人福利小说| 亚洲国产欧美人成| 国产老妇女一区| 男人狂女人下面高潮的视频| 有码 亚洲区| 欧美区成人在线视频| 成年女人毛片免费观看观看9| 亚洲国产欧美人成| 全区人妻精品视频| 亚洲av电影不卡..在线观看| a级毛片a级免费在线| 亚洲精品一卡2卡三卡4卡5卡| 日本黄色片子视频| 麻豆国产97在线/欧美| 欧美区成人在线视频| 如何舔出高潮| 亚洲最大成人av| 麻豆久久精品国产亚洲av| 国内久久婷婷六月综合欲色啪| 日本免费a在线| 久久精品人妻少妇| 精品久久久久久久久久免费视频| 免费搜索国产男女视频| 村上凉子中文字幕在线| 精品国内亚洲2022精品成人|