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

    Genetic dissection and genomic prediction for pork cuts and carcass morphology traits in pig

    2023-12-18 08:49:20LeiXieJiangtaoQinLinRaoDengshuaiCuiXiTangLiqingChenShijunXiaoZhiyanZhangandLushengHuang

    Lei Xie, Jiangtao Qin, Lin Rao, Dengshuai Cui, Xi Tang, Liqing Chen, Shijun Xiao, Zhiyan Zhang* and Lusheng Huang*

    Abstract

    Keywords Carcass morphology traits, Genomic selection, Genotype imputation, GWAS, Pork cuts

    Background

    Pig carcass cutting is a process of decomposing the postmortem carcass into various cuts with different sizes and weights according to the tissue structure of different anatomical parts, followed by trimming, cooling, packaging, and preservation.The economic value of pork cuts varies depending on their quantity and quality.Different pork cuts also require diverse cooking and processing methods [1-3].In recent years, the outbreak of African swine fever in China has issued many policies restricting the transportation of live pigs to prevent the spread of the virus [4, 5], which has presented a new opportunity for the development of chilled meat.Furthermore, due to the rise in living standards and the fast-paced lifestyle,consumers have shifted their pork consumption habits from purchasing hot carcasses for direct cutting and selling to opting for pre-packaged chilled meat that suits their preferences [1, 6].This further led to the widespread acceptance and adaptation of chilled meat by the majority of consumers.Consequently, many pig companies are rapidly deploying slaughterhouses and expanding their slaughter-processing capabilities within their production chain to optimize the carcass economic value.Mapping and identifying quantitative trait loci (QTLs) for pork cuts will help to breed merit pigs with higher proportion of expensive pre-cut products to increase the overall value of cuts.To the best of our knowledge, there is a lack of reports on QTLs and causal genes that affect pork cuts, as well as investigations into genomic selection or the evaluation of prediction accuracy for pork cuts.The identification of QTLs and investigation of the genetic mechanisms of pork cut attributes serve as the foundation for enhancing the economic value of pork cuts by improving the accuracy of genomic selective breeding.

    In this study, 17 pork cuts and 12 carcass morphology traits were measured on 2,012 pigs from four populations genotyped using the CC1 PorcineSNP50K BeadChip(CC1 Chip) [7, 8].The aim of this study was to identify QTLs affecting proportion of pork cuts to evaluate the accuracy of selection and the feasibility of industrial application.We employed imputation-based wholegenome sequence (WGS) association analysis to uncover potential causal mutations and major genes affecting pork cuts, comparing it with haplotype-based CC1 Chip genotyping data association analysis [9-11].These results are essential for pig companies who aim to enhance their advantage in the consumer market, core competitiveness, and brand value.Moreover, genetic dissection of pork cuts is vital for understanding carcass composition,which provides critical reference for studying regulatory mechanisms of skeletal and muscle growth and development in different parts of pigs.

    Materials and methods

    Animals, feeding and sampling

    A total of 2,012 pigs were randomly sampled from Muyuan Food Co., Ltd.(Henan, China) for pork cut evaluation, as described by Xei et al.[12].The experimental pigs including 265 Landrace (LR, 95 sows and 170 barrows),698 Yorkshire (YK, 435 sows and 263 barrows), 689 Landrace × Yorkshire hybrid (LY, 402 sows and 287 barrows),and 258 Duroc × Landrace × Yorkshire hybrid (DLY, 115 sows and 143 barrows).All pigs were raised under consistent feeding environments and nutritional conditions,and they were provided with the same commercial diets and had unrestricted access to water.More details of breeding environment and pedigree family structure were described in our previous study [12].Each time approximately 100 pigs were randomly selected from 500 to 1,000 market-aged pigs for slaughter testing.A total of 22 batches of pigs were measured for pork cuts and carcass morphology traits (Table S1).These pigs were uniformly slaughtered centrally, following the specifications described in the Operating Procedures of Livestock and Poultry Slaughtering - Pig (GB/T 17236-2019) [13], at an average age of 180 d.

    Phenotypic determination

    Twelve carcass morphology traits were measured for all individuals, including carcass straight length (SL),oblique length (OL), thoracic number (THN), lumbar number (LUN), thoracic length (THL), lumbar length(LUL), single lumbar length (SLUL), shoulder backfat depth (SBD), 6th_7th rib backfat depth (RBD), waist backfat depth (WBD), hip backfat depth (HBD), and the mean of backfat depth (MBD).Additionally, the carcass was cut into 17 pork cuts as shown in Fig.1, and their weight was measured, including three primal cuts (shoulder cut (SC), middle cut (MC), leg cut (LC)), and 14 subprimal cuts (boneless Boston shoulder (BBS), boneless picnic shoulder (BPS), front ribs (FR), fore leg bones(FLB), scapula bone (SB), loin (LO), belly (BE), ribs (RI),chine bones (CB), back fat (BF), boneless leg (BL), tenderloin (TL), hind leg bones (HLB), tail and pelvis bone(TPB)).The determination methods and processes were described in the previously published study [12].Each pork cut was carefully weighed and measured by the investigators.The proportion of each pork cut was determined through the division of the weight of pork cut by the weight of the entire carcass.

    Genotyping

    Genomic DNA was extracted from the muscle tissue of each animal using the routine phenol/chloroform extraction method.Individuals were genotyped using the CC1 PorcineSNP50 BeadChip (51,368 SNPs) [7, 8] according to the manufacturer’s protocol.The marker density and accuracy of the CC1 Chip were described in our prior study [7, 8].Thresholds of individual call rates > 90%,SNPs call rates < 95%, minor allele frequency (MAF) < 5%,and Hardy-Weinberg disequilibrium (P< 10-5) were filtered out using PLINK (v1.90b6.24) [14].After quality control, 40,016 SNPs and 2,012 animals were retained for further analysis.

    Fig.1 Standardized pork cuts and their corresponding pork cuts in a commercial pig carcass

    Imputation of whole-genome sequence variants

    Genotype imputation for the experimental population was performed using IMPUTE5 [15] from a highquality haplotype reference panel including 42,620,918 variants as described by Tong et al.[16].The haplotype reference panel included whole-genome sequencing data of 1,096 samples from 43 pig breeds (n≥ 3)with an average sequencing depth of 17.1 X.The detailed imputation process was described in our previous study [17].Variants were called using GATK following the best practice flowchart and were quality controlled by following criteria: (1) SNP: QD < 2.0,QUAL < 30.0, MQ < 40.0, SOR > 3.0, FS > 60.0, MQRank-Sum < -2.5, ReadPosRankSum < -8.0; (2) INDEL:QD < 2.0, QUAL < 30.0, MQ < 40.0, FS > 200.0, Read-PosRankSum < -20.0.SNPs in the target panel were further filtered with call rate < 95%, or minor allele frequency (MAF) < 5%, or Hardy-Weinberg disequilibrium (P< 10E -5) by PLINK (v1.90b6.24) [14].The haplotypes of the target panel (Sscrofa 11.1) were constructed by SHAPEIT4.2 [18] and PHASEBOOK [19].Then, genotype imputation was performed between the target and reference panels by IMPUTE5 with default parameters [15].The imputation accuracy was evaluated by an internal cross-validation solution of IMPUTE5.Specifically, the genotypes of one locus in all individuals in the target panel were masked at a time,and then the masked genotypes were imputed with the haplotype information from the reference panel.The genotypic concordance rate and squared correlation(R2) between original genotypes from the target panel and imputed genotypes were calculated as imputation accuracies.The accuracies (MeanR2/concordance rate)of the imputed genotypes for the experimental population were 0.89/99.16%, which implied a high quality of the imputed genotypes.

    Imputation-based of whole-genome sequence GWAS(IGWAS)

    Single locus association analysis was conducted using the GEMMA software (version 0.98.1) [20] with a linear mixed model (LMM) that accounts for SNP-based population structure and relatedness between individuals.

    Haplotype-based CC1 Chip genotyping data GWAS(HGWAS)

    Haplotypes of the SNP genotypes were constructed by PHASEBOOK [19].It assumes that all haplotypes in the population can be traced back to a predetermined number (K = 10) of ancestral haplotypes [26].Then a hidden Markov model was employed to infer the ancestral haplotypes inherited by each individual at each locus [21].To detect the association between phenotypes and the haplotype status, a linear mixed framework was used similar to single locus association with a difference in the incidence matrices.In this model,Xis the incidence matrices of the ancestral haplotypes rather than SNP genotypes.The haplotype effects were fitted as random effects.Gis the kinship matrix calculated from the SNP genotypes using VanRaden’s method.

    Statistical models to genome prediction

    Two genomic selection models were implemented to evaluated the genomic accuracy of pork cuts.(1)Genomic best linear unbiased prediction (GBLUP) [21],which is the most widely used model in genome breeding practice.The mixed linear model is as follows:

    where y is the vector of phenotype, μ is the overall mean, α is the fixed effect including sex, populations and slaughter batches, a is the vector of genomic breeding values of all individuals, e is the vector of residuals, 1 is a vector of ones, W is the indicator matrix of α , and Z is the indicator matrix of a.Assume that e follows a normal distribution of N(0,Iσ2e), andafollows a normal distribution of N(0,Gσ2a).Where σ2ais the additive genetic variance,Gis the kinship matrix obtained from genotype data (included CC1 PorcineSNP50 BeadChip genotype and genome-wide imputed SNPs), which was calculated using VanRaden’s method [21], and the detailed calculation method can be found in Yang et al.[27].Then a is solved from the mixed model equations (MME) [28].In this study, the MME formula is solved by using GCTA software [29], and the estimated genome breeding value of the individual is ^a.

    (2) Bayesian sparse linear mixed model (BSLMM),which assumes that the effects of markers follow a mixture of two normal distributions [30].It assumes that all markers have at least a small effect, but some proportion of markers have an additional large effect.The model consists of a standard linear mixed model, with one random effects term, and with sparsity inducing priors on the regression coefficients, corresponding formula is:

    where y is the vector of the corrected phenotype, μ is the phenotype mean, α is a vector of the fixed effect including sex, populations and slaughter batches, W is the corresponding indicator matrix for α, ε is the residual effect following a multivariate normal distribution, τ-1is the variance of the residual errors, Inis an n-vector of 1s,nis the number of phenotypic individuals.Z is the genotype indicator matrix; ~β is the SNP substitution effect vector come from a mixture of two normal distributions:

    where σ2a/pτ is the variance for the SNPs with large effects, σ2b/pτ is the variance for the SNPs with minor effects, p is the number of SNPs, and π denotes the proportion of SNPs with large effects.SNP effect ~β was estimated by GEMMA software (version 0.98.1) [20] uses the Markov chain Monte Carlo (MCMC) algorithm and the

    Evaluation of the accuracy of genomic prediction

    The accuracy of genomic predictions was evaluated using the fivefold cross-validation method and leave-one-out method.In the fivefold cross-validation method, the population (combined population, YK population or LY population) were divided into five equal groups.For each test, one group of individuals served as the validation dataset, while the other four groups constituted the reference dataset.The test was repeated until all individuals had predicted GEBV, and then the prediction accuracy was calculated.In leave-one-out method, the main step is to take one individual out as the verification group each time, and the remaining individuals as the reference group.The test was repeated to circularly predict the GEBV of all individuals, and calculate the prediction accuracy.

    The prediction accuracy was calculated using the formula proposed by Hayes et al.[31], the formula is as follows:

    where A is the prediction accuracy, yvalis the adjusted phenotype of each animal, GEBV is the genomic estimated breeding values, and h2is the heritability of the trait.Estimates of heritability (Table S2) for all traits refer to our previous studies [32].Apand Awdenote, respectively, the prediction accuracy of GEBV for the proportion and weight of pork cuts.To further investigate the genomic prediction accuracy impacted by pre-selection of SNPs, we perform GWAS analysis on the reference dataset and selected SNPs which significantly associated with the phenotype to predict the GEBV of individuals in the validation dataset.In the GWAS based on genotype imputation data, we selected the SNPs withP-values < 0.01 to predict GEBV.Considering that the SNPs of microarray genotyping are much less than the imputation data, we selected SNPs withP-value < 0.05 to predict GEBV in SNP Chip data.

    Results

    Summary of HGWAS

    In haplotype-based association studies, we identified a total of 14 QTLs significantly associated with pork cuts and 14 QTLs significantly associated with carcass morphology traits (Table 1).In shoulder cuts, we found three QTLs associated with the proportion of BBS and FLB(Table 1), with the most significant SNP (rs0700815,P= 4.03 × 10-9) associated with the proportion of FLB located at 31,161,760 bp ofSus scrofachromosome (SSC)7.This QTL region contains genes (GRM4,HMGA1,SMIM29,NUDT3andPPARD) associated with body height and limb bone length [33-35].In middle cuts, we identified 6 QTLs significantly associated with the weight and proportion of RI, BF, and MC (Table 1), with the most significant SNP (rs0702042,P= 1.05 × 10-16) associated with the RI proportion located at 97,732,109 bp of SSC7.This QTL region contains theVRTNgene, which has been identified and functionally validated as a causative gene affecting the number of thoracic vertebrae and ribs [36].In leg cuts, we found three QTLs significantly associated with the weight and proportion of HLB and LC (Table 1), with the most significant SNP (rs1705050,P= 6.25 × 10-11) associated with the HLB weight located 188,108 bp downstream ofBMP2gene on SSC17.In carcass morphology traits, we identified 14 QTLs significantly associated with carcass length (SL, OL), length and number of vertebrae (THL, LUL, THN, LUN, SLUL).The two major candidate genes identified in carcass morphology traits affecting carcass length and vertebral length wereVRTNandBMP2.

    Additionally, we detected new QTLs significantly associated with pork cuts, such as a QTL on SSC1 significantly associated with the weight and proportion of BF, with the most significant SNP (rs0700815,P= 1.57 × 10-7) located at 161,408,832 bp and a QTL on SSC5 associated with LC weight, with the most significant SNP (rs0501529,P= 9.35 × 10-7) located at the position of 81,315,221 bp, 460,749 bp away fromIGF1gene.

    Summary of imputation-based IGWAS

    Based on imputed genotype data, we identified a total of 167 QTLs significantly related to pork cuts and carcass morphology traits (Table S3).The majority of QTLs identified by HGWAS were also validated in the IGWAS,comprising 54 QTLs associated with weight of pork cuts,8 QTLs associated with carcass weight, 58 QTLs associated with proportion of pork cuts, and 47 QTLs associated with carcass morphology traits.

    In shoulder cuts, a total of 25 QTLs were identified for the weight of pork cuts and 26 QTLs for the proportion of pork cuts (Table 2, Fig.1, and Table S3).Notably, the largest number of QTLs affecting the weight and proportion of the BBS was observed, with a total of 22 QTLs.The most significant SNP was located at 11,938,089 bp on SSC14, and it was significantly associated with both the weight and proportion of BBS,withP-values of 2.75 × 10-9and 1.58 × 10-9, respectively.This SNP was located in the intronic region of theELP3gene.In middle cuts, a total of 15 and 19 QTLs were identified affecting the weight and proportion of the pork cuts respectively (Table 2, Fig.1, TableS3).The QTLs significantly associated with the weight and proportion of CB were the most.The most significant SNP (rs17_15644200) affecting CB weight was located at 15,644,200 bp on SSC17 with aP-value of 1.70 × 10-9.And rs17_15384749 (15,384,749 bp), located near rs17_15644200, showed a significant association with CB proportion, with aP-value of 2.18 × 10-7.Both SNPs are situated upstream of theBMP2gene.Similarly,the SNPs at positions 97, 130, 183 bp and 97,576,486 bp on SSC7 are top SNPs affecting CB weight and proportion, withP-values of 2.88 × 10-7and 3.92 × 10-9,respectively.These SNPs are located upstream of theVRTNgene at positions 484,524 bp and 38,221 bp.Additionally, two QTLs (97,578,564 - 97,112,240 bp and 97,596,043 - 96,354,619 bp) containing causative gene ofVRTNaffecting vertebra number also significantly affected the weight and proportion of RI.In leg cuts, a total of 14 and 13 QTLs were identified affecting the weight and proportion of the pork cuts respectively.Among them, the greatest number of QTLs that affect the weight and proportion of TL were identified.The most significant SNP affecting the weight of TL was located at 68,490,542 bp on SSC10, with aP-value of 8.07 × 10-8,located in the intronic region of theWDR37.Furthermore, the QTLs significantly associated with HLB weight and proportion were located on SSC17 at 14,621,182-19,590,143 bp and 149,33,905-17,042,539 bp, coveringBMP2.

    Table 1 Significant loci associated with pork cuts and carcass morphology traits by haplotypes-based GWAS

    In carcass morphology traits, 27 QTLs associated with carcass length and vertebral length, 4 QTLs associated with vertebral number, and 16 QTLs associated with the thickness of backfat were detected (Table 3 and Table S3).Among them,VRTNon SSC7 andBMP2on SSC17 were found to be the major QTLs affecting carcass length,vertebrae length, and number of vertebrae (Fig.2).The QTL near toVRTNwas significantly associated with various traits such as carcass SL, OL, THL, LUL, THN,and LUN and QTL nearBMP2was also significantly associated with SL, OL, THL, LUL, and SLUL.Interestingly,VRTNwas found to affect carcass length and total vertebral length by increasing the number of vertebrae,while theBMP2may affect these traits by affecting the length of every vertebra.In backfat thickness traits, the most significant SNP was located at 12,758,893 bp on SSC7 with 38,175 bp upstream of theATXN1gene, which was significantly associated with MBD, with aP-value of 4.05 × 10-8(Table 3).Furthermore, two QTLs affecting the MBD were identified in the region of 159,644-161,160 kb of SSC1 and 7,347-7,356 kb of SSC2, which affect RBD and WBD (Table 3).The most significant SNPs in these two QTLs were rs1_161160798 (SSC1:161,160,798 bp) and rs2_7347710 (SSC2: 7,347,710 bp),located at 386,674 bp downstream ofMC4Rand 158,720 bp downstream ofBATF2gene, and with theP-values of 3.45 × 10-7and 1.23 × 10-7, respectively.

    Table 2 (continued)

    Accuracy of genomic predictions

    The accuracy of GEBV for all traits using SNP Chip data were presented in Table 4.In pork cuts, the highest prediction accuracy was RI (Ap= 0.693, Aw= 0.664), followed by BPS (Ap= 0.665, Aw= 0.640), and the lowest prediction accuracy was TPB (Ap= 0.342, Aw= 0.438).In carcass morphology traits, the highest prediction accuracy was THN (A = 0.882), followed by LHN(A = 0.749), and the lowest prediction accuracy was LUN(A = 0.373) (Table 4).Additionally, pork cuts and carcass morphology traits with the highest prediction accuracy using the GBLUP model were SB (Ap= 0.586, Aw= 0.554)and THL (A = 0.579), respectively (Table 4).Importantly,the accuracy of prediction using the BSLMM model was significantly higher than that of the GBLUP model(P= 9.54 × 10-8) (Fig.3a), with THN showing the greatest improvement of 0.333.Additionally, we found that the prediction accuracy of the leave-one-out method was significantly higher than that of the fivefold cross-validation method (P= 1.27 × 10-10) (Fig.3b).

    Different populations, marker densities and pre-selecting markers

    Fig.2 GWAS results of length-related carcass morphology traits.(left) Manhattan plots for carcass morphology traits with the data after imputation.(right) Quantile-quantile plots (Q-Q plots) for carcass morphology traits.In the Manhattan plots, the y-axis and x-axis represent the -log10(P-value)of the SNPs and the genomic positions separated by chromosomes, respectively.The tomato puree points represent SNPs that exceeded the genome-wide significance threshold (-log10(5 × 10-8)).The green points represent SNPs that exceeded the suggestive significance threshold (-log10(1 × 10-6)).In Q-Q plots, the y-axis and x-axis represent the expected and observed -log10(P-value), respectively

    Table 3 Significant loci associated with carcass morphology traits by imputation-based GWAS

    Table 4 (continued)

    Fig.3 Boxplots comparing the prediction accuracy of GEBV across different models, validation methods, and SNP datasets.a Comparison of the accuracy of predicting GEBV using the GBLUP and BSLMM models based on the CC1 Chip genotyping data.b Comparison of the accuracy of predicting GEBV using the fivefold cross-validation and leave-one-out method based on the CC1 Chip genotyping data.c Comparison of the accuracy of predicting GEBV using the CC1 Chip genotyping data and genotype imputation data based on GBLUP models.CC1 Chip represents CC1 Chip genotyping data.Imputation represents imputation-based of whole-genome sequence

    We found that the prediction accuracy based on the CC1 Chip genotype data was significantly higher than that based on sequence imputation data by GBLUP model(P= 6.16 × 10-5, Fig.3c).This shows that the accuracy of the CC1 chip data for genome selection of pork cuts and carcass morphology traits is better.We propose two potential explanations for this result.Firstly, the CC1 Chip, developed collaboratively by the National Laboratory of Pig Genetic Improvement and Breeding Technology and over 12 universities and research institutes in China, includes causal loci that influence body length and weight.Secondly, the poor prediction accuracy of GEBVs based on genotype imputation data may be attributed to the fact that over 98% of imputed genotypes were not associated (P> 0.05) with the phenotype, and these loci may be unfavorable to the prediction of GEBV.Previous studies have found that the accuracy of GEBV prediction can be improved by excluding markers that have no effect on traits or have inconsistent effects among different populations [37-39].Therefore, we pre-selected a set of SNPs to predict GEBV in WGS genotype imputation data.The prediction accuracy of CC1 Chip data was still significantly higher than that of pre-selected genotype imputation data, but the difference in prediction accuracy became smaller (Fig.4a), with aP-value of 0.027, in the combined populations.However, the prediction accuracy of genotype imputation data was significantly higher than that of chip data in the YK and LY populations(Fig.4a), with significantP-values of 9.88 × 10-24and 1.01 × 10-24, respectively.Similarly, we chose the GWAS significant loci based on CC1 chip data for genomic prediction.The genomic prediction based on the CC1 Chip data showed that the accuracy of GEBVs for different traits in the combined populations using GWAS significant loci was lower than that using all SNPs (Fig.4b),with aP-value of 8.76 × 10-5.However, the prediction accuracy in the YK populations and LY populations was the opposite (Fig.4b).Furthermore, we compared the prediction accuracy under pre-selection strategy of SNP Chip data and imputation data, we found that the prediction accuracy of imputation-based data was significantly higher than that of the CC1 Chip-based data (Fig.5a).The results indicate that the selection of GWAS significant loci for GEBV prediction has substantially improved accuracy in single-breed populations, whether using CC1 Chip data or genotype imputation data.However, in the combined population, the prediction accuracy of GEBVs using all markers from CC1 Chip data outperformed others.Also, the prediction accuracy of GEBVs varies significantly across populations when using different datasets(Fig.5b).

    In summary, when predicting GEBVs using genomewide data, it is advisable to exclude non-relevant loci, also known as pre-selection markers, through GWAS analysis.Different populations may require different strategies for genomic selection.

    Discussion

    Candidate genes affecting body size

    We identified three candidate genes associated with skeletal development, namelyVRTN,BMP2, andHMGA1.A causal mutation (g.19034 A > C) inVRTNwas found to be significantly correlated with thoracic vertebra number in our previous studies, and was confirmed by a series of biochemical experiments [36].In this study,QTLs were also identified in theVRTN, which was significantly associated with the weight and proportion of RI and CB, SL, OL, THL and LUL.Li et al.[40] found that the rs320706814 SNP located approximately 123 kb upstream of theBMP2was the strongest candidate affecting carcass length.However, this study found that the QTL upstream of theBMP2was associated with weight and proportion of FLB, HLB and RI, SL, OL, THL,LUL and SLUL.And, Zhang et al.[35] identifiedHMGA1andPPARDas candidate for limb bone length in pigs in the Large White × Minzhu intercross population.Furthermore, other studies have reported thatHMGA1is a strong candidate gene affecting pig body size [35, 41, 42].This study found that a QTL in the intron region of theHMGA1gene was significantly associated with the proportion of FLB.Overall,VRTN,BMP2, andHMGA1are prominent candidate genes influencing pig body size and play crucial roles in bone development.

    Effects of marker preselection, marker density,and reference population size on genomic prediction

    Based on previous research, we know that several factors can influence the accuracy of predicting genomic estimated breeding values (GEBVs).These include the selection and size of the reference population [43, 44],marker density [45, 46], pre-selection of markers [37],prediction models [47-49], and heritability of traits [50,51].We compared the effect of different populations on GEBV prediction accuracy and observed significantly higher accuracy in the combined populations when using CC1 Chip data compared to the YK and LY populations(Fig.5b).This may be due to the limited size of the YK and LY populations, which reduces the accuracy of GEBV prediction.However, using GWAS significant loci for GEBV prediction resulted in significant improvement in accuracy for the YK and LY populations, although it remained lower than that of the combined population.Apart from the reference population size, the variation in linkage disequilibrium between markers in combined populations and single-breed populations also affects prediction accuracy.In the combined population, linkage disequilibrium blocks formed between markers are smaller.Thus,assuming a specific marker has an effect in the combined population, it is more likely due to its higher linkage disequilibrium with the QTL, rather than longer linkage blocks within a single breed.Previous research by Roos et al.[52] also showed that the accuracy of genome prediction is the highest when multiple populations are combined to form a training set, but a higher labeling density is also required.Higher marker density can improve prediction accuracy to some extent, but not all markers will have an impact on traits.In our study, we found that the accuracy of GEBV predictions using genotype imputation data was lower than that based on CC1 Chip genotyping data.However, when using GWAS significant loci to predict GEBVs of different traits, the accuracy of genotype imputation data significantly improved, and in the singlebreed population, the accuracy of genotype imputation data was significantly higher than that of CC1 Chip data.It can be seen that while increasing the marker density,we also need to pre-select the markers to improve GEBVs prediction accuracy [37-39].

    Fig.5 Boxplot comparing the prediction accuracy of GEBV in different populations.Comparison of the accuracy of predicting GEBV using the CC1 chip-based GWAS significant loci data and imputation-based GWAS significant loci in different populations

    Feasibility of genome-based selection for pork cuts

    As we all know, the most important thing in animal breeding is to select elite individuals and those are identified as the candidates with high EBVs.One of the widely used molecular breeding methods is marker-assisted selection, which involves identifying QTLs associated with traits of interest and then using models incorporating these QTLs to predict EBV in individuals [50].In this study, QTLs related to pork cuts were identified, which has important reference value for breeding pork cuts using marker-assisted selection.However,marker-assisted selection has been gradually replaced by molecular breeding methods based on genomic selection in recent years [53-55].Genomic selection requires establishing a reference population containing phenotype and genotype individuals, evaluating the effect value of each marker on the target phenotype using a suitable model, and then genotyping the individuals that need to be predicted.The GEBVs of each individual are calculated using the estimated marker effect value of the reference population, and individuals are selected and retained based on their GEBVs ranking [56].This method improves the accuracy of selective breeding and shortens the generation interval.It is especially effective for difficult-to-measure phenotypes and phenotypes with low heritability [57, 58].In our previous study, we found that most of the pork cuts were medium to high heritability traits.This suggests that breeding for pork cuts using genomic selection may have higher predictive accuracy.In this study, we predicted the GEBVs of pork cuts weight and proportion and found that the prediction accuracy of pork cuts was similar to that of carcass morphology traits, and the accuracy ranged from 0.342 to 0.693.The prediction accuracy of some pork cuts can even reach above 0.65, such as the proportion of RI and BPS.In addition, the pork cuts are the traits of pigs after slaughter, and it is still challenging to predict the weight and proportion of pork cuts through live bodies.Therefore, the use of genomic selection would be a practical way to select elite pigs for pork cuts early in life.

    Conclusion

    In this study, we identified 14 QTLs and 112 QTLs associated with 17 pork cuts, as well as candidate genes, using HGWAS and IGWAS for the first time.Our results suggest the independent regulation of skeletal development by several genes across different body parts.Specifically,we identifiedHMGA1as a candidate gene that affects the size of the fore leg bones,VRTNas a causal gene that affects the number of vertebral and rib bones andBMP2as candidate gene that affects the size of both hind leg bones and fore leg bones, as well as the length of a single vertebral bone.The QTLs and candidate genes we identified have important implications for marker-assisted selection and genome selection.Moreover, we conducted genomic selection of pork cuts and carcass morphology traits in different populations.We found that the prediction accuracy of GEBVs for pork cuts ranged from 0.342 to 0.693, and that the predictive accuracy of several traits,including ribs, boneless picnic shoulder, tenderloin, hind leg bones, and scapula bones, exceeded 0.6.We also found that genomic selection strategy of using BSLMM model, with higher density of effective markers and preselecting markers can improve the accuracy of GEBVs.Furthermore, we constructed the first reference populations for genome selection of pork cuts in pigs.These reference populations contain the genetic information of main commercial breeds of Landrace, Yorkshire, and Duroc, which can be directly used for genome selection for most of the commercial pig companies.Overall, our study provides valuable insights into the genetics of pork cuts in pigs and lays a foundation for improving the effi-ciency of pig breeding programs.

    Abbreviations

    BBS Boneless boston shoulder

    BE Belly

    BF Back fat

    BL Boneless leg

    BPS Boneless picnic shoulder

    CB Chine bones

    DLY Duroc × Landrace × Yorkshire hybrid

    EBV Estimated breed value

    FLB Fore leg bones

    FR Front rib

    GEBV Genomic estimated breed value

    GWAS Genome-wide association study

    HBD Hip backfat depth

    HGWAS Haplotype-based of CC1 Chip genotyping data GWAS

    HLB Hind leg bones

    IGWAS Imputation-based of whole-genome sequence GWAS

    LC Leg cut

    LR Landrace

    LO Loin

    LUL Lumbar length

    LUN Lumbar number

    LY Landrace and Yorkshire hybrid

    MBD Mean of backfat depth

    MC Middle cut

    MCP Meat cut proportion

    OL Oblique length

    QTL Quantitative trait locus

    RBD 6th_7th rib backfat depth

    RI Ribs

    SB Scapula bone

    SBD Shoulder backfat depth

    SC Shoulder cut

    SL Straight length

    SLUL Single lumbar length

    SNP Single nucleotide polymorphism

    SSC Sus scrofa chromosome

    THL Thoracic length

    THN Thoracic number

    TL Tenderloin

    TPB Tail and pelvis bone

    WBD Waist backfat depth

    YK Yorkshire

    Supplementary Information

    The online version contains supplementary material available at https:// doi.org/ 10.1186/ s40104- 023- 00914-4.

    Additional file 1: Table S1.Summary information for four populations.

    Additional file 2: Table S2.Estimates of heritabilities with their standard error (SE) for pork cuts and carcass morphology traits in the combined populations.

    Additional file 3: Table S3.Significant loci associated with pork cuts and carcass morphology traits by GWAS after genotype imputation.

    Acknowledgements

    The authors would like to acknowledge National Natural Science Foundation of China [grant number 32160782] for their financial support.Special thanks were given to Muyuan Food Co., Ltd.(Henan, China) and Longda Muyuan Meat Food Co., Ltd.(Henan, China) for providing the swine and meat samples.Authors’ contributions

    LH conceived and designed the study and revised the manuscript.LX performed most of the analysis, and was a major contributor in writing the manuscript.ZZ analyzed the data and wrote and revised the manuscript.JQ assisted in sample collection and participated in discussion.LR, DC, XT and LC assisted in sample collection.SX, participated in the designed the study.LH,ZZ, and LX participated in the discussion of the results.All authors read and approved the final manuscript.

    Funding

    This work was supported by National Natural Science Foundation of China[grant number 32160782].

    Availability of data and materials

    The data that support the findings of this study are available from the corresponding author upon reasonable request.

    Declarations

    Ethics approval and consent to participate

    All procedures involving animals followed the guidelines for the care and use of experimental animals (GB/T 27416-2014, Laboratory animal institutionsgeneral requirements for quality and competence) approved by the National Standard of the People’s Republic of China.The ethics committee of Jiangxi Agricultural University specially approved this study.

    Consent for publication

    Not applicable.

    Competing interests

    The authors declare that they have no conflict of interest.

    Received: 11 April 2023 Accepted: 2 July 2023

    我要看日韩黄色一级片| 在线播放国产精品三级| 国产精品亚洲美女久久久| 国产极品精品免费视频能看的| 久久精品国产自在天天线| 国产精品自产拍在线观看55亚洲| 欧美xxxx性猛交bbbb| 成人综合一区亚洲| 身体一侧抽搐| 色综合亚洲欧美另类图片| 偷拍熟女少妇极品色| 可以在线观看的亚洲视频| 日本一本二区三区精品| 国产熟女欧美一区二区| 日日撸夜夜添| 久久午夜福利片| 校园春色视频在线观看| 免费看光身美女| 欧美一区二区精品小视频在线| 欧美黑人欧美精品刺激| 久久国内精品自在自线图片| 国产 一区 欧美 日韩| 国产黄片美女视频| 久久精品夜夜夜夜夜久久蜜豆| 超碰av人人做人人爽久久| 成人性生交大片免费视频hd| 美女被艹到高潮喷水动态| 亚洲一区二区三区色噜噜| 久久精品国产清高在天天线| 午夜老司机福利剧场| 亚洲不卡免费看| 久久精品国产亚洲av涩爱 | 亚洲国产精品合色在线| 欧美bdsm另类| 亚洲国产精品成人综合色| 成年女人毛片免费观看观看9| 尾随美女入室| 亚洲乱码一区二区免费版| 综合色av麻豆| 能在线免费观看的黄片| 小蜜桃在线观看免费完整版高清| 亚洲专区中文字幕在线| 淫秽高清视频在线观看| 麻豆精品久久久久久蜜桃| a级一级毛片免费在线观看| 欧美三级亚洲精品| 国产男人的电影天堂91| 午夜福利视频1000在线观看| 国产色爽女视频免费观看| 在线观看免费视频日本深夜| 亚洲av免费在线观看| 男插女下体视频免费在线播放| 久久久久精品国产欧美久久久| 一个人看的www免费观看视频| 少妇丰满av| 免费观看人在逋| 性欧美人与动物交配| 国产精品不卡视频一区二区| 亚洲自拍偷在线| 深爱激情五月婷婷| 中出人妻视频一区二区| 看十八女毛片水多多多| 久久久成人免费电影| 国产精品久久久久久av不卡| 日本在线视频免费播放| 婷婷精品国产亚洲av| 男女做爰动态图高潮gif福利片| 91久久精品国产一区二区成人| 99热这里只有精品一区| 欧美日韩综合久久久久久 | 18+在线观看网站| 嫩草影视91久久| a级毛片免费高清观看在线播放| 国产精品一区www在线观看 | 免费高清视频大片| 能在线免费观看的黄片| 色5月婷婷丁香| 亚洲乱码一区二区免费版| 观看美女的网站| 久久草成人影院| 看免费成人av毛片| 久久久成人免费电影| 女人十人毛片免费观看3o分钟| 女的被弄到高潮叫床怎么办 | 春色校园在线视频观看| 久久久久久久亚洲中文字幕| 久久久久久九九精品二区国产| 男女边吃奶边做爰视频| 成年免费大片在线观看| 亚洲国产精品合色在线| 国产精品免费一区二区三区在线| 亚洲av免费高清在线观看| 啦啦啦观看免费观看视频高清| 精品一区二区三区视频在线观看免费| 精品久久国产蜜桃| 啦啦啦啦在线视频资源| 亚洲av二区三区四区| 99热网站在线观看| 国产淫片久久久久久久久| 日日摸夜夜添夜夜添小说| 久久精品人妻少妇| 搡老妇女老女人老熟妇| 尤物成人国产欧美一区二区三区| 亚洲一区二区三区色噜噜| 亚洲精品粉嫩美女一区| 日本撒尿小便嘘嘘汇集6| 国产精品亚洲美女久久久| 欧美日韩精品成人综合77777| 色播亚洲综合网| 国产精品永久免费网站| 午夜福利在线观看免费完整高清在 | 欧美人与善性xxx| av国产免费在线观看| 少妇猛男粗大的猛烈进出视频 | 亚洲成a人片在线一区二区| 亚洲av不卡在线观看| 联通29元200g的流量卡| 国产亚洲精品av在线| 黄色日韩在线| 村上凉子中文字幕在线| 男人舔奶头视频| 精品午夜福利在线看| 日日干狠狠操夜夜爽| 99riav亚洲国产免费| 少妇人妻精品综合一区二区 | av在线蜜桃| 国产免费av片在线观看野外av| 搞女人的毛片| 亚洲av二区三区四区| 国产男靠女视频免费网站| 我的女老师完整版在线观看| 偷拍熟女少妇极品色| 亚洲av免费高清在线观看| aaaaa片日本免费| 观看免费一级毛片| 久久久久久久久久成人| 嫩草影视91久久| 国国产精品蜜臀av免费| 中出人妻视频一区二区| 国产一级毛片七仙女欲春2| 色哟哟·www| 免费看av在线观看网站| 三级男女做爰猛烈吃奶摸视频| 国产亚洲91精品色在线| 三级毛片av免费| 两性午夜刺激爽爽歪歪视频在线观看| 日日撸夜夜添| 日本撒尿小便嘘嘘汇集6| 国产精品一区二区三区四区免费观看 | 99热网站在线观看| 非洲黑人性xxxx精品又粗又长| 日本熟妇午夜| 精品乱码久久久久久99久播| 日本黄色片子视频| 国产精品久久久久久久电影| 我要看日韩黄色一级片| 级片在线观看| 日日啪夜夜撸| 中亚洲国语对白在线视频| 天堂网av新在线| 此物有八面人人有两片| 嫩草影院入口| 赤兔流量卡办理| 久久精品国产自在天天线| 午夜福利视频1000在线观看| 伦精品一区二区三区| 一进一出好大好爽视频| 很黄的视频免费| 99久久九九国产精品国产免费| 亚州av有码| 午夜激情福利司机影院| 一个人看的www免费观看视频| 国产亚洲91精品色在线| 午夜福利成人在线免费观看| 成人av一区二区三区在线看| 精品免费久久久久久久清纯| 日韩大尺度精品在线看网址| 五月伊人婷婷丁香| 精品一区二区三区视频在线| 色综合色国产| 色哟哟·www| 精品日产1卡2卡| 一区二区三区激情视频| 一级毛片久久久久久久久女| 国产精品一区二区免费欧美| 久久欧美精品欧美久久欧美| videossex国产| 免费一级毛片在线播放高清视频| 少妇人妻一区二区三区视频| 久久午夜亚洲精品久久| 国产免费男女视频| 长腿黑丝高跟| 亚洲人成网站在线播| 欧美最新免费一区二区三区| 淫秽高清视频在线观看| 精品久久久久久久久久久久久| 亚洲专区中文字幕在线| 97热精品久久久久久| 色综合婷婷激情| 99久久精品热视频| 成年免费大片在线观看| 免费在线观看影片大全网站| 日韩av在线大香蕉| 成人精品一区二区免费| 成人国产麻豆网| 一进一出抽搐动态| 女同久久另类99精品国产91| 俄罗斯特黄特色一大片| 精品久久久久久成人av| 国产精品久久久久久亚洲av鲁大| 亚洲va在线va天堂va国产| 亚洲精华国产精华液的使用体验 | 亚洲成人中文字幕在线播放| 国产成人aa在线观看| 亚洲色图av天堂| 亚州av有码| 日本黄色片子视频| 一进一出抽搐gif免费好疼| 最近中文字幕高清免费大全6 | 久久婷婷人人爽人人干人人爱| 三级毛片av免费| 三级男女做爰猛烈吃奶摸视频| 露出奶头的视频| 高清在线国产一区| 日韩欧美精品v在线| 日本欧美国产在线视频| 亚洲av二区三区四区| 亚洲精品国产成人久久av| 别揉我奶头 嗯啊视频| 男人舔奶头视频| 国产成人影院久久av| 伦精品一区二区三区| 亚洲熟妇中文字幕五十中出| 九九爱精品视频在线观看| 欧美+亚洲+日韩+国产| 久久天躁狠狠躁夜夜2o2o| 午夜日韩欧美国产| 日本精品一区二区三区蜜桃| 国产精品永久免费网站| av国产免费在线观看| 国产精品美女特级片免费视频播放器| 日本黄色视频三级网站网址| 色播亚洲综合网| 国产精品野战在线观看| 国产91精品成人一区二区三区| 免费看av在线观看网站| 亚洲 国产 在线| 日韩欧美精品免费久久| 少妇熟女aⅴ在线视频| 变态另类丝袜制服| 午夜久久久久精精品| 无遮挡黄片免费观看| 亚洲欧美日韩高清专用| 麻豆成人av在线观看| 久久久精品欧美日韩精品| 最近最新中文字幕大全电影3| 久久人人精品亚洲av| 麻豆国产av国片精品| 亚洲色图av天堂| 国产精品自产拍在线观看55亚洲| 99久久精品一区二区三区| 九九久久精品国产亚洲av麻豆| 在线看三级毛片| 亚洲一区二区三区色噜噜| 丝袜美腿在线中文| 国产在线精品亚洲第一网站| 人妻夜夜爽99麻豆av| 国产91精品成人一区二区三区| 男女之事视频高清在线观看| 97超级碰碰碰精品色视频在线观看| 少妇熟女aⅴ在线视频| 色吧在线观看| 久久国内精品自在自线图片| a级一级毛片免费在线观看| 大型黄色视频在线免费观看| 亚洲精品粉嫩美女一区| 一本一本综合久久| 国产色爽女视频免费观看| 午夜福利高清视频| 观看美女的网站| 久久久久久久久久久丰满 | 亚洲熟妇中文字幕五十中出| 国产精品美女特级片免费视频播放器| 韩国av一区二区三区四区| 国产 一区精品| 一个人观看的视频www高清免费观看| 九九爱精品视频在线观看| 永久网站在线| 欧美性猛交黑人性爽| 中出人妻视频一区二区| 熟女电影av网| 国产一区二区在线av高清观看| 色视频www国产| h日本视频在线播放| 五月玫瑰六月丁香| 美女高潮的动态| 亚洲欧美日韩东京热| 久久久久久久午夜电影| 联通29元200g的流量卡| 久久中文看片网| 亚洲avbb在线观看| 欧美bdsm另类| 人妻久久中文字幕网| 99国产极品粉嫩在线观看| 亚洲国产色片| 国产女主播在线喷水免费视频网站 | 欧美色视频一区免费| 亚洲经典国产精华液单| 天天躁日日操中文字幕| 欧美激情久久久久久爽电影| а√天堂www在线а√下载| 中文资源天堂在线| 亚洲性夜色夜夜综合| 免费在线观看日本一区| 久久中文看片网| 日本在线视频免费播放| 日韩欧美一区二区三区在线观看| 男人狂女人下面高潮的视频| 国产 一区精品| 99热这里只有是精品在线观看| 免费在线观看影片大全网站| av女优亚洲男人天堂| 亚洲国产高清在线一区二区三| 嫩草影院入口| 三级男女做爰猛烈吃奶摸视频| 69av精品久久久久久| 一个人看视频在线观看www免费| 亚洲av日韩精品久久久久久密| 精品久久久久久久久久免费视频| 深夜精品福利| 老师上课跳d突然被开到最大视频| 少妇被粗大猛烈的视频| 亚洲精品456在线播放app | 精品久久久久久久久久久久久| 久久精品夜夜夜夜夜久久蜜豆| 99久国产av精品| 日本与韩国留学比较| 午夜激情欧美在线| av中文乱码字幕在线| 国产麻豆成人av免费视频| 免费在线观看影片大全网站| 99精品久久久久人妻精品| 成年女人看的毛片在线观看| 久久九九热精品免费| 中文亚洲av片在线观看爽| 性插视频无遮挡在线免费观看| 国产精品国产高清国产av| 又黄又爽又免费观看的视频| 男女之事视频高清在线观看| 黄色配什么色好看| 狠狠狠狠99中文字幕| 久久久久久久久久成人| 欧美成人性av电影在线观看| 国内精品久久久久久久电影| 人人妻,人人澡人人爽秒播| 1024手机看黄色片| 成人毛片a级毛片在线播放| 99在线视频只有这里精品首页| 成年女人看的毛片在线观看| 精品一区二区三区视频在线| 日韩一本色道免费dvd| 国产真实乱freesex| 亚洲一区二区三区色噜噜| 亚洲av免费高清在线观看| 日本-黄色视频高清免费观看| 一本久久中文字幕| 麻豆久久精品国产亚洲av| 亚洲av免费高清在线观看| 午夜福利在线观看免费完整高清在 | 国产私拍福利视频在线观看| 久久久久久久亚洲中文字幕| 国产黄a三级三级三级人| 观看免费一级毛片| 国产欧美日韩精品亚洲av| 日韩大尺度精品在线看网址| 男女边吃奶边做爰视频| 在线播放无遮挡| 欧美性猛交╳xxx乱大交人| 日韩在线高清观看一区二区三区 | 大又大粗又爽又黄少妇毛片口| 亚洲熟妇中文字幕五十中出| 天堂av国产一区二区熟女人妻| 国产精品日韩av在线免费观看| 天堂av国产一区二区熟女人妻| 久久久久久久午夜电影| 91狼人影院| 深爱激情五月婷婷| 黄色一级大片看看| www.www免费av| 麻豆成人av在线观看| 18禁裸乳无遮挡免费网站照片| 全区人妻精品视频| 国产真实乱freesex| 在线观看av片永久免费下载| 桃色一区二区三区在线观看| 五月伊人婷婷丁香| 久久精品国产鲁丝片午夜精品 | 免费av观看视频| 日韩强制内射视频| 日韩国内少妇激情av| 成人午夜高清在线视频| 国产精品1区2区在线观看.| 亚洲精品久久国产高清桃花| 精品一区二区三区人妻视频| 国产美女午夜福利| 亚洲成人精品中文字幕电影| 午夜福利成人在线免费观看| 久久精品国产清高在天天线| 99热精品在线国产| 黄色视频,在线免费观看| 欧美激情国产日韩精品一区| 欧美激情在线99| 免费av不卡在线播放| 最后的刺客免费高清国语| 动漫黄色视频在线观看| 久99久视频精品免费| 91午夜精品亚洲一区二区三区 | 亚洲综合色惰| 99热精品在线国产| av在线天堂中文字幕| 熟妇人妻久久中文字幕3abv| 99久久无色码亚洲精品果冻| 亚洲国产高清在线一区二区三| 久久久色成人| av天堂在线播放| 免费看av在线观看网站| 久久这里只有精品中国| 久久久精品大字幕| 亚洲精品国产成人久久av| 伦精品一区二区三区| 久久6这里有精品| 在线观看66精品国产| 亚洲欧美清纯卡通| 久久热精品热| 蜜桃亚洲精品一区二区三区| 久久精品国产亚洲av香蕉五月| 久久人妻av系列| 99在线人妻在线中文字幕| 国产精品伦人一区二区| 搡女人真爽免费视频火全软件 | 热99re8久久精品国产| 国产乱人伦免费视频| 九九热线精品视视频播放| 婷婷丁香在线五月| 亚洲乱码一区二区免费版| 精品一区二区三区视频在线| 全区人妻精品视频| 欧美zozozo另类| 中文资源天堂在线| 亚洲性夜色夜夜综合| 欧美最黄视频在线播放免费| 在线观看一区二区三区| 欧美极品一区二区三区四区| 91av网一区二区| 欧美丝袜亚洲另类 | 啪啪无遮挡十八禁网站| 色5月婷婷丁香| 五月玫瑰六月丁香| 中文字幕人妻熟人妻熟丝袜美| 少妇人妻一区二区三区视频| 热99re8久久精品国产| 我的老师免费观看完整版| 国产精华一区二区三区| 又紧又爽又黄一区二区| 欧美在线一区亚洲| 国产成人一区二区在线| 琪琪午夜伦伦电影理论片6080| 国产色婷婷99| 日韩欧美三级三区| 国产精品永久免费网站| 国产一区二区在线观看日韩| 夜夜夜夜夜久久久久| 神马国产精品三级电影在线观看| 成人性生交大片免费视频hd| 日本-黄色视频高清免费观看| 夜夜爽天天搞| 韩国av一区二区三区四区| 欧美激情在线99| 亚洲国产高清在线一区二区三| 51国产日韩欧美| 亚洲精品日韩av片在线观看| 国产色爽女视频免费观看| a级毛片免费高清观看在线播放| 伊人久久精品亚洲午夜| 国产亚洲精品久久久com| 日日啪夜夜撸| 久久精品国产鲁丝片午夜精品 | 亚洲最大成人中文| 亚洲va在线va天堂va国产| 亚洲精品亚洲一区二区| 老司机福利观看| 亚洲精品色激情综合| 91在线观看av| 久久午夜福利片| 日本免费一区二区三区高清不卡| 免费看日本二区| 我的女老师完整版在线观看| 亚洲成av人片在线播放无| 毛片一级片免费看久久久久 | 最近最新免费中文字幕在线| 成人毛片a级毛片在线播放| 老司机午夜福利在线观看视频| 亚洲av.av天堂| 中文字幕av在线有码专区| 日韩av在线大香蕉| 91久久精品电影网| 欧美极品一区二区三区四区| 长腿黑丝高跟| 看免费成人av毛片| 国产淫片久久久久久久久| 亚洲精品日韩av片在线观看| 精品久久久久久久久久免费视频| 在线观看美女被高潮喷水网站| 麻豆成人av在线观看| 午夜激情福利司机影院| 久久99热这里只有精品18| 天堂av国产一区二区熟女人妻| 狠狠狠狠99中文字幕| 999久久久精品免费观看国产| 韩国av在线不卡| 免费av毛片视频| 成人精品一区二区免费| 中亚洲国语对白在线视频| 91午夜精品亚洲一区二区三区 | 天堂动漫精品| 国产精品乱码一区二三区的特点| 真人一进一出gif抽搐免费| 无遮挡黄片免费观看| av在线老鸭窝| 一个人免费在线观看电影| 日本-黄色视频高清免费观看| 欧美3d第一页| 亚洲国产欧洲综合997久久,| 久久精品国产亚洲av涩爱 | 欧美zozozo另类| 91久久精品电影网| 欧美高清性xxxxhd video| 亚洲国产欧美人成| 国语自产精品视频在线第100页| 国产精品美女特级片免费视频播放器| 久久亚洲精品不卡| 美女免费视频网站| 亚洲人成伊人成综合网2020| 九九久久精品国产亚洲av麻豆| 看黄色毛片网站| 毛片一级片免费看久久久久 | 久久九九热精品免费| 国产真实伦视频高清在线观看 | 亚洲精品国产成人久久av| av福利片在线观看| 99热这里只有是精品在线观看| 国产精品免费一区二区三区在线| 伦精品一区二区三区| 日本欧美国产在线视频| 精品久久久久久久久亚洲 | 久久久久久伊人网av| 成人美女网站在线观看视频| 亚洲精品色激情综合| 日日摸夜夜添夜夜添小说| 极品教师在线免费播放| 不卡视频在线观看欧美| 免费观看人在逋| 国产免费男女视频| 亚洲精华国产精华精| 身体一侧抽搐| 美女cb高潮喷水在线观看| 精品一区二区三区人妻视频| 国产美女午夜福利| 亚洲精品亚洲一区二区| 亚洲成a人片在线一区二区| 3wmmmm亚洲av在线观看| 亚洲久久久久久中文字幕| 丰满的人妻完整版| 日韩大尺度精品在线看网址| 露出奶头的视频| 色综合婷婷激情| 免费大片18禁| 天天一区二区日本电影三级| 色精品久久人妻99蜜桃| 亚洲第一区二区三区不卡| 国产精品一区www在线观看 | 美女 人体艺术 gogo| 美女大奶头视频| 国产精品永久免费网站| 哪里可以看免费的av片| 欧美高清成人免费视频www| 亚洲av熟女| 免费观看在线日韩| 国产精品美女特级片免费视频播放器| 国产精品福利在线免费观看| 内地一区二区视频在线| 精品午夜福利在线看| netflix在线观看网站| 岛国在线免费视频观看| 九九在线视频观看精品| 尾随美女入室| 免费大片18禁| 欧美一区二区精品小视频在线| 99久久成人亚洲精品观看| 一级毛片久久久久久久久女| 91久久精品电影网| 一本精品99久久精品77| 国产精品国产三级国产av玫瑰| 不卡视频在线观看欧美| 欧美xxxx性猛交bbbb| 精品乱码久久久久久99久播| 国产精品一区二区三区四区免费观看 | 亚洲成a人片在线一区二区| 国产亚洲91精品色在线| 亚洲专区国产一区二区| 久久精品国产99精品国产亚洲性色| 男女边吃奶边做爰视频| 夜夜看夜夜爽夜夜摸| 精品99又大又爽又粗少妇毛片 | 精品不卡国产一区二区三区|