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

    Genetic basis of negative heterosis for growth traits in chickens revealed by genome-wide gene expression pattern analysis

    2021-12-17 11:49:36ChunningMaiChaoliangWenZhiyuanXuGuiyunXuSiruiChenJiangxiaZhengCongjiaoSunandNingYang

    Chunning Mai,Chaoliang Wen,Zhiyuan Xu,Guiyun Xu,Sirui Chen,Jiangxia Zheng,Congjiao Sunand Ning Yang

    Abstract

    Background:Heterosis is an important biological phenomenon that has been extensively utilized in agricultural breeding.However,negative heterosis is also pervasively observed in nature,which can cause unfavorable impacts on production performance.Compared with systematic studies of positive heterosis,the phenomenon of negative heterosis has been largely ignored in genetic studies and breeding programs,and the genetic mechanism of this phenomenon has not been thoroughly elucidated to date.Here,we used chickens,the most common agricultural animals worldwide,to determine the genetic and molecular mechanisms of negative heterosis.

    Results:We performed reciprocal crossing experiments with two distinct chicken lines and found that the body weight presented widely negative heterosis in the early growth of chickens.Negative heterosis of carcass traits was more common than positive heterosis,especially breast muscle mass,which was over?40% in reciprocal progenies.Genome-wide gene expression pattern analyses of breast muscle tissues revealed that nonadditivity,including dominance and overdominace,was the major gene inheritance pattern.Nonadditive genes,including a substantial number of genes encoding ATPase and NADH dehydrogenase,accounted for more than 68% of differentially expressed genes in reciprocal crosses(4257 of 5587 and 3617 of 5243,respectively).Moreover,nonadditive genes were significantly associated with the biological process of oxidative phosphorylation,which is the major metabolic pathway for energy release and animal growth and development.The detection of ATP content and ATPase activity for purebred and crossbred progenies further confirmed that chickens with lower muscle yield had lower ATP concentrations but higher hydrolysis activity,which supported the important role of oxidative phosphorylation in negative heterosis for growth traits in chickens.

    Conclusions:These findings revealed that nonadditive genes and their related oxidative phosphorylation were the major genetic and molecular factors in the negative heterosis of growth in chickens,which would be beneficial to future breeding strategies.

    Keywords:Chicken,Growth,Gene expression patterns,Heterosis,Oxidative phosphorylation

    Background

    Heterosis,first proposed by Shull in 1908[1],is defined as the deviation between F1reciprocal crosses and their parental lines mean[2].Heterosis has become a routine strategy for livestock and crop breeding and has driven great improvements in performance or livability over the last century.The fundamental mechanism underlying heterosis will determine whether it can be manipulated for the benefit of agriculture and biotechnology.Numerous studies have attempted to explain the genetic mechanism of heterosis,and three classic quantitative genetic hypotheses have been proposed:dominance[3,4],overdominance[1,5],and epistasis[6,7].However,these three models were mainly theoretical and could not provide a full explanation for the molecular basis and physiological causes of heterosis[8,9].

    At the molecular level,variation in gene expression is thought to constitute a significant source of phenotypic diversity[10].Investigation of differentially expressed genes(DEGs)between crossbred and their parent lines might contribute to improving our understanding of the molecular basis for heterosis.In particular,the gene expression patterns involved in some metabolic pathways are obviously correlated with heterotic phenotypes.For example,several studies on heterosis in rice revealed that differentially expressed genes between hybrids and their parents were involved in energy metabolism,which contributed in a significant way to the increased yield of hybrids[11,12].Fujimoto et al.[13]demonstrated that the higher photosynthetic efficiency of Arabidopsis hybrids was obviously associated with the early increase in the activity of genes involved in chlorophyll biosynthesis and photosynthesis,which contributed to increased heterotic biomass.Similar studies were also reported in animals.Hedgecock et al.[14]conducted transcriptomic analysis in Crassostrea gigas and found that nonadditive genes and their related protein metabolism played important roles in growth heterosis.

    In addition to positive heterosis,negative heterosis is also observed in many farm animals[15–19]and plants[20–23].Negative heterosis can cause unfavorable impacts on production performance in agriculture,such as negative heterosis of body weight in quails[15,16],carcass traits in beefs[17],hybrid necrosis in plants[20,21],and hybrid weakness of shoot dry weight[22]and height[23]in rice.The poultry industry has a long history of using crosses between different populations to take advantage of strain complementarity.In theory,the magnitude of heterosis is inversely correlated to the extent of genetic similarity between parental lines,and interspecific crosses show greater heterosis than intraspecific crosses[24].However,abundant evidence has revealed the existence of negative heterosis for growth traits in chickens when the genetic background of the parents varies greatly.Williams et al.[25]used the high-and low-body-weight chicken lines from Virginia Tech for heterosis analysis,and approximately?24% and?16% heterosis was observed for the body weights of reciprocal crosses at 4 and 8 weeks of age,respectively.Similar results were reported by Jull and Quinn[26],Maw[27],Liu et al.[28],and Sutherland et al.[29].Negative heterosis for carcass performance,especially for muscle mass in reciprocal crosses,was reported by Sun et al.[30]when broiler and layer chickens were used as parents.

    However,the genetic and molecular basis of heterosis for growth traits in chickens is still a mystery.Chickens are the most common and widespread domestic animals worldwide,as well as a great source of meat for humans.Revealing the genetic mechanisms of negative heterosis for growth traits will improve yield to meet the chicken meat demands of humans worldwide.Breast muscle is the largest proportion of body weight and is an important indicator of the growth rate in chickens[31];thus,investigating the negative heterosis of breast muscle mass can be a breakthrough to explore this unclear phenomenon.The development of transcriptome sequencing technologies has allowed unbiased and reproducible sequencing of whole transcriptomes,which are valuable for characterizing the patterns of gene expression and have been used to unravel the mechanisms of heterosis[32–34].In the current study,Cornish(C)and Rhode Island White(R)were used as parental lines to produce the F1generation.Cornish,as a standard broiler breed,has been selected for growth and,in particular,for muscle growth.Rhode Island White,a layer breed,has been intensively selected for egg production.The RNA-sequencing strategy was used to identify the transcriptomic differences in the breast muscle of reciprocal crosses and their parental lines.The objective of this study was to provide new insight into the molecular basis of negative heterosis for growth performance in chickens.

    Methods

    Experimental populations

    Two domesticated chicken breeds,Cornish(meat-type chicken,C line)and Rhode Island White(egg-type chicken,R line),from Beijing Huadu Yukou Poultry Industry Co.,Ltd.were employed as parents in this study to produce purebred progenies and reciprocal crosses(Fig.1a).The C line has been selected for 42-day body weight for 7 generations,whereas the R line has been intensively selected for total egg production to 300 days of age for 15 generations.We selected 10 males and 120 females from the 8th generation of the C line and 10 males and 80 females from the 16th generation of the R line as parents according to the following criteria:(i)cocks in each line with similar body weight and good semen quality and(ii)hens in each line with similar body weight and high egg production.These chickens were housed with individual cages in the same poultry facility.Each male(both C and R)was mated with 6 C and 4 R females by artificial insemination.The eggs were collected and recorded daily.Finally,a total of 632 chicks(347 for females and 285 for males)with clear pedigree information were hatched on the same day and used for subsequent studies.

    Fig.1 Body weight for purebred and crossbred progenies and heterosis of body weight for reciprocal crosses.A schematic diagram of two pure lines,Rhode Island White(R,layers)and Cornish(C,broilers),used as parental lines to produce the F1 generation(CC,CR,RC and RR),is shown in the upper left corner of the figure.Body weight of females(a)and males(b)from hatching to 8 weeks of age.Heterosis as a percentage(H%)of body weight for females(c)and males(d)in reciprocal crosses.The values of H% of body weight are shown below the figure.Values presented in blue and green indicate that H% was highly significant(P<0.01).Values presented in black indicate that H% was not statistically significant(P>0.05)

    Phenotypic measurement and sample collection

    At hatching,chicks were identified as males and females by vent sexing and then reared in separated cages under the same environment with free access to feed and water.The hatched chicks were wing-banded for individual identification.The four genetic combinations were reared in different cages to eliminate size disparities and reduce competition.The body weights were measured weekly from hatch to 8 weeks of age(Table S1).The length of the left shank and sternum for chickens were measured at 3,6,and 8weeks of age(Table S2).At 6 weeks of age,we randomly selected 64 female chickens(14,14,16 and 20 for CC,RR,CR and RC,respectively)and 44 male chickens(10,9,15 and 10 for CC,RR,CR and RC,respectively)from different half-sib families.These chickens were euthanized by cervical dislocation,and the following carcass traits were measured by an electronic balance:slaughter weight(measured after bloodletting),both left and right wing weight,breast muscle weight(pectoralis major and minor),drumstick weight(bone and muscle),and drumstick bone weight.The drumstick muscle weight was calculated as the drumstick weight subtracted from the drumstick bone weight.Heterosis as a percentage(H%)of the abovementioned traits was calculated according to the following equation:

    where F1iis the phenotype of individual i from reciprocal crosses;N is the number of birds in RC or CR.We obtained the P-value using the pt.function in the R program(https://www.r-project.org/)according to the tvalue and the degrees of freedom.H% was considered significant and highly significant if P-value<0.05,and P-value<0.01,respectively.Meanwhile,six female offspring(from 4–6 half-sib families)of each group,except for five female offspring in the RC group,were selected,and the left pectoralis major muscle of these chickens was isolated for subsequent RNA sequencing.

    RNA extraction and sequencing

    Total RNA was extracted using the TRIzol?Reagent(Invitrogen,USA)according to the manufacturer’s instructions and then dissolved in DEPC-treated water.To ensure that RNA was isolated successfully,the extracted RNA was first evaluated by 1% agarose gel electrophoresis.Then,the RNA purity,concentration and integrity of all eligible RNA extraction were determined by a NanoPhotometer? spectrophotometer(IMPL EN,CA,USA),a Qubit?2.0 Fluorimeter(Life Technologies,CA,USA),and an RNA Nano 6000 Assay Kit from the Bioanalyzer 2100 system(Agilent Technologies,CA,USA),respectively.The samples with an RNA integrity number value that greater than 7.0 were considered as high-quality RNA samples.A total of 23 samples were qualified for RNA sequencing library construction.Approximately 3 μg of RNA per sample was subjected to RNA-seq library construction using the NEBNext? UltraTM RNA Library Prep Kit(Illumina,USA)according to the manufacturer’s guide.After PCR amplification and purification,150 bp paired-end sequencing was performed on the Illumina Hiseq X Ten platform(Illumina Inc.,San Diego,CA,USA)and generated in nearly 750.13 million raw reads.

    Quality control and mapping

    To minimize mapping errors,reads that met the following parameter were removed:a)containing adaptors;b)with more than 10% unknown nucleotides;c)with more than 50% low-quality bases(Qphred≤20).The chicken reference genome(galGal5)and gene model annotation files were downloaded from the Ensembl database(ftp://ftp.ensembl.org/pub/release-91/).After quality control,over 721.42 million high quality reads with Q20>95%(Table S3)were aligned to the chicken reference genome using Hisat2(v2.0.5)[36].Approximately 76% of the high quality reads in each sample were mapped to the reference genome.Over 80% of reads were assigned to exonic regions,approximately 4% were assigned to intronic regions,and 16% were assigned to intergenic regions.

    Differential gene expression analyses

    The mapped reads of each sample were assembled by StringTie(v.1.3.3b)[37].The function of novel genes was annotated based on the Pfam database(v.31.0)[38].Then,the gene count matrix table was generated by featureCounts(v1.5.0-p3)[39].FPKM(fragments per kilobase million)values were extracted from the StringTie outputs.To enhance the statistical power for DEGs,the genes with an average FPKM<1 were removed.Meanwhile,the sex-linked genes were removed from the following analysis.After these steps,11,050 genes were filtered,and the remaining 11,544 genes were used for differential expression analysis between two purebred lines(CC vs.RR)and between reciprocal crosses and purebred lines(CR vs.CC,CR vs.RR,RC vs.CC and RC vs.RR)using the DESeq2 package(v.1.16.1)[40]in R project.We presented DESeq2,a method for differential analysis of count data,using the empirical Bayes shrinkage method to estimate dispersions and fold changes.The P-value was calculated by the Wald test.To control the false discovery rate,the resulting P-values were adjusted for multiple testing using the Benjamini-Hochberg method.Genes with an adjusted P-value<0.05 were considered differentially expressed genes in the corresponding comparison.

    Evaluation of differential inheritance patterns

    We used the average FPKM value of each group and the adjusted P-value to evaluate different inheritance patterns of genes(Table S4)[32].These genes were further classified into three inheritance patterns:additivity,dominance and overdominance,based on the level of gene expression exhibited by reciprocal crosses and parental lines.In brief,additivity(I and XII)occurred when the gene expression was significantly different between the two parental lines(adjusted P-value<0.05),and the gene expression of reciprocal crosses(CR or RC)was higher than one parental line but lower than the other parental line.Gene expression within CR/RC that was not significantly different from one parental line but significantly higher(or lower)than the other parental line was regarded as dominance(II,IV,IX,and XI).Gene expression within CR/RC that was significantly higher(or lower)than both parental lines(CC and RR)was considered overdominance(V,VI,VIII,III,VII,and X).

    To confirm the reliable of gene expression patterns reveled by RNA-seq,we performed quantitative realtime PCR(qRT-PCR)experiments.qRT-PCR reactions were performed with three technical replicates for each individual.The details of qRT-PCR and related results have been previously described[34].

    GO enrichment and KEGG pathway analyses

    To investigate the biological function of nonadditive genes involved,we performed functional enrichment analysis,including Gene Ontology(GO)categories and Kyoto Encyclopedia of Genes and Genomes(KEGG)pathways,using the Clusterprofile package[41]in the R project.The GO terms and KEGG pathways with FDR<0.05(BH method)were considered significant.

    ATP content assay

    ATP content was detected using an ATP assay kit(S0026B,Beyotime Biotechnology,China)as described in a previous study[42].The method is based on the theory that luciferase catalyzes luciferin to form fluorescence,which requires energy provided by ATP.Thus,the emitted fluorescence intensity is linearly related to the ATP concentration.Briefly,tissue samples(20 mg)were homogenized on ice with 150 μL of ice-cold assay buffer.It was then centrifuged at 12,000 r/min for 10 min at 4°C to remove insoluble materials,and the supernatant was collected.An aliquot(100μL)of ATP detection working solution was added to each well of a white 96-well plate.After incubation for 3min at room temperature,50μL of supernatant was added to the wells.Luminescence was measured by a fluorescence microplate reader.

    ATPase activity assay

    ATPase activity was assessed using an ATPase activity assay kit(MAK113,Millipore Sigma,St.Louis,MO,USA)according to the manufacturer’s instructions.ATPase hydrolyzes ATP into ADP and free phosphate.Free phosphate causes the malachite green reagent to form a stable dark green colorimetric product that is proportional to the ATPase activity.In brief,tissue samples(20mg)were homogenized on ice with 200μL of ice-cold assay buffer.They were then centrifuged at 14,000r/min for 10min at 4°C to remove insoluble materials,and the supernatant was collected.An aliquot(30μL)of the reaction mixture solution was added to each well of a 96-well flat-bottom plate and incubated for 30min at room temperature.Then,200μL of reagent was added to each well and incubated for an additional 30 min at room temperature to terminate the enzyme reaction.Finally,the absorbance was determined at 620 nm for all samples.

    Statistical analysis

    Differences in breast muscle weight,ATP content and ATPase activity among parental lines and reciprocal crosses were assessed using ANOVA followed by Tukey’s HSD test in the R program.The results were considered to be statistically significant when the adjusted P-value was less than 0.05.

    Results

    Negative Heterosis of body weights and carcass traits

    As described in Fig.1a,we chose the C and R breeds to produce purebred(CC and RR)and reciprocal crossbred progenies(RC and CR).The body weight of each progeny was measured weekly from hatching to 8 weeks of age.The traits of shank length and sternum length were measured at 3,6,and 8 weeks of age.The correlation among body weights at different ages for females and males varied from 0.115 to 0.991 and from 0.012 to 0.991,respectively(Fig.S1).As shown in the dynamic growth of parental lines and reciprocal crosses,the body weights of females and males in CR and RC from 2 to 8 weeks of age were lower than the average of CC and RR(Fig.1a,b and Table S1),although reciprocal crosses exhibited slightly positive heterosis for the length of shank and sternum at 6 and 8 weeks of age(Table S5).The degree of heterosis for body weight is displayed in Fig.1c and d in terms of heterosis as a percentage.The H% of body weight varied from?21.63% to 7.29% and from?16.62% to 7.26% for females and males,respectively.The negative heterosis of females and males reached a maximum value between the fifth and sixth weeks of age.In females,compared with CR(range of?13.16% to 7.29%),the H% was smaller in RC,which varied from?21.63% to?1.43%.In males,compared with CR(range of?14.56% to 7.26%),the H% was smaller in RC,which varied from?16.62% to 4.82%.In CR(females)and RC(males),the H% of body weight decreased from hatching to 5weeks of age and showed a slight increasing trend from 5 to 8weeks of age,while in RC(females)and CR(males),the inflection point was 6weeks of age.

    Given that fast growing broilers are mostly marketed at 6 weeks of age,we randomly slaughtered 108 chickens from four groups for carcass composition analysis at 42 days of age.As shown in Fig.2 and Table S6,most carcass traits,including slaughter weight,breast muscle weight,drumstick weight,and drumstick muscle weight,showed extremely significant negative heterosis(P<0.01)in reciprocal crosses of females and males.Among these various carcass traits,the negative heterosis of breast muscle weight was the largest,e.g.,?42.35%(CR)and?49.93%(RC)in females and?40.29%(CR)and?40.75%(RC)in males.Meanwhile,the correlation of body weight and breast muscle weight at 6weeks of age was 0.98 for both females and males(Fig.S2).

    Fig.2 Heterosis of carcass performance for purebred and crossbred progenies at 6 weeks of age.The breast muscle weight,drumstick bone weight,drumstick weight(bone and muscle)and wing weight were measured on both sides.The drumstick muscle weight was calculated as the drumstick weight subtracted from the drumstick bone weight.a-b Female.c-d Male.For(a)and(c),the dashed red line represents the midparent value.For(b)and(d),ns,*and**indicate that the heterosis as a percentage(H%)was not statistically significant(P>0.05),significant(P<0.05)and highly significant(P<0.01),respectively

    Inheritance of gene expression in reciprocal crosses

    As noted above,the negative heterosis of body weight and carcass traits were widespread in the present study.The more fundamental question is why the reciprocal progenies exhibited this phenomenon.Thus,the transcriptional data of breast muscle tissues for the four groups were used to analyze the differences in gene expression between parental lines and reciprocal crosses.A principal component analysis(PCA)was performed to visualize the differences in gene expression.The PCA plot showed that the four groups were obviously separated from each other(Fig.3a),indicating that there were visible differences in gene expression between the two parental lines and between the parental lines and reciprocal crosses.

    A total of 6253 DEGs between the two parental lines and between the reciprocal crosses and parental lines were identified(Fig.3b and Fig.S3),e.g.,5147(RR vs.CC),3205(CR vs.CC),1628(CR vs.RR),1184(RC vs.RR)and 4470(RC vs.CC).These DEGs were divided into 12 types(I,II,III,IV,V,VI,VII,VIII,IX,X,XI and XII;for details,see Table S4)based on the level of gene expression exhibited by reciprocal crosses and parental lines.The number of the 12 type genes in the CR and RC groups is shown in Fig.3c and Table S7.The 12 types were further classified into 3 main inheritance patterns:additivity(I,XII),dominance(II,IV,IX and XI),and overdominance(III,V,VI,VII,VIII and X).The number of dominant genes was 3324 and 3851 in CR and RC,respectively.The number of overdominant genes was 293 and 406 in RC and CR,respectively.Nonadditivity,including dominance and overdominace,was the major gene inheritance pattern.Nonadditive genes accounted for 68.99% and 76.20% of DEGs in CR and RC,respectively(Fig.3d).

    Fig.3 Analysis of gene inheritance patterns.a Principal component analysis of the reciprocal crosses(CR,RC)and the parental lines(RR,CC).b The number of DEGs among F1progenies.c Inheritance patterns of DEGs between reciprocal crosses and parental lines.DEGs were divided into 12 types,e.g.,class I,II,III,IV,V,VI,VII,VIII,IX,X,XI and XII,and further classified into three inheritance patterns:additivity(class I and XII),dominance(class II,IV,IX,and XI)and overdominance(class III,V,VI,VII,VIII,and X),based on the level of gene expression exhibited by reciprocal crosses and parental lines.Additivity,dominance,and overdominance are presented in blue,orange,and purple,respectively.Each class was accompanied by diagrams representing the relative expression levels of the maternal line(left dot),F1(middle dot),and paternal line(right dot).The number of DEGs in each class is shown above this class(green numbers,represented as the RC group)and below(blue numbers,represented as the CR group).d The proportion of additive,dominant and overdominant genes in DEGs

    Nonadditive inheritance is related to oxidative phosphorylation

    Apart from the focus on gene expression patterns in reciprocal crosses,we are more interested in the biological processes that additive and nonadditive genes are related to.We tested for enrichment of these genes against GO and KEGG pathways to detect the metabolic pathways involved.The functional enrichment analyses showed no significant GO terms or pathways detected in additive genes of the RC group,although one KEGG pathway,ribosome,was significantly enriched in additive genes of the CR group.However,as shown in Fig.4a,the dominant genes of CR and RC were both significantly enriched in 46 GO terms,including 2 GO terms of molecular function,32 GO terms of cell composition and 12 GO terms of biological process.The majority of these GO terms were associated with mitochondrial components and energy metabolism.Furthermore,the KEGG pathway analysis showed that one shared pathway,oxidative phosphorylation,was significantly enriched in the dominant genes of the CR and RC groups(Fig.4b).Additionally,overdominant genes of both RC and CR were also enriched in GO terms related to mitochondrial components and energy metabolism(Table S8).The oxidative phosphorylation pathway was also detected in the overdominant genes of the CR and RC groups(Fig.S4).

    Fig.4 Functional enrichment analysis for dominant genes.a Significant GO terms of dominant genes in reciprocal crosses.Each dot represents a GO term,and the size of a dot represents the number of genes enriched in the GO terms.The shade of the colored dots indicates the level of significance of the GO terms.The names of GO terms in purple,blue and brown represent the GO terms that belonged to molecular function,cell composition and biological process,respectively.b and c KEGG pathway analysis for dominant genes in the CR and RC groups,respectively.Each dot represents a KEGG pathway,and the size of a dot represents the number of genes enriched in the pathway.The color of a dot represents the KEGG classification in the pathway.The dashed red lines indicate significance levels(adjusted P-value<0.05).The dots that passed dashed red lines are regarded as significant pathways

    Given that nonadditive genes were significantly enriched in the pathway of oxidative phosphorylation,we further analyzed the nonadditive genes in the CR and RC groups.The number of nonadditive genes enriched in oxidative phosphorylation was 33 and 59 in the CR and RC groups,respectively(Fig.5a).These genes were related to NADH dehydrogenase,cytochrome c reductase,cytochrome c oxidase,ATP synthase,ATPase and succinate dehydrogenase.Among those,31 shared genes were detected in the CR and RC groups(Fig.5a and Table S9).Considering that nonadditive genes contained 10 different types,we further analyzed which type was important to the process of oxidative phosphorylation and found that types IV and II were the major gene expression patterns in CR and RC,respectively(Fig.5b).Type IV in CR and type II in RC accounted for 77.42%(24 of 31)and 61.29%(19 of 31)of shared genes,respectively.It is worth noting that the expression level of these genes in reciprocal progenies biased to the R line,and the expression level in the R line was significantly higher than that in the C line.

    Fig.5 Analyses of nonadditive genes enriched in oxidative phosphorylation of reciprocal crosses.a Overlap of nonadditive genes enriched in oxidative phosphorylation in the CR and RC groups.Each dot represents one gene.The names of the enzymes encoded by nonadditive genes are listed on the right.b Heatmap of shared gene expression levels in reciprocal crosses and their parents.The types of gene expression patterns in the CR and RC are shown on the left side.Schematic diagrams of the expression patterns and the number of genes are shown on the right.The color of the gene names represents the encoded enzyme,which is the same as the enzyme in plot(a)

    ATP content and ATPase activity detection

    The gene expression pattern results showed that nonadditive genes were related to the biological process of oxidative phosphorylation,implying that energy metabolism plays a vital role in negative heterosis of breast muscle.To further confirm the relationship between oxidative phosphorylation and negative heterosis for growth traits,we detected the ATP concentration and hydrolysis activity of breast muscle tissues for the CC,RR,CR,and RC groups.ATPases are a group of enzymes that catalyze the hydrolysis of ATP to form ADP.The detection of ATP content and ATPase activity revealed that the group with lower growth traits had lower ATP content but higher ATPase activity(Fig.6a-d).As shown in Fig.6a,breast muscle mass was significantly higher in CC than in RR,RC,and CR.The same trend was observed in body weight(Fig.6b)and ATP content(Fig.6c),and an opposite trend for ATPase activity is presented in Fig.6d.

    Fig.6 ATP content and ATPase activity detection.Difference analysis of breast muscle weight(a),body weight(b),ATP content(c),and ATPase activity(d)among RR,RC,CR,and CC.For a-d,each dot represents a sample.The central red dot represents the mean value of the corresponding group.***,**and*indicate adjusted P-values less than 0.001,0.01,and 0.05,respectively

    Discussion

    The utilization of heterosis has contributed tremendously to the increased productivity in many domesticated animals and crops for decades.In terms of the calculation formula,H% can be a positive or negative sign.Compared with the extensive studies on positive heterosis[11,13,43–45],the phenomenon of negative heterosis is overlooked in breeding programs and genetic studies,even though it exists widely in nature.In the present study,we observed that negative heterosis of body weight and carcass traits in juvenile chickens was more common than positive heterosis.This phenomenon was also reported by Williams et al.[25],Liu et al.[28]and Sutherland et al.[46].Among the carcass characteristics,we found that meat production displayed the largest negative heterosis in reciprocal crosses of females and males.The H% of breast and drumstick muscle weight was over?40% and?20%,respectively.These results were consistent with previous research showing that the negative heterosis of breast muscle weight in crosses was the largest among carcass traits when using broilers and layers as parents[30].Positive or negative heterosis does not imply superiority or inferiority since it depends on the trait’s biological significance and production preference[25,34].In livestock production,the negative heterosis of growth and meat yield was unfavorable since it reduced the edible carcass portions.We characterized the transcriptome profiles of breast muscle in reciprocal crosses and the parental lines herein to reveal the potential mechanisms of negative heterosis for growth traits in chickens.

    A large number of DEGs between reciprocal crosses and the parental lines were identified.The number of DEGs between two parental lines was greater than that between reciprocal crosses and their parental lines.This result was consistent with a previous study[11]and indicated that the genetic difference between two parental lines was larger than that between reciprocal crosses and their parental lines.Nonadditive genetic variance can result from a nonlinear phenotypic effect of alleles at one locus,as in the case of dominant or recessive allele pairs in classical genetics.Thus,the nonadditive expression pattern is critically important to the formation of heterosis[2,14].Recently,gene expression pattern analysis of chicken liver tissues revealed that overdominant genes related to lipid metabolism played a central role in the heterosis of fat deposition[34].Wu et al.[33]reported that dominant genes involved in carbohydrate metabolism were associated with heterosis for body weight in Drosophila melanogaster.In the present study,we classified these DEGs between reciprocal crosses and their parental lines into additivity,dominance and overdominance.Our results revealed that nonadditivity,including dominance and overdominance,was the major gene expression pattern in reciprocal crosses.Similar results were observed in Arabidopsis[2],Crassostrea gigas[14]and chickens[34].Previous reports in Medicago sativa[47]and Larix kaempferi[48]showed that the proportion of nonadditive genes in heterotic hybrids was higher than that in nonheterotic hybrids.It should be noted that nonadditive genes accounted for 76% of DEGs in the RC group,which was more than that observed in the CR group(69%),and the degree of negative heterosis for growth traits in the RC group was higher than that in the CR group.These results implied that the magnitude of the heterotic response was related to the proportion of genes with nonadditive expression.

    To better understand the molecular basis of negative heterosis, functional enrichment analysis was performed to gain insight into the biological relevance of nonadditive inheritance in reciprocal crosses.We found that the process of oxidative phosphorylation was significantly enriched in nonadditive genes of reciprocal crosses,indicating the special and crucial roles of energy metabolism in the negative heterosis of growth traits.Several previous studies have described the correlation of oxidative phosphorylation with heterosis in corn[49],wheat[50]and rice[12].Seymour et al.[2]found that the growth-related traits of Arabidopsis hybrids were associated with energy production via oxidative phosphorylation.This association was also reported in animals.McDaniel and Grimwood[51]demonstrated that heterosis of body weight in Drosophila melanogaster was correlated with oxidative phosphorylation efficiency.To validate the role of oxidative phosphorylation in negative heterosis of muscle yield in chickens,we further detected the ATP content and ATPase activity of breast muscle tissues for reciprocal crosses and parental lines.ATPase,as an essential enzyme in energy metabolism,catalyzes the hydrolysis of ATP to form ADP and harnesses the energy released from the breakdown of the phosphate bond to perform other cellular reactions.Our results showed that chickens with lower breast muscle weight had lower ATP content but higher ATPase activity,suggesting that chickens with higher ATP consumption had lower meat production.This finding corroborated that energy metabolism contributed strongly to negative heterosis and might help provide effective strategies for reducing the rate of ATP hydrolysis to improve muscle yield in crossbreds.Since the expression level of nonadditive genes involved in oxidative phosphorylation of reciprocal progenies biased to the R line (egg-type chicken),the expression level in the R line was significantly higher than that in the C line.Thus,the objective of reducing the rate of ATP hydrolysis might be achieved by decreasing the difference in parental weights or increasing the proportion of broiler parentage in the crossbred population.

    Growth is a complex polygenetic trait.To identify significant genes underlying the observed negative heterosis,we extracted the nonadditive genes detected in the process of oxidative phosphorylation in reciprocal crosses.A total of 31 shared genes were detected in reciprocal crosses.These genes encoding NADH dehydrogenase,cytochrome c reductase,cytochrome c oxidase,ATP synthase,ATPase and succinate dehydrogenase were all reported to be involved in the regulation of muscle growth and development,such as ATP5C1[52],ATP5G3[53],ATP5H[54,55],ATP5J[53],ATP6AP1[56],COX6B1[57],NDUFA1[52],NDUFA4[52],NDUFA5[52],NDUFA6[52],NDUFV2[53,58],NDUFS6[53],UQCR10[52],UQCR11[52],UQCRFS1[57],SDHA[52,53,57]and SDHB[52,53,58].Among these shared genes in reciprocal crosses,more than 60% of nonadditive genes exhibited a similar expression pattern to the layer line.The growth rate and body weight of layer chickens are considerably lower than those of broilers.It might be the large disparity of growth between layers and broilers and the differences in resource allocations that led to negative heterosis of growth traits in crossbred progenies.The Galgal5 may not be optimal chicken genome reference due to GRCg6a is available now,but should be sufficient to draw a conclusion that the important role of nonadditive genes and their related oxidative phosphorylation in negative heterosis for growth traits in chickens,since we confirmed the results by the detection of ATP content and ATPase activity.However,Gene expression is a dynamic process[59],and our research focused on gene expression analysis in juvenile chickens.We expect to determine whether the contributions of nonadditive genes would persist over time and to what degree they would impact the heterosis of growth traits in future experiments.In addition,the negative heterosis of growth traits in males was similar to that observed in females.However,sex-linked factors[60,61],such as hormones,may influence the growth rate.Therefore,further experiments should be performed to confirm that nonadditive genes and their related oxidative phosphorylation are also the major genetic and molecular factors in the negative heterosis of growth in males.

    Conclusions

    Our research focused on the phenomenon of heterosis in chickens and found that negative heterosis of growth traits was more common than positive heterosis,especially for muscle yield.Whole genome-wide gene expression pattern analysis showed that nonadditivity was the major mode of gene action in crossbred chickens.Nonadditive genes related to the biological process of oxidative phosphorylation played a critical role in the formation of negative heterosis for growth traits.Chickens with higher ATP consumption had lower muscle production.Our study revealed fundamental mechanisms of negative heterosis for growth traits in chickens and has important implications for muscle yield improvement.

    Supplementary Information

    The online version contains supplementary material available at https://doi.org/10.1186/s40104-021-00574-2.

    Additional file 1.Figure S1.Correlation among body weights of females and males at different ages.Figure S2.Correlation among body weight and muscle mass at 6 weeks of age.Figure S3.Volcano plot of differentially expressed genes between reciprocal crosses and parental lines.Figure S4.KEGG pathway analysis of overdominant genes in reciprocal crosses.

    Additional file 2:Table S1.Descriptive statistics for body weight of females and males from hatch to 8weeks of age.Table S2.Descriptive statistics for the left shank and sternum length of females and males during the experiment.Table S3.Summary statistics for transcriptome sequencing data.Table S4.Classification of different expression patterns of genes.Table S5.Descriptive statistics for heterosis of the length of shank and sternum for females and males during the experiment.Table S6.Heterosis of carcass performance for F1progenies at 6weeks of age.Table S7.Gene expression patterns of differentially expressed genes in reciprocal crosses.Table S8.GO terms of the top 20 overdominant genes in the RC and CR groups.Table S9.Detailed information on the overlap of nonadditive genes in RC and CR groups.

    Abbreviations

    C:Cornish;CC:The purebred progenies that used Cornish as parents;CR:The crossbred progenies that used Cornish as the paternal line and Rhode Island White as the maternal line;DEGs:Differentially expressed genes;FPKM:Fragments per kilobase million;GO:Gene Ontology;H%:Heterosis as a percentage;KEGG:Kyoto Encyclopedia of Genes and Genomes;R:Rhode Island White;RC:The crossbred progenies that used Rhode Island White as the paternal line and Cornish as the maternal line;RR:The purebred progenies that used Cornish as parents;PCA:Principal component analysis;qRT-PCR:Quantitative real-time PCR

    Acknowledgments

    We thank Dr.Guiqin Wu and Beijing Huadu Yukou Poultry Industry Co.,Ltd.,for providing the experimental chickens.

    Authors’contributions

    NY,CS and CW were involved in the conception of the study as well as the study design.CM and CW carried out the experiments.All authors participated in phenotypic and sample collection.CM,CW,NY and CS conducted bioinformatics and statistical analysis.CM detected the ATP content and ATPase activity of parental lines and reciprocal crosses.CM and CW designed all of the figures and wrote the manuscript.NY and CS were responsible for critical revisions of the manuscript drafts.All authors read and approved the final manuscript.

    Funding

    This work was supported by the National Natural Science Foundation of China(No.31930105),China Agriculture Research Systems(CARS-40)and China Postdoctoral Science Foundation(No.2020 M680028).

    Availability of data and materials

    The RNA sequencing data are available from the Sequence Read Archive(https://www.ncbi.nlm.nih.gov/sra) with BioProject number PRJNA524721.

    Declarations

    Ethics approval and consent to participate

    This study was carried out obeying the Guidelines for Experimental Animals established by the Animal Care and Use Committee of China Agricultural University.

    Consent for publication

    Not applicable.

    Competing interests

    The authors declare that they have no conflict of interest.

    Received:22 October 2020 Accepted:21 February 2021

    搡老熟女国产l中国老女人| cao死你这个sao货| 侵犯人妻中文字幕一二三四区| 亚洲综合色网址| 国产成人av教育| 欧洲精品卡2卡3卡4卡5卡区| 电影成人av| 欧美黑人精品巨大| 亚洲精品中文字幕在线视频| 久久国产精品影院| 丁香六月欧美| 自拍欧美九色日韩亚洲蝌蚪91| 最新在线观看一区二区三区| 男人的好看免费观看在线视频 | 如日韩欧美国产精品一区二区三区| 免费观看人在逋| 欧美日韩中文字幕国产精品一区二区三区 | 99精品久久久久人妻精品| 老司机午夜十八禁免费视频| 国产1区2区3区精品| 天堂俺去俺来也www色官网| 欧美精品高潮呻吟av久久| 久久精品国产亚洲av香蕉五月 | 在线天堂中文资源库| 亚洲伊人色综图| www.自偷自拍.com| 99久久人妻综合| 在线视频色国产色| 一a级毛片在线观看| 久久精品国产清高在天天线| 亚洲精品美女久久久久99蜜臀| 国产成人精品在线电影| 少妇裸体淫交视频免费看高清 | videos熟女内射| 久久久精品免费免费高清| 岛国毛片在线播放| 日韩欧美免费精品| 中国美女看黄片| 国产成人一区二区三区免费视频网站| 国产精品亚洲av一区麻豆| 熟女少妇亚洲综合色aaa.| 色尼玛亚洲综合影院| 免费av中文字幕在线| e午夜精品久久久久久久| 99热国产这里只有精品6| 性少妇av在线| 欧洲精品卡2卡3卡4卡5卡区| 高潮久久久久久久久久久不卡| 成人精品一区二区免费| av电影中文网址| 侵犯人妻中文字幕一二三四区| 天天躁夜夜躁狠狠躁躁| 国产又色又爽无遮挡免费看| 欧美激情 高清一区二区三区| 精品久久久久久久毛片微露脸| 精品福利永久在线观看| 久久午夜亚洲精品久久| 18禁裸乳无遮挡动漫免费视频| 久久ye,这里只有精品| 免费在线观看日本一区| 又黄又粗又硬又大视频| 黄网站色视频无遮挡免费观看| av中文乱码字幕在线| 国产成人系列免费观看| 国产精品久久久av美女十八| 美女高潮喷水抽搐中文字幕| 免费女性裸体啪啪无遮挡网站| 国产成人精品在线电影| 最新的欧美精品一区二区| 国产成人一区二区三区免费视频网站| 久久久久精品国产欧美久久久| 成人精品一区二区免费| 9热在线视频观看99| 久久久久久人人人人人| 久久狼人影院| 亚洲av日韩精品久久久久久密| 亚洲五月天丁香| x7x7x7水蜜桃| 黄色a级毛片大全视频| 精品久久久精品久久久| 欧美亚洲日本最大视频资源| 巨乳人妻的诱惑在线观看| 99在线人妻在线中文字幕 | 亚洲一区二区三区欧美精品| 欧美日韩一级在线毛片| 波多野结衣一区麻豆| 多毛熟女@视频| 在线观看免费高清a一片| 免费在线观看视频国产中文字幕亚洲| www日本在线高清视频| 美女视频免费永久观看网站| 精品一品国产午夜福利视频| 国产一区二区激情短视频| 国产精品 国内视频| 人人妻人人澡人人看| 欧美大码av| 亚洲中文日韩欧美视频| 美女国产高潮福利片在线看| 日本五十路高清| 中文字幕人妻丝袜制服| 精品一区二区三卡| 久久久久精品人妻al黑| 久久久久久久精品吃奶| 麻豆av在线久日| 两性夫妻黄色片| 国产精品 国内视频| 久久久久久久久久久久大奶| 视频在线观看一区二区三区| 怎么达到女性高潮| 999久久久精品免费观看国产| 在线永久观看黄色视频| 最新的欧美精品一区二区| 亚洲精品乱久久久久久| 手机成人av网站| 欧美在线黄色| 亚洲五月婷婷丁香| av国产精品久久久久影院| 国产精品久久久久久精品古装| 美女高潮喷水抽搐中文字幕| 午夜精品久久久久久毛片777| 怎么达到女性高潮| 亚洲精品在线观看二区| 天堂中文最新版在线下载| 国产精品亚洲av一区麻豆| 亚洲专区字幕在线| av在线播放免费不卡| 在线看a的网站| 中文字幕色久视频| 国产精品99久久99久久久不卡| 久久亚洲精品不卡| 嫁个100分男人电影在线观看| 国产高清国产精品国产三级| 亚洲精品中文字幕在线视频| 久久精品人人爽人人爽视色| 国产成人欧美在线观看 | 在线天堂中文资源库| 欧美日韩一级在线毛片| 国产精品久久久久成人av| 波多野结衣av一区二区av| 黄频高清免费视频| 久久ye,这里只有精品| 少妇猛男粗大的猛烈进出视频| 欧美亚洲日本最大视频资源| 国产一卡二卡三卡精品| 啪啪无遮挡十八禁网站| svipshipincom国产片| 欧美精品亚洲一区二区| 欧美黑人欧美精品刺激| 夜夜夜夜夜久久久久| 搡老熟女国产l中国老女人| 欧美精品高潮呻吟av久久| 大型av网站在线播放| 大片电影免费在线观看免费| 日本黄色视频三级网站网址 | 精品少妇久久久久久888优播| 亚洲精品国产区一区二| 国内毛片毛片毛片毛片毛片| 亚洲国产中文字幕在线视频| 不卡一级毛片| 大陆偷拍与自拍| 少妇的丰满在线观看| 欧美精品啪啪一区二区三区| 久久天堂一区二区三区四区| 在线观看免费午夜福利视频| 亚洲国产精品合色在线| 久久中文看片网| 中文欧美无线码| 亚洲 欧美一区二区三区| 亚洲全国av大片| 在线天堂中文资源库| 国产麻豆69| 韩国av一区二区三区四区| 免费少妇av软件| 国产精品久久视频播放| 黑人操中国人逼视频| 久久狼人影院| 国产精品秋霞免费鲁丝片| 亚洲国产精品合色在线| videosex国产| 国产片内射在线| 999精品在线视频| 国内毛片毛片毛片毛片毛片| 怎么达到女性高潮| 日韩欧美一区视频在线观看| 亚洲人成伊人成综合网2020| 国产亚洲一区二区精品| 亚洲精品中文字幕一二三四区| 成人影院久久| 日韩成人在线观看一区二区三区| a在线观看视频网站| 亚洲精品久久午夜乱码| 亚洲中文日韩欧美视频| 亚洲五月婷婷丁香| 最新美女视频免费是黄的| 国产精品98久久久久久宅男小说| 国产精品久久久av美女十八| 夜夜爽天天搞| 精品熟女少妇八av免费久了| 久久精品国产99精品国产亚洲性色 | 热99国产精品久久久久久7| 午夜福利乱码中文字幕| 久久精品成人免费网站| 免费观看a级毛片全部| 法律面前人人平等表现在哪些方面| av不卡在线播放| 欧美精品高潮呻吟av久久| 少妇裸体淫交视频免费看高清 | 他把我摸到了高潮在线观看| aaaaa片日本免费| 精品福利观看| 亚洲九九香蕉| 50天的宝宝边吃奶边哭怎么回事| 一区福利在线观看| 日本撒尿小便嘘嘘汇集6| 国产精品亚洲一级av第二区| 久久精品国产亚洲av高清一级| 精品一区二区三区av网在线观看| 亚洲成人免费av在线播放| 在线观看免费日韩欧美大片| 老司机深夜福利视频在线观看| 亚洲,欧美精品.| 夜夜夜夜夜久久久久| 天堂俺去俺来也www色官网| 久久人妻福利社区极品人妻图片| 精品一区二区三区av网在线观看| 黄色女人牲交| 国产区一区二久久| 欧美日韩亚洲高清精品| 欧美亚洲 丝袜 人妻 在线| 王馨瑶露胸无遮挡在线观看| 久久天堂一区二区三区四区| 亚洲熟女精品中文字幕| 国产精品98久久久久久宅男小说| 一夜夜www| 国产精品永久免费网站| 好男人电影高清在线观看| 狠狠婷婷综合久久久久久88av| 精品久久久久久久毛片微露脸| 国产成人精品久久二区二区免费| 日日爽夜夜爽网站| 欧美成人午夜精品| 午夜免费观看网址| 国产精华一区二区三区| 精品一区二区三区av网在线观看| 女人被躁到高潮嗷嗷叫费观| 久久精品人人爽人人爽视色| 日本撒尿小便嘘嘘汇集6| 女性被躁到高潮视频| 亚洲久久久国产精品| 久久久国产一区二区| 欧美+亚洲+日韩+国产| 久久精品91无色码中文字幕| 免费一级毛片在线播放高清视频 | 91精品国产国语对白视频| 国产精品免费视频内射| 黄色女人牲交| 18在线观看网站| 亚洲av电影在线进入| 欧美日韩av久久| 捣出白浆h1v1| 国产区一区二久久| 97人妻天天添夜夜摸| 久久狼人影院| 亚洲五月婷婷丁香| 50天的宝宝边吃奶边哭怎么回事| 免费看a级黄色片| 成在线人永久免费视频| 国产亚洲欧美精品永久| 高清欧美精品videossex| 国产真人三级小视频在线观看| 嫁个100分男人电影在线观看| 亚洲午夜理论影院| 日韩精品免费视频一区二区三区| 嫁个100分男人电影在线观看| 日日摸夜夜添夜夜添小说| 国产一区在线观看成人免费| 在线观看免费日韩欧美大片| 欧美大码av| 久久人人97超碰香蕉20202| 天天躁狠狠躁夜夜躁狠狠躁| 亚洲成人手机| 亚洲av成人一区二区三| 亚洲精品一卡2卡三卡4卡5卡| 曰老女人黄片| 国产主播在线观看一区二区| 91成年电影在线观看| 十八禁人妻一区二区| 国内久久婷婷六月综合欲色啪| 色老头精品视频在线观看| 欧美最黄视频在线播放免费 | av欧美777| 在线天堂中文资源库| 制服诱惑二区| 黑人欧美特级aaaaaa片| 日本黄色视频三级网站网址 | 国产成人影院久久av| 国产精品九九99| a在线观看视频网站| 啦啦啦视频在线资源免费观看| 人人妻,人人澡人人爽秒播| 久久99一区二区三区| 久久热在线av| 国产深夜福利视频在线观看| 18禁美女被吸乳视频| 国产精品一区二区在线不卡| 黄色女人牲交| 成人18禁在线播放| 午夜福利一区二区在线看| 精品久久蜜臀av无| 久久国产乱子伦精品免费另类| 国产精品美女特级片免费视频播放器 | 精品国产乱码久久久久久男人| 黄色怎么调成土黄色| 欧美黑人精品巨大| 色综合欧美亚洲国产小说| av中文乱码字幕在线| 免费人成视频x8x8入口观看| 欧美中文综合在线视频| 国产精品电影一区二区三区 | 桃红色精品国产亚洲av| 99久久99久久久精品蜜桃| 一本一本久久a久久精品综合妖精| 日本撒尿小便嘘嘘汇集6| av免费在线观看网站| 国产亚洲欧美98| aaaaa片日本免费| 午夜精品国产一区二区电影| 亚洲色图 男人天堂 中文字幕| 9热在线视频观看99| 日韩欧美在线二视频 | 国产精品免费大片| 日韩精品免费视频一区二区三区| 日本撒尿小便嘘嘘汇集6| av免费在线观看网站| 日日摸夜夜添夜夜添小说| 一级毛片女人18水好多| 国产深夜福利视频在线观看| 免费观看人在逋| 黑人猛操日本美女一级片| 少妇粗大呻吟视频| 一夜夜www| 脱女人内裤的视频| 国产又爽黄色视频| 18在线观看网站| 美女国产高潮福利片在线看| 欧美精品人与动牲交sv欧美| 97人妻天天添夜夜摸| 男女下面插进去视频免费观看| 欧美人与性动交α欧美软件| 欧美日韩乱码在线| 超碰成人久久| 69av精品久久久久久| 日韩免费高清中文字幕av| 国产精品自产拍在线观看55亚洲 | 国产av一区二区精品久久| 美女扒开内裤让男人捅视频| 少妇被粗大的猛进出69影院| 色尼玛亚洲综合影院| 成年动漫av网址| 校园春色视频在线观看| 久久精品国产亚洲av香蕉五月 | 老汉色av国产亚洲站长工具| 99精品在免费线老司机午夜| 亚洲精品国产精品久久久不卡| 精品国内亚洲2022精品成人 | 国产精品永久免费网站| 宅男免费午夜| 午夜免费鲁丝| 精品视频人人做人人爽| 少妇裸体淫交视频免费看高清 | 中国美女看黄片| 一级毛片精品| 精品一品国产午夜福利视频| 日韩一卡2卡3卡4卡2021年| 欧美日韩黄片免| 最近最新中文字幕大全免费视频| 久久久久久久久免费视频了| 在线视频色国产色| 国产极品粉嫩免费观看在线| 午夜激情av网站| 免费高清在线观看日韩| 国产精品 欧美亚洲| 亚洲伊人色综图| 在线观看免费日韩欧美大片| 欧美日韩国产mv在线观看视频| 日韩精品免费视频一区二区三区| 久久久国产成人免费| 人人妻人人添人人爽欧美一区卜| 久久久精品区二区三区| 俄罗斯特黄特色一大片| 成人18禁高潮啪啪吃奶动态图| 亚洲国产毛片av蜜桃av| 男人的好看免费观看在线视频 | 桃红色精品国产亚洲av| 不卡一级毛片| 久久久精品国产亚洲av高清涩受| 91成人精品电影| 成人国语在线视频| 国产精品香港三级国产av潘金莲| 亚洲欧洲精品一区二区精品久久久| 啦啦啦免费观看视频1| 熟女少妇亚洲综合色aaa.| 久久精品国产亚洲av香蕉五月 | 亚洲精品粉嫩美女一区| 看黄色毛片网站| 成人亚洲精品一区在线观看| 搡老熟女国产l中国老女人| 国产亚洲精品久久久久5区| 免费日韩欧美在线观看| 免费人成视频x8x8入口观看| 国产视频一区二区在线看| 成人手机av| 亚洲国产欧美日韩在线播放| 老司机午夜十八禁免费视频| 丁香欧美五月| 欧美国产精品一级二级三级| 久久中文看片网| 精品少妇久久久久久888优播| 19禁男女啪啪无遮挡网站| 久久国产亚洲av麻豆专区| 露出奶头的视频| 一本大道久久a久久精品| 亚洲国产精品一区二区三区在线| 少妇 在线观看| 在线观看免费日韩欧美大片| 好男人电影高清在线观看| av片东京热男人的天堂| 天天影视国产精品| 亚洲欧洲精品一区二区精品久久久| 99热只有精品国产| 在线观看免费高清a一片| 国产精品自产拍在线观看55亚洲 | 高清毛片免费观看视频网站 | 国产三级黄色录像| 精品乱码久久久久久99久播| av片东京热男人的天堂| 黄色 视频免费看| 国产精品 国内视频| 高清欧美精品videossex| 欧美日韩中文字幕国产精品一区二区三区 | 亚洲在线自拍视频| 在线天堂中文资源库| 色播在线永久视频| 日本a在线网址| 一级黄色大片毛片| 午夜精品久久久久久毛片777| 99久久人妻综合| 老司机午夜福利在线观看视频| 精品人妻在线不人妻| 成人影院久久| 热99久久久久精品小说推荐| 亚洲avbb在线观看| 不卡一级毛片| 久久久久久久国产电影| 亚洲九九香蕉| 看片在线看免费视频| 少妇被粗大的猛进出69影院| 99热网站在线观看| 欧美亚洲日本最大视频资源| 视频在线观看一区二区三区| 欧美精品啪啪一区二区三区| 亚洲性夜色夜夜综合| 女人爽到高潮嗷嗷叫在线视频| 中文字幕人妻丝袜制服| 国产视频一区二区在线看| 久久香蕉激情| 国产一区有黄有色的免费视频| 免费久久久久久久精品成人欧美视频| 久久久精品免费免费高清| 高清av免费在线| 亚洲av成人不卡在线观看播放网| 国产黄色免费在线视频| 亚洲精品久久午夜乱码| 80岁老熟妇乱子伦牲交| 成人三级做爰电影| 波多野结衣av一区二区av| 国产麻豆69| 精品电影一区二区在线| 村上凉子中文字幕在线| 婷婷丁香在线五月| 国产精品免费一区二区三区在线 | 久久久久久久国产电影| 精品午夜福利视频在线观看一区| 国产在线观看jvid| 免费在线观看亚洲国产| bbb黄色大片| 成人永久免费在线观看视频| 亚洲精品美女久久久久99蜜臀| 正在播放国产对白刺激| 亚洲欧美色中文字幕在线| 久久久国产成人精品二区 | 日韩三级视频一区二区三区| 欧美 亚洲 国产 日韩一| a在线观看视频网站| 自线自在国产av| 亚洲综合色网址| 美女高潮喷水抽搐中文字幕| av电影中文网址| 国产黄色免费在线视频| 国产在线一区二区三区精| 美女扒开内裤让男人捅视频| 啦啦啦 在线观看视频| 亚洲专区字幕在线| 国产精品偷伦视频观看了| 嫩草影视91久久| 亚洲国产中文字幕在线视频| 亚洲精品国产精品久久久不卡| 欧美国产精品va在线观看不卡| 又大又爽又粗| av中文乱码字幕在线| 丰满人妻熟妇乱又伦精品不卡| 国产精品九九99| 亚洲一区中文字幕在线| 欧美色视频一区免费| 亚洲美女黄片视频| 久久国产精品人妻蜜桃| 大陆偷拍与自拍| 黑人巨大精品欧美一区二区mp4| bbb黄色大片| 亚洲精品久久成人aⅴ小说| 日日夜夜操网爽| 大型av网站在线播放| 悠悠久久av| 中出人妻视频一区二区| 美女高潮喷水抽搐中文字幕| 亚洲自偷自拍图片 自拍| 久久性视频一级片| 欧美精品高潮呻吟av久久| 久久人人爽av亚洲精品天堂| 午夜亚洲福利在线播放| 999精品在线视频| 久久久久国产一级毛片高清牌| 国产一卡二卡三卡精品| 美国免费a级毛片| 黄片大片在线免费观看| 亚洲成人免费av在线播放| 久久国产精品人妻蜜桃| 欧美日韩乱码在线| 身体一侧抽搐| 国产亚洲欧美精品永久| 国产精品九九99| 精品福利观看| 久久这里只有精品19| 国产蜜桃级精品一区二区三区 | 亚洲av欧美aⅴ国产| 国产欧美日韩一区二区三区在线| 法律面前人人平等表现在哪些方面| 国产精品九九99| 午夜激情av网站| 日本wwww免费看| 热99久久久久精品小说推荐| 天天躁日日躁夜夜躁夜夜| 丝袜人妻中文字幕| 色在线成人网| 高潮久久久久久久久久久不卡| 久久精品国产亚洲av高清一级| 悠悠久久av| 欧美中文综合在线视频| 国产精品一区二区在线观看99| 黑人巨大精品欧美一区二区mp4| 操美女的视频在线观看| 亚洲成人国产一区在线观看| 欧美乱码精品一区二区三区| 两性夫妻黄色片| 国产一区在线观看成人免费| 免费一级毛片在线播放高清视频 | 亚洲熟妇中文字幕五十中出 | 久久久久久免费高清国产稀缺| 日本vs欧美在线观看视频| 亚洲国产精品sss在线观看 | 国产日韩一区二区三区精品不卡| 久久人妻福利社区极品人妻图片| 亚洲精品中文字幕在线视频| 熟女少妇亚洲综合色aaa.| 久久精品人人爽人人爽视色| 一边摸一边做爽爽视频免费| 国产亚洲欧美98| 天堂动漫精品| 欧美日韩精品网址| aaaaa片日本免费| 亚洲男人天堂网一区| 超碰成人久久| 日本一区二区免费在线视频| 欧美午夜高清在线| 老司机深夜福利视频在线观看| 午夜精品久久久久久毛片777| 一区二区日韩欧美中文字幕| 王馨瑶露胸无遮挡在线观看| 两性夫妻黄色片| 一区福利在线观看| 亚洲欧美一区二区三区久久| 青草久久国产| 一区福利在线观看| 亚洲欧美一区二区三区久久| tocl精华| 一a级毛片在线观看| 午夜日韩欧美国产| 十八禁人妻一区二区| 国产又爽黄色视频| 热re99久久国产66热| 久久精品亚洲av国产电影网| 丁香欧美五月| 一区二区三区国产精品乱码| 大香蕉久久成人网| 国产精品香港三级国产av潘金莲| 精品国产国语对白av| 亚洲熟女精品中文字幕| 久久精品国产亚洲av香蕉五月 | 丝袜人妻中文字幕| 在线观看一区二区三区激情| 黄色毛片三级朝国网站| 最新的欧美精品一区二区| 欧美在线黄色| netflix在线观看网站| 两性夫妻黄色片|