Large white pig genome prediction method based on improved genome feature principal component regression

Abstract

The invention discloses a large white pig genome prediction method based on improved genome feature principal component regression, and the method defined as GFPCR divides a traditional PCR method into two parts, carries out simple regression on a genome feature part, and carries out PCR analysis on residual SNPs. When the genomic features contain more cause mutation, the model can obtain higher genomic prediction accuracy. Meanwhile, the GFPCR model genome prediction accuracy is higher than that of a PCR model when the day-age character reaches 100 Kg, and the GFPCR further improves the genome prediction accuracy due to the fact that the proportion of mutation in preselected genome characteristics is large.

Description

Technical Field

[0001] The present invention relates to the technical field of genome prediction, and in particular to a method for predicting the genome of Large White pigs based on improved principal component regression of genome features. Background Art

[0002] Improving the accuracy of genomic selection has always been one of the research hotspots in the field of genomic prediction technology. Although genomic selection can use all markers in the whole genome without the need to identify quantitative trait loci (QTL), integrating biological information and applying it to genomic selection may further improve the accuracy of genomic prediction for traits with smaller population sizes or that are more difficult to obtain records. In addition, when genomic features contain more causal mutations, the genomic prediction accuracy of traditional methods such as GBLUP (genomic Best Linear Unbiased Prediction) and PCR (Principal Component Regression) will be affected. Therefore, a method is needed to overcome the problem of decreased genomic prediction accuracy of traditional methods due to the presence of more causal mutations in genomic features. Summary of the invention

[0003] The present invention aims to improve the accuracy of genome prediction and provides a Large White pig genome prediction method based on improved genome feature principal component regression.

[0004] To achieve this object, the present invention adopts the following technical solutions:

[0005] A method for predicting the genome of Large White pigs based on improved principal component regression of genome features is provided, the steps comprising:

[0006] S1, screening genome prediction sample individuals, and performing correction phenotype calculation on the screened sample individuals;

[0007] S2, constructing a genome-wide association analysis model combining genome and pedigree information, and using the genome-wide association analysis model to analyze and obtain genome features containing genome and pedigree information;

[0008] S3, using the improved genomic feature principal component regression method to build a genomic prediction model, and using the genomic features analyzed in step S2 as model input, predicting and outputting genomic prediction results.

[0009] Preferably, the screening of individual samples is based on the condition that the traits meet the standards, and the traits include age at reaching 100 kg body weight, back fat thickness corrected for reaching 100 kg body weight, and total litter size.

[0010] Preferably, the trait reaching standard condition includes genotype data quality control condition and / or pedigree quality control condition, and the genotype data quality control condition is: individual genotype detection rate>90%, and / or marker genotype detection rate>90%, and / or minimum allele frequency>0.05;

[0011] The pedigree quality control condition is: whether the kinship of the individual pair is greater than 0.7, if so, retain the two individual samples in the individual pair, if not, delete the individual pair.

[0012] Preferably, the method for calculating the corrected phenotype for the two traits of age at reaching 100 kg body weight and corrected backfat thickness at reaching 100 kg body weight is expressed by the following expression (1):

[0013] y ijkas =BL i +G j +TYS k +a a +e ijkas (1-1)

[0014] In expression (1-1), y ijkas Corrected backfat thickness for age at 100 kg or body weight at 100 kg; BL i is the strain effect; G j is the gender effect; TYS k To determine the year-season effect; a a For the value of a sport, obey in is the additive genetic variance, A is the pedigree relationship matrix; e ijkas is a random residual, subject to distributed, is the remaining residual variance.

[0015] Preferably, the method for correcting the phenotypic calculation of the total litter size trait is expressed by the following expression (1-2):

[0016] y ijkas =BL i +G j +TYS k +a a +Pe a +e ijkas (1-2)

[0017] In expression (1-2), y ijkas is the total litter size; BL i is the strain effect; P j Parity effect; TYS k To determine the year-season effect; a aFor the value of a sport, obey in is the additive genetic variance, A is the pedigree relationship matrix; Pe a is the permanent environmental effect of individual a, subject to distributed; represents the permanent environmental effect variance; e ijks is a random residual, subject to distributed, represents the remaining residual variance.

[0018] Preferably, in step S2, the constructed genome-wide association analysis model is expressed by the following expression: the genome-wide association analysis model for growth traits is (1-1), and the genome-wide association analysis model for reproductive traits is (1-2). The only difference is that the kinship matrix H used by ssGWAS combines genomic and pedigree information, and

[0019]

[0020] In expression (1-3), A 11 is the submatrix of the pedigree kinship matrix A that associates individuals with genotype information, A 22 is an individual without genotype in matrix A, 12 or A 21 is the submatrix of the A matrix associating the relationship between genotyped individuals and non-genotyped animals, and G is the genomic kinship matrix, which is calculated by the following formula:

[0021] G=0.95G0+0.05A 22 (1-4)

[0022] Among them, G0 is the correction matrix of G, making it consistent with A 22 Compatible;

[0023] The inverse matrix H of the H matrix -1 The calculation is as follows:

[0024]

[0025] Preferably, in step S3, constructing the genome prediction model is expressed by the following expression (1-13):

[0026] y c =1μ+Sd+Zg+e (1-13)

[0027] In expression (1-13), y cis the phenotypic value after the age at which the pig reaches 100 kg body weight, or the corrected backfat thickness or total litter size at 100 kg body weight, μ is the mean of the analysis population, 1 is a vector with 1 elements, S is the genomic feature matrix, d is the corresponding effect vector of the genomic feature, Z is the design matrix constructed according to the genotype value, g is the individual genomic breeding value vector, and Where G is the genomic kinship matrix, is the additive genetic variance, e is the random residual, and distribution, where I is the identity matrix, is the residual variance.

[0028] Preferably, in step S3, the process of performing genome prediction by the genome prediction model comprises the steps of:

[0029] S31, input the genome feature matrix S and the genome feature corresponding effect vector d into the genome prediction model expressed by expression (1-13), and the model predicts the output y c value;

[0030] S32, add genomic features as covariates from y c After the value is split, the remaining polygenes and the remaining residuals constitute the corrected phenotype y * ;

[0031] S33, for the correction phenotype y * Singular value decomposition was performed to obtain the remaining multi-genomic feature effects as the model prediction results.

[0032] Preferably, in step S32, the genomic features are used as covariates from y c The method of dividing the value is:

[0033] y*=y c -1μ-sd

[0034] Preferably, the correction phenotype y formed in step S32 * The following expression (1-14) expresses:

[0035] y * =Zg+e=1μ+Xβ+e (1-14)

[0036] In expression (1-14), μ is the population mean, 1 is a vector with 1 elements, β is the residual multi-genomic feature effect (p×1), and its variance is I represents the identity matrix, represents the variance of the remaining multi-genome characteristic effect, p represents the number of selected remaining multi-genotype group characteristic effects, X is (n×p) the remaining multi-genotype group characteristic effect matrix, n represents the number of individuals with chips, e is the random residual, and obeys distributed.

[0037] The GF_PCR method provided by the present invention divides the traditional PCR method into two parts, performs a simple regression on the genomic feature part, and performs PCR analysis on the remaining SNPs. When the genomic features contain more causative mutations, the model can obtain higher genomic prediction accuracy. In addition, the genomic prediction accuracy of the GF_PCR model for the 100Kg body weight age trait is higher than that of the PCR model, because the proportion of causative mutations in the pre-selected genomic features is large, and GF_PCR further improves the genomic prediction accuracy. In addition, the heritability of the corrected back fat thickness of Large White pigs up to 100kg body weight is the highest, and its genomic prediction accuracy is also the highest. The GBLUP genomic prediction accuracy of the three traits of Large White pigs is higher than that of the traditional method, because the reference population in this embodiment is larger. BRIEF DESCRIPTION OF THE DRAWINGS

[0038] In order to more clearly illustrate the technical solutions of the embodiments of the present invention, the following briefly introduces the drawings required for use in the embodiments of the present invention. Obviously, the drawings described below are only some embodiments of the present invention, and for ordinary technicians in this field, other drawings can be obtained based on these drawings without creative work.

[0039] Figure 1 This is a diagram of the implementation steps of the Large White Pig genome prediction method based on improved genome feature principal component regression provided by an embodiment of the present invention;

[0040] Figure 2 This is a comparison chart of the genomic prediction accuracy of different methods for the age-related traits of Large White pigs reaching 100 kg in weight;

[0041] Figure 3 This is a comparison chart of the genomic prediction accuracy of different methods for the backfat thickness trait of Large White pigs with a body weight of 100 kg;

[0042] Figure 4 This is a comparison chart of the accuracy of genomic prediction of total litter size in Large White pigs using different methods;

[0043] Figure 5 This is the relationship between the prediction ability and the number of principal components drawn by the 5-fold cross validation (PCR_cv) and HAT (Hat) methods for the 100 kg body weight age trait of Large White pigs;

[0044] Figure 6 This is a graph showing the relationship between the prediction ability and the number of principal components of the 5-fold cross validation (PCR_cv) and HAT (Hat) methods for the backfat thickness trait of Large White pigs up to 100 kg body weight;

[0045] Figure 7 This is the relationship between the prediction ability and the number of principal components drawn by the 5-fold cross validation (PCR_cv) and HAT (Hat) methods for the total litter size trait of Large White pigs;

[0046] Figure 8 This is a comparison chart of the number of principal components of the 5-fold cross validation (PCR_cv) and HAT (Hat, red) methods for different traits of Large White pigs;

[0047] Fig. 9 This is a graph showing the relationship between the prediction ability and the number of principal components of the GF_PCR 5-fold cross validation and HAT (Hat) methods for the 100 kg body weight age trait of Large White pigs;

[0048] Fig.10 This is a graph showing the relationship between the prediction ability and the number of principal components of the trait of backfat thickness corrected for Large White pigs up to 100 kg body weight using GF_PCR 5-fold cross validation and HAT (Hat);

[0049] Fig.11 This is the relationship between the prediction ability and the number of principal components of the total litter size trait of Large White pigs drawn by GF_PCR 5-fold cross validation and HAT (Hat);

[0050] Fig.12 This is the relationship between the prediction ability and the number of principal components of the GF_PCR 5-fold cross validation and HAT (Hat) methods for the 100 kg body weight age trait of Large White pigs;

[0051] Fig.13 This is a graph showing the relationship between the prediction ability and the number of principal components of the trait of backfat thickness corrected for Large White pigs up to 100 kg body weight using GF_PCR 5-fold cross validation and HAT (Hat);

[0052] Fig.14 This is the relationship between the prediction ability and the number of principal components of the total litter size trait of Large White pigs drawn by GF_PCR 5-fold cross validation and HAT (Hat);

[0053] Fig.15 This is the relationship diagram of the number of principal components of the GF_PCR 5-fold cross validation and HAT (Hat, red) methods for different traits of Large White pigs;

[0054] Fig.16 This is a schematic diagram of the unbiased genomic prediction of different methods for the age-related traits of Large White pigs reaching 100 kg body weight;

[0055] Fig.17It is a schematic diagram of the unbiased genomic prediction of different methods for the corrected backfat thickness trait of Large White pigs up to 100 kg body weight;

[0056] Fig.18 It is a schematic diagram of the unbiasedness of genomic prediction of total litter size trait of Large White pigs by different methods;

[0057] Fig.19 This is a schematic diagram of the time-consuming statistics of different traits of Large White pigs based on the marker effect method (HAT method, red). DETAILED DESCRIPTION

[0058] The technical solution of the present invention is further described below with reference to the accompanying drawings and through specific implementation methods.

[0059] Among them, the drawings are only used for illustrative explanations, and they only represent schematic diagrams rather than actual pictures, and should not be understood as limitations on this patent; in order to better illustrate the embodiments of the present invention, some parts of the drawings may be omitted, enlarged or reduced, and do not represent the size of the actual product; for those skilled in the art, it is understandable that some well-known structures and their descriptions in the drawings may be omitted.

[0060] The same or similar numbers in the drawings of the embodiments of the present invention correspond to the same or similar parts; in the description of the present invention, it should be understood that if the terms "upper", "lower", "left", "right", "inner", "outer", etc. appear, the orientation or position relationship indicated is based on the orientation or position relationship shown in the drawings, which is only for the convenience of describing the present invention and simplifying the description, rather than indicating or implying that the device or element referred to must have a specific orientation, be constructed and operated in a specific orientation. Therefore, the terms describing the position relationship in the drawings are only used for illustrative purposes and cannot be understood as a limitation on this patent. For ordinary technicians in this field, the specific meanings of the above terms can be understood according to specific circumstances.

[0061] In the description of the present invention, unless otherwise clearly specified and limited, if the term "connection" or the like appears to indicate the connection relationship between components, the term should be understood in a broad sense, for example, it can be a fixed connection, a detachable connection, or an integral connection; it can be a mechanical connection or an electrical connection; it can be a direct connection or an indirect connection through an intermediate medium, it can be the internal connection of two components or the interaction relationship between two components. For ordinary technicians in this field, the specific meanings of the above terms in the present invention can be understood according to specific circumstances.

[0062] In the Large White pig genome prediction method based on improved genomic feature principal component regression provided by the embodiment of the present invention, genomic selection is performed on the growth and reproductive traits of the pig, and the estimation accuracy and unbiasedness of three different heritability traits GEBV (genomic estimated breeding value) of age at 100 kg body weight, back fat thickness corrected at 100 kg body weight and total litter size are compared; the accuracy and unbiasedness of GBLUP (genomic feature BLUP, genomic feature best linear unbiased prediction), BayesA (Bayes A), BayesB (Bayes B), BayesC (Bayes C), Bayes RidgeRegression (Bayes Ridge Regression), Bayes Lasso (Bayes Least Absolute Shrinkage and Selection Operator), BayesR (Bayes R), PCR (Principal component regression) and their derived new methods are compared.

[0063] 1.1 Experimental materials and methods

[0064] 1.1.1 Test materials

[0065] For reaching 100Kg body weight age and reaching 100Kg body weight correction back fat thickness two characters, the present embodiment selects 3290 large white pigs of Shanghai Xinnong, wherein 1836 are from Canadian system, 1454 are from French system, all carry out generation performance measurement when about 100Kg, obtain live body weight and live body back fat thickness, and it is corrected to reach 100Kg body weight age (labeled as "AGE100") and reach 100Kg body weight correction back fat thickness (labeled as "BF100"). Keep reaching 100Kg body weight age 130-220 days, reach 100Kg body weight correction back fat thickness 5-23mm record. Keep reaching 100Kg body weight age and reaching 100Kg body weight correction back fat thickness all have phenotype individual. Canadian system large white and French system large white population are raised in two different fields respectively, and two fields nutrition management level is basically the same, and two strains have no affinity. The pedigrees of the two lines consist of a total of 4,719 Large White pigs. The pedigree records were from 2012 to 2020, and the average traceable generations of individuals with chips and phenotypes are 6.12 generations.

[0066] For the total litter size trait, this example selects 2043 Large White pigs from a Shanghai company, of which 1097 are from the Canadian line and 946 are from the French line. Delete the records of the same individual with the same parity. The two groups of Canadian Large White and French Large White are raised in two different farms, with basically the same nutritional management level in the two farms, and the two strains are not related. The pedigrees of the two strains consist of a total of 3237 Large White pigs, and the pedigree records are recorded from 2012 to 2020. The average number of generations that can be traced back for individuals with chips and phenotypes is 5.48 generations.

[0067] 1.1.2 Genotype data and quality control

[0068] GeneSeek GGP Porcine HD chip was used to genotype the large white pigs. PLINK V1.90 software was used for quality control, and the quality control conditions were as follows: individual genotype detection rate>90%, marker genotype detection rate>90%, minimum allele frequency>0.05, and only autosomal SNP markers were retained. After quality control, 36430 SNPs remained for the two traits of age at 100 kg body weight and back fat thickness corrected for 100 kg body weight, and 3290 large white pigs. 36451 SNPs remained for the total litter size trait, and 2043 large white pigs. Beagle5.1 was used to fill the chip data, and the filling method was briefly described as follows: the SNP genotypes of the individual missing sites were estimated by filling each other within the group. In addition, the purpose of filling the chip data is to estimate the missing genotypes in the chip and then conduct whole genome association analysis and genome prediction.

[0069] 1.1.3 Pedigree compilation and quality control

[0070] The records of individuals born before their parents were deleted, and the pedigree was tested and corrected using the BLUPF90 family program. In addition, the IBS (Identical by state) relationship matrix was constructed using PLINK v1.90 software to check whether the individual pairs with a kinship greater than 0.7 in the IBS relationship matrix were related. If not, the two individuals were deleted.

[0071] 1.1.4 Corrected phenotype calculation

[0072] In order to avoid the problem of repeated use of information leading to increased computational time, the sum of the estimated breeding values ​​(EBV) based on the pedigree and the residuals was used as the corrected phenotype (y c ), and the corrected phenotype is used as the response variable for genome prediction. In this embodiment, the EBV estimation model for the two traits of age at reaching 100 kg body weight and back fat thickness corrected at reaching 100 kg body weight is expressed by the following formula (1-1):

[0073] y ijkas =BL i +G j +TYS k +a a +e ijkas (1-1)

[0074] In expression (1-1), y ijkas Corrected backfat thickness for age at 100 kg or body weight at 100 kg; EL i is the strain effect; G j is the gender effect; TYS k To determine the annual season effect, the seasons are divided according to the month, that is, December to February of the previous year is winter, March to May is spring, June to August is summer, and September to November is autumn; a a For the value of a sport, obey in is the additive genetic variance, A is the pedigree relationship matrix; e ijkas is a random residual, subject to distributed.

[0075] In expression (1-1), BL i The calculation method is briefly described as follows:

[0076] Because there are two different strains in the population, a matrix is ​​constructed with one strain as 1 and the other strain as 2.

[0077] G j The calculation method is briefly described as follows:

[0078] Because all individuals have only two categories, male and female, a matrix is ​​constructed with male as 1 and female as 2.

[0079] TYS k The calculation method is briefly described as follows:

[0080] Because the seasons are divided into four seasons and 5.5 years, there are 22 levels in total.

[0081] a a The calculation method is briefly described as follows:

[0082] is the individual additive effect and also the individual breeding value, which is obtained by the BLUP method.

[0083] The calculation method is briefly described as follows:

[0084] The additive genetic variance was obtained by constructing a mixed linear model and then using the restricted maximum likelihood method.

[0085] The calculation method of A is briefly described as follows:

[0086] Through the pedigree information of the group, the relationship between individuals is determined and a block matrix is ​​constructed.

[0087] e ijks The calculation method is briefly described as follows:

[0088] The remaining residual effect is the part after deducting the line effect, sex effect, year-season effect and additive breeding value from the phenotypic value.

[0089] The calculation method is briefly described as follows:

[0090] The remaining residual variance is also obtained by constructing a mixed linear model and then using the restricted maximum likelihood method.

[0091] In addition, i, j, k, a, and s in expression (1-1) represent the phenotypic values ​​of the i-th variety, j-th sex, k-th season, a-th individual, and s-th individual, respectively.

[0092] The EBV estimation model for total litter size is:

[0093] y ijkas =BL i +P j +TYS k +Za a +Pe a +e ijkas (1-2)

[0094] In expression (1-2), y ijkas is the total litter size; BL i is the strain effect; P j Parity effect; TYS k To determine the annual season effect, the seasons are divided according to the month, that is, December to February of the previous year is winter, March to May is spring, June to August is summer, and September to November is autumn; a a For the value of a sport, obey in is the additive genetic variance, A is the pedigree relationship matrix; Pe a is the permanent environmental effect of individual a, subject to distributed; represents the permanent environmental effect variance; e ijks is a random residual, subject to distributed.

[0095] In expression (1-2), BL i TYS k , Z, a a 、e ijks , The calculation method of BL is the same as that in expression (1) i TYS k , Z, a a 、e ijks , The calculation method is the same. j The calculation method is: divide the data into 7 levels according to the actual number of fetuses. a The calculation method is as follows: the permanent environmental effect is based on the actual individual's phenotypic data of seven parities as repeated data; The calculation method of is: The variance of permanent environmental effect is also obtained by constructing a mixed linear model and then using the restricted maximum likelihood method.

[0096] The corrected phenotype calculation process expressed by expressions (1-1)-(1-2) is calculated by the DMUAI program in the DMU software, and the genetic parameters of the three traits (including additive variance, permanent environmental variance, residual variance and heritability) and the heritability of growth traits are obtained at the same time. Heritability of reproductive traits in represents additive genetic variance.

[0097] 1.1.5 Genome-wide association analysis

[0098] According to the chip SNPs obtained after quality control of growth traits and reproductive traits, PLINK v1.90 software was used for further pruning, with 50 SNPs as a sliding window, 5 SNPs each time, and only one SNP was retained for the paired sites with LD (r2) greater than 0.9. After pruning, 26552 SNPs remained for the growth genotype and 26502 SNPs remained for the reproductive genotype. LD (r2) is an indicator to measure the degree of linkage between SNP markers.

[0099] In order to make better use of all phenotypic data, a one-step GWAS (genome-wide association study) was used for analysis (one-step GWAS is defined as ssGWAS). A total of 13,978 growth trait phenotype records and 21,343 reproductive trait records were included. Phenotypic and pedigree data were renumbered using RENUMF90 software, genetic parameters and heritability were estimated by AIRENLF90 and BLUPF90, and SNP effects were calculated using POLTGSF90.

[0100] The genome-wide association analysis model for growth traits is (1-1), and the genome-wide association analysis model for reproductive traits is (1-2). The only difference is that the kinship matrix H used in ssGWAS combines genomic and pedigree information and assumes

[0101]

[0102] Among them, A 11 is the submatrix of the pedigree kinship matrix A that associates individuals with genotype information, A 22 is an individual without genotype in matrix A, 12 or A 21 is the submatrix of the A matrix associating the relationship between genotyped individuals and non-genotyped animals, and G is the genomic kinship matrix, which is calculated by the following formula:

[0103] G=0.95G0+0.05A 22 (1-4)

[0104] Among them, G0 is the correction matrix of G, making it consistent with A 22 Compatible.

[0105] The inverse matrix H of the H matrix -1 The calculation is as follows:

[0106]

[0107] It should be noted here that the expression of G is as follows:

[0108] where q j is the frequency of the smallest allele, Z is the genotype of the individual (constructed according to the number of genotypes with the smallest allele frequency, represented by 0, 1, and 2 respectively), and m represents the total number of SNPs.

[0109] It should also be noted that in expression (1-4), the technical effect that the value "0.95" can bring is that the proportion of G is 0.95; the technical effect that the value "0.05" can bring is that the proportion of A is 0.05.

[0110] 1.1.6 Genomic prediction model

[0111] 1.1.6.1 Principal Component Regression Model

[0112] The traditional PCR model (principal component regression, PCR) is shown in the following expression (1-6):

[0113] y c =1μ+Xβ+e (1-6)

[0114] Among them, y cis the age at which 100 kg body weight is reached, or the back fat thickness corrected for reaching 100 kg body weight, or the total litter size corrected phenotypic value, μ is the average value of the corrected phenotype of the analysis population, 1 is a vector with 1 elements, β is the (m×1) SNPs effect vector, and its variance I represents the identity matrix, represents the variance of the SNP effect, m represents the total number of SNPs, X is the n×m SNPs matrix, n represents the number of individuals with chips, e is the random residual, and obeys distributed.

[0115] The PCR method requires The conventional method of PCR is to use singular value decomposition (SVD) to decompose X. n×m To break it down:

[0116]

[0117] In expression (1-7), U n×n is a standard orthogonal matrix, containing the eigenvectors of the genome kinship matrix G, thus containing the kinship information between individuals, and U T U=UU T =I;Δ n×n is a real matrix with non-negative diagonal elements, and its diagonal elements are called singular values; V m×n The elements in the matrix contain the linkage disequilibrium structure information between SNPs, and V T V = I. The matrix V contains X T Character vector of X, the matrix containing the LD between SNPs as its r 2 The symbol for square root.

[0118] X T X performs SVD by multiplying the following matrices:

[0119] X T X=(UΔV T ) T (UΔV T )=VΔ T U T UΔV T =VΔ T ΔV T =VDV T (1-8)

[0120] XX T =(UΔV T )(UΔV T ) T =UΔVT VΔ T U T =UΔΔ T U T =UDU T (1-9)

[0121] Where D = Δ T Δ, and Δ T Δ=ΔΔ T

[0122] The component score matrix T is:

[0123]

[0124] Principal component regression is to use the component score matrix T to perform multiple regression on y, that is, y = Tb + e, and the estimated value of its score coefficient is It can be expressed as:

[0125]

[0126] Since the columns of the score matrix are orthogonal, we can get (T T T) -1 =(V T X T XV) -1 =D -1 .therefore, Since D is a diagonal element, the calculation is relatively simple. Substituting T=XV into (1-11) we get:

[0127]

[0128] According to the original model y c =1μ+Xβ+e, based on the least squares estimation, we can get β=(X T X) -1 X T y. Therefore, the effect β of each SNP is expressed as GEBV was calculated as GEBV = Xβ.

[0129] The PCR model was used to predict the genome using a self-written program.

[0130] 1.1.6.2 New method for principal component regression of genomic features

[0131] The genomic feature principal component regression (GF_PCR) model uses genomic features as prior information. The genomic features are used as covariates or fixed effects of the model, and the model is:

[0132] yc =1μ+Sd+g+e (1-13)

[0133] In expression (1-13), y c is the phenotypic value after the age at which a body weight of 100 kg is reached, or the corrected backfat thickness or total litter size is reached after the body weight of 100 kg is reached, μ is the average value of the corrected phenotype of the population, 1 is a vector with 1 elements, S is a genomic feature matrix, which in this embodiment is composed of SNPs explaining different variation ratios identified by GWAS (Genome-wide association study), d is the corresponding effect vector of the genomic feature, g is the individual genomic breeding value vector, and obeys Where G is the genomic kinship matrix, is the additive genetic variance, e is the random residual, and distribution, where I is the identity matrix, is the remaining residual effect variance.

[0134] The genomic feature matrix S in expression (1-13) is expressed as follows:

[0135] Genomic signature matrix obtained by principal component analysis of SNP genotypes.

[0136] The corresponding effect vector d of the genomic feature is expressed as follows:

[0137] Effect sizes relative to genomic features, a column vector.

[0138] The individual genomic breeding value vector g is expressed as follows:

[0139] is the individual additive effect and also the individual breeding value, which is obtained by the BLUP method and is a column vector.

[0140] The genomic kinship matrix G is expressed as follows:

[0141] The relationship matrix of the genome constructed by SNP markers.

[0142] The calculation method is briefly described as follows:

[0143] Additive genetic variance was obtained by restricted maximum likelihood method after constructing mixed linear models.

[0144] After using the genomic features as covariates to partition the corrected phenotypic values, the remaining polygenes and the remaining residuals constitute the new corrected phenotype y * ,Right now:

[0145] y * =Zg+e=1μ+Xβ+e (1-14)

[0146] The method of using genomic features as covariates to decompose the corrected phenotypic values ​​into two parts: residual polygenic and residual variance.

[0147] For the new corrected phenotype y * The principal component regression model in 1.1.6.1 was used to perform singular value decomposition and obtain the SNPs effect.

[0148] In this embodiment, the results obtained by the whole genome association analysis of 1.1.5 define the genomic features. In order to avoid reusing the information of the validation population, only the reference population is used for GWAS. The proportion of variance explained by each SNP is sorted from large to small, and the SNPs sets that can explain 1%, 2% and 3% of the total variation are used as genomic features. In order to avoid reusing SNPs information, in (1-13), the SNPs used to construct the G array are not included in the SNPs set of genomic features.

[0149] The GF_PCR model was used to predict the genome using a self-written program.

[0150] 1.1.6.3GBLUP model

[0151] The GBLUP model (Genomic best linear unbiased prediction) was used to predict genomic estimated breeding values. The growth trait model is as follows:

[0152] y=Xβ+Zg+e (1-15)

[0153] In expression (1-15), y is the phenotypic record (age at 100 kg body weight or corrected back fat thickness at 100 kg body weight), X is the fixed effect design matrix, β is the strain, sex, and year-season effect, and the seasonal division is the same as 1.1.4, Z is the design matrix based on genotype values, g is the individual genomic breeding value vector, and obeys Where G is the genomic kinship matrix, is the additive genetic variance, e is the random residual, and distributed.

[0154] The reproductive traits use the repeatability model, which is as follows:

[0155] y=Xβ+Zg+Pe+e (1-16)

[0156] In expression (1-16), y is the total litter size phenotype record, X is the fixed effect design matrix, β is the strain, parity, and year-season effect, and the seasonal division is the same as 1.1.4, Z is the design matrix based on genotype values, g is the individual genome breeding value vector, and obeys Where G is the genomic kinship matrix, is the additive genetic variance, Pe is the permanent environmental effect, and distribution, e is the random residual, obey distributed.

[0157] The G array construction method is:

[0158]

[0159] In expression (1-17), m represents the total number of SNPs.

[0160] The Z matrix is ​​the genome marker information matrix, expressed as:

[0161]

[0162] Among them, Z ij represents the jth marker of the ith individual, p j In the expression (1-18), the genotype is 11 for the homozygous genotype of the smallest allele, 22 for the homozygous genotype of the relative allele, and 12 for the heterozygous genotype. Another way to express it is 0, 2, and 1.

[0163] The GBLUP model genome predicted individual GEBVs and were calculated by the DMUAI program in the DMU software.

[0164] 1.1.6.4ssGBLUP model

[0165] The ssGBLUP model expression is the same as that of GBLUP, except that g is subject to Where H is defined as:

[0166]

[0167] In expression (1-19), A 11 , A 12 and A 22 G is a submatrix based on the pedigree relationship matrix A. Subscripts 1 and 2 represent whether the individual has a genotype, respectively. Subscript 1 indicates that the individual has a genotype, and subscript 2 indicates that the individual does not have a genotype. w =(1-w)G a +wA 22 In this embodiment, w = 0.1. a represents the genome relationship matrix, A′ represents the inverse matrix of the pedigree relationship matrix of individuals without genotypes. 12 Represents the transposed matrix of the pedigree relationship matrix between individuals with genotypes and individuals without genotypes.

[0168] The individual GEBV predicted by the ssGBLUP model genome was calculated by the DMUAI program in the DMU software.

[0169] 1.1.6.5G ROH BLUP Model

[0170] G ROH The BLUP model expression is the same as GBLUP, except that g obeys Among them G ROH is the ROH-based kinship matrix. The ROH-based kinship matrix is ​​obtained by multiplying the matrix obtained by the segIBD function in the optiSel package of the R language by 2.

[0171] 1.1.6.6 Bayes Ridge Regression, BayesA, BayesB, BayesC, and Bayesian Lasso Models

[0172] The five methods of Bayes Ridge Regression, BayesA, BayesB, BayesC and Bayesian Lasso are used to perform genome prediction using the BGLR package. The BGLR package is an R package for fitting Bayesian generalized linear regression functions that supports continuous, binary and ordered response types.

[0173] For the two growth traits of age at 100 kg body weight or back fat thickness corrected at 100 kg body weight, the models of these five methods are similar to model (1-15), where y, β, and e are defined the same as in model (1-15), but g represents an M×1 SNP effect vector, and Z is an N×M SNP matrix with genotype encoding of 0=AA, 1=AB, and 2=BB. j The prior for and the prior probability π.

[0174] For the total litter size reproductive trait, the models for the five methods were:

[0175] y c =1μ+Zg+e (1-20)

[0176] In expression (1-20), y c is the phenotypic value corrected for total litter size, μ is the population mean, 1 is a vector with 1 elements, e is the same as defined in model (1-15), g represents the M×1 SNP effect vector, and Z is the N×M SNP matrix, whose genotype encoding is 0=AA, 1=AB, 2=BB. j The prior for and the prior probability π.

[0177] In the Bayes Ridge Regression (BRR) model, all SNPs have an effect and obey π is equal to 0. The prior distribution of is a standardized inverse chi-square distribution (x -2 ), whose hyperparameter is the degrees of freedom df β =5, scale parameter S β is the discrete function of the original observation in the sample, according to the model prior R 2 Assignment. In the BGLR package, using default parameters, R 2 = 0.5, i.e. the value of the hyperparameter is set to 50% of the variation as prior for linear prediction and 50% as residual in the model.

[0178] In the BayesA model, all SNPs have an effect and obey π is equal to 0. represents the effect variance of the j-th SNP, with v degrees of freedom and a scale parameter of The inverse chi-square distribution of .

[0179] The BayesB model has a different prior assumption from the BayesA model. Its SNPs effects obey the two-component joint distribution. It assumes that a few SNPs have effects and the effects obey Most SNPs have no effect, i.e.:

[0180]

[0181] In expression (1-21), j = 1, ..., p, π represents the proportion of SNPs with no effect,

[0182] In the BayesC model, the effects of SNPs have the same variance, represents the effect variance of the j-th SNP, with v degrees of freedom and a scale parameter of The inverse chi-square distribution of . Therefore, the SNPs effect with probability (1-π) follows the multivariate mixed Student t distribution The parameter π follows the π~beta(p0,π0) prior distribution with p0>0 and π0∈[0,1], where In the BGLR package, the default parameters are used, setting the parameters p0=0.5, π0=10.

[0183] The Bayesian Lasso model is similar to model (1-15). The prior distribution of the SNPs effect is assumed to be a double exponential distribution with parameter λ 2 Unknown, subject to prior distribution λ 2~gamma(r,s). In the BGLR package, the default parameter r=1.1 is used, and the scale parameter s is based on the model's prior R 2 calculate.

[0184] For the above five Bayes methods, Gibbs sampling iterations were set to 25,000, the first 5,000 were set to burn-in, and thin was set to 20. GEBV was calculated as GEBV i =∑Z ij g j .

[0185] 1.1.6.7 BayesR model

[0186] The BayesR model looks like this:

[0187] y c =1μ+Zg+e (1-22)

[0188] In expression (1-22), y c is the phenotypic value of the age at reaching 100 kg body weight, the back fat thickness corrected for reaching 100 kg body weight, and the total litter size after correction, μ is the population mean, 1 is a vector with 1 elements, g represents the M×1 SNPs effect vector, and the SNPs effect obeys Effect variance of SNPs There are 4 possible values,

[0189] Right now

[0190] Therefore, the effect of each SNP follows a mixed distribution composed of four normal distributions, which is in is the total genetic variance, Z is the N×M SNPs matrix (Erbe et al., 2012), e is the random residual, and follows distributed.

[0191] BayesR used Gibbs sampling to construct the posterior distribution, the MCMC chain length was set to 30,000, and the first 10,000 times were set as burn-in, and every 20 samples were used to estimate the SNPs effects and variance components. GEBV was calculated as GEBV i =∑Z ij g j .

[0192] The BayesR model for genome prediction was implemented using the BayesR Fortran program.

[0193] 1.1.7 Cross-validation of genome prediction

[0194] For the growth trait data set (age at reaching 100 kg body weight, corrected back fat thickness at reaching 100 kg body weight), the group was divided into a grandparent group and an offspring group with 2020.01.01 as the boundary, including 2648 Large White pigs in the grandparent group and 642 Large White pigs in the offspring group, and there was no parent-child relationship among the individuals in the offspring group; for the total litter size trait, the group was divided into a grandparent group and an offspring group with 2019.01.01 as the boundary, including 1514 Large White pigs in the grandparent group and 529 Large White pigs in the offspring group, and there was no parent-child relationship among the individuals in the offspring group.

[0195] Eight methods, including GBLUP, Bayes Ridge Regression, BayesA, BayesB, BayesC, Bayesian Lasso, PCR and GF_PCR, were used to predict the three traits of age at 100 kg body weight, back fat thickness corrected at 100 kg body weight and total litter size. The accuracy and unbiasedness were verified by cross-validation. In order to simulate the actual production situation as much as possible, the offspring population was randomly divided into 5 parts during cross-validation, 4 of which and the grandparent population constituted the reference population, and the other was used as the validation population to obtain the GEBV estimate of each individual.

[0196] Accuracy and unbiasedness were used as indicators of genomic prediction: accuracy was the correlation coefficient between the corrected phenotypic value and the individual estimated GEBV; unbiasedness was the regression coefficient of the corrected phenotypic value on the individual estimated GEBV. The cross-validation process for each trait of each GEBV estimation method was repeated 20 times, and the average of the 20 correlation values ​​was taken as the final accuracy.

[0197] In addition, for the PCR method and GF_PCR method, the cross-validation method and HAT method were used to calculate the accuracy of GEBV. The phenotypic value predicted by the HAT method was:

[0198]

[0199] in, Obtained from (1-12). Therefore, the HAT matrix can be defined as

[0200] H=TD -1 T T (1-24)

[0201] Let the residual r=yy=(IH)y, we can get:

[0202]

[0203] Among them, h jj is the jth diagonal element of the HAT matrix H. Let the total variance Since y has been standardized, SS = n-1. Therefore, the prediction reliability is:

[0204]

[0205] The HAT method can replace cross-validation to obtain genomic prediction accuracy and can determine the choice of the optimal number of PCs.

[0206] 1.2 Test results

[0207] 1.2.1 Basic descriptive statistics of phenotypic data

[0208] The basic descriptive statistics of the two strains of Large White pigs reaching 100 kg in age and 100 kg in corrected back fat thickness are shown in Table 1-1, and the basic descriptive statistics of the total number of piglets born by the two strains of Large White pigs are shown in Table 1-2.

[0209] Table 1-1 Basic statistics of age at reaching 100 kg body weight and corrected back fat thickness at reaching 100 kg body weight

[0210]

[0211] Table 1-2 Basic statistics of total litter size

[0212]

[0213] In the third column of Table 1-1, "day" means "day" and "mm" means millimeter.

[0214] 1.2.2 Genetic parameter estimation

[0215] The variance components and heritability of the age at 100kg body weight, backfat thickness at 100kg body weight correction and total litter size of Large White pigs are shown in Table 1-3. The heritability of the three traits of Large White pigs ranges from 0.111 to 0.491, among which the heritability of total litter size is the lowest, at 0.111. The heritability of the age at 100kg body weight and backfat thickness at 100kg body weight correction of Large White pigs is similar.

[0216] Table 1-3 Variance components and heritability of the age at which Large White pigs reach 100 kg body weight, corrected backfat thickness at 100 kg body weight, and total litter size

[0217]

[0218] 1.2.3 Genomic prediction accuracy

[0219] 1.2.3.1 Accuracy of genome prediction by GBLUP, Bayes and PCR methods

[0220] For the three traits of age at 100 kg body weight, corrected backfat thickness at 100 kg body weight, and total litter size, the genomic prediction accuracy obtained by PCR was the highest ( Figure 2 , Figure 3 , Figure 4 ).

[0221] For the two growth traits of age at 100kg body weight and backfat thickness corrected for 100kg body weight, the genomic prediction accuracy of the five Bayes methods, Bayes RidgeRegression, BayesA, BayesB, Bayes, and Bayesian Lasso, was similar, and the accuracy was lower among the eight genomic prediction methods. For the age at 100kg body weight trait, the genomic prediction accuracy of the BayesR method was slightly higher than that of the GBLUP method, while for the backfat thickness corrected for 100kg body weight trait, the genomic prediction accuracy of the BayesR method was slightly lower than that of the GBLUP method. Overall, the genomic prediction accuracy of backfat thickness corrected for 100kg body weight was higher than that of age at 100kg body weight ( Figure 2 , Figure 3 , Appendix Table B-1, Appendix Table B-2).

[0222] Table B-1 Genomic prediction accuracy of different methods for the age of Large White pigs reaching 100 kg body weight

[0223]

[0224]

[0225] Table B-2 Genomic prediction accuracy of different methods for the corrected backfat thickness trait of Large White pigs reaching 100 kg body weight

[0226]

[0227]

[0228] In the total litter size trait, the genomic prediction accuracy of the five Bayes methods, Bayes Ridge Regression, BayesA, BayesB, Bayes, and BayesianLasso, is also similar, but its prediction accuracy is higher than that of the GBLUP method among the eight genomic prediction methods, but its average prediction accuracy is lower than that of the PCR method, and its variance is greater than that of the GBLUP method and the PCR method. The genomic prediction accuracy of the BayesR method is much lower than that of the other eight genomic prediction methods ( Figure 4 , Appendix Table B-3).

[0229] Table B-3 Accuracy of genomic prediction of total litter size in Large White pigs by different methods

[0230]

[0231]

[0232] In the PCR method, for the three traits of Large White pigs, age at 100 kg body weight, corrected backfat thickness at 100 kg body weight, and total litter size, the accuracy first increased, then stabilized, and finally decreased with the increase of the number of principal components, regardless of the cross-validation method or the HAT method ( Figure 5 , Figure 6 , Figure 7 ). For the two traits of age at 100 kg body weight and total litter size in Large White pigs, the accuracy of the HAT method was higher than that of the 5-fold cross-validation method, especially for the age at 100 kg body weight trait; for the trait of corrected backfat thickness at 100 kg body weight, the accuracy of the HAT method was similar to that of the 5-fold cross-validation method, especially when the number of principal components used was greater than 500.

[0233] The genomic prediction accuracy of the HAT method is not lower than that of the 5-fold cross-validation method. When the genomic prediction accuracy reaches the highest, the number of principal components selected by the HAT method and the 5-fold cross-validation method is basically the same in the two traits of age at 100 kg body weight and total litter size. In the trait of back fat thickness corrected for 100 kg body weight, the number of principal components selected by the HAT method is slightly higher than that of the 5-fold cross-validation method, which is 598. However, in the 5-fold cross-validation method, when the number of principal components is selected from 300 to 600, the genomic prediction accuracy is between 0.798 and 0.807 ( Figure 8 ).

[0234] 1.2.3.2 Genomic feature PCR new method and GFBLUP (Genomic feature BLUP, genomic feature best linear unbiased prediction) method genomic prediction accuracy

[0235] The three traits of age at 100 kg body weight, back fat thickness at 100 kg body weight correction and total litter size in Large White pigs were predicted using a new genomic feature PCR (GF_PCR) and GFBLUP method and their accuracy was obtained ( Figure 2 , Figure 3 , Figure 4 ).

[0236] For the age trait of Large White pigs reaching 100 kg in weight, according to the variation explanation ratio obtained from the GWAS analysis results, SNPs with a variation explanation ratio of 1%, 2% or 3% were selected as genomic features. The genomic prediction accuracy of the new GF_PCR method was higher than that of the PCR method and the GFBLUP method. There was no difference in the genomic prediction accuracy when the genomic features were used as fixed effects or covariates. When the SNPs that explained 2% of the total variation were used as fixed effects or covariates, the genomic prediction accuracy of the new GF_PCR method was the highest; when the SNPs that explained 3% of the total variation were used as fixed effects or covariates, the genomic prediction accuracy of the new GF_PCR method decreased ( Figure 2 , Appendix Table B-4).

[0237] Table B-4 Genomic prediction accuracy of different variation ratios of traits of Large White pigs reaching 100 kg body weight at age (number of principal components)

[0238]

[0239] For the trait of backfat thickness corrected for Large White pigs up to 100 kg body weight, SNPs with a variation explanation ratio of 1%, 2% or 3% were selected as genomic features according to the variation explanation ratio obtained from the GWAS analysis results. When genomic features were used as fixed effects, SNPs that explained 1% of the overall variation were selected, and the genomic prediction accuracy of the new GF_PCR method was basically the same as that of the PCR method, both of which were higher than the GFBLUP method. As the number of SNPs explaining the overall variation ratio increased, the accuracy of the new GF_PCR method decreased, but its genomic prediction accuracy was still higher than that of the GFBLUP method. When genomic features were used as covariates, the genomic prediction accuracy of the new GF_PCR method was higher than that of the genomic features as fixed effects. When SNPs that explained 1% or 2% of the overall variation were selected, the genomic prediction accuracy of the new GF_PCR method was higher than that of the PCR method; when SNPs that explained 3% of the overall variation were selected, the genomic prediction accuracy of the new GF_PCR method was basically the same as that of the PCR method; as the number of SNPs explaining the overall variation ratio increased, the accuracy of the new GF_PCR method decreased, and this trend was the same as when genomic features were used as fixed effects ( Figure 3 , Appendix Table B-5). Table B-5 Genomic prediction accuracy of different variation ratios of the corrected backfat thickness trait of Large White pigs reaching 100 kg body weight (number of principal components)

[0240]

[0241] For the total litter size trait of Large White pigs, according to the variation explanation ratio obtained from the GWAS analysis results, SNPs with a variation explanation ratio of 1%, 2% or 3% were selected as genomic features. The genomic prediction accuracy of the new GF_PCR method was slightly lower than that of the PCR method, and the GFBLUP method had the highest genomic prediction accuracy. When genomic features were used as fixed effects or covariates, the accuracy of the new GF_PCR method decreased with the increase of the SNPs explaining the overall variation ratio. The genomic prediction accuracy trends of the two methods were the same ( Figure 4 , Appendix Table B-6).

[0242] Table B-6 Genomic prediction accuracy of different variation ratios of total litter size in Large White pigs (number of principal components)

[0243]

[0244]

[0245] In the PCR part of the new GF_PCR method, the accuracy of the three traits of Large White pigs reaching 100 kg body weight, 100 kg corrected back fat thickness and total litter size, whether it is the cross-validation method or the HAT method, with the increase of the number of principal components, first increased, then stabilized, and finally decreased ( Figures 9 to 14 ). When genomic features were used as fixed effects, in the two traits of age at 100 kg body weight and total litter size, when the number of principal components was less than 1000, the accuracy of the HAT method was basically the same as that of the 5-fold cross-validation method. In the trait of backfat thickness corrected at 100 kg body weight, the accuracy of the HAT method was lower than that of the 5-fold cross-validation method. When genomic features were used as covariates, the accuracy of the HAT method was basically the same as that of the 5-fold cross-validation method in the three traits.

[0246] When the genomic prediction accuracy reached the highest and the genomic features were used as fixed effects, the number of principal components selected by the HAT method and the 5-fold cross-validation method was basically the same in the three traits; when the genomic features were used as covariates, the number of principal components selected by the HAT method and the 5-fold cross-validation method was basically the same in the two traits of age at 100 kg body weight and total litter size. In the trait of back fat thickness corrected for 100 kg body weight, the number of principal components selected by the HAT method was slightly higher than that of the 5-fold cross-validation method, both of which were 598. However, in the 5-fold cross-validation method, when the number of principal components was 450 to 600, the genomic prediction accuracy was between 0.808 and 0.810. As the proportion of variation that SNPs can explain increases, the standard deviation of the number of principal components required to achieve the maximum genomic prediction accuracy for the three traits in 20 5-fold cross-validation methods gradually decreases ( Fig.15 ).

[0247] 1.2.4 Unbiased Genomic Prediction

[0248] 1.2.4.1 Unbiased genome prediction using GBLUP, GFBLUP, Bayes and PCR methods

[0249] For the two growth traits of age at 100 kg and backfat thickness at 100 kg, the five Bayes methods, Bayes RidgeRegression, BayesA, BayesB, Bayes and Bayesian Lasso, had similar unbiased genomic predictions, and the bias was larger among the eight genomic prediction methods. For age at 100 kg, the PCR method had the best unbiased genomic prediction, and its regression coefficient was closest to 1, followed by the GBLUP method, GFBLUP method and BayesR method ( Fig.16 , Appendix Table B-7); in the trait of back fat thickness corrected for body weight up to 100 kg, the unbiasedness of genomic prediction by BayesR method was closest to 1, followed by GBLUP method, GFBLUP method and PCR method ( Fig.17 , Appendix Table B-8).

[0250] Table B-7 Unbiased prediction of different methods of genomic traits for the age of Large White pigs reaching 100 kg body weight

[0251]

[0252]

[0253] Table B-8 Unbiased prediction of the trait of corrected backfat thickness in Large White pigs up to 100 kg by different genomic methods

[0254]

[0255]

[0256] In the total litter size trait, the PCR method has the best genomic unbiasedness, which is 1.192±0.045, followed by the GFBLUP method. The five Bayes methods, Bayes Ridge Regression, BayesA, BayesB, Bayes and Bayesian Lasso, have similar genomic prediction unbiasedness, but their prediction unbiasedness range is also larger among the eight genomic prediction methods, and their regression coefficients are between 0.703 and 0.775. The unbiasedness of the GBLUP method is slightly lower than the above six methods, and its standard deviation is smaller, which is 0.013. The unbiasedness of the BayesR method genomic prediction is much lower than the unbiasedness of the other eight genomic prediction methods ( Figure 5-18, Appendix Table B-9).

[0257] Table B-9 Unbiased genomic prediction of total litter size in Large White pigs by different methods

[0258]

[0259]

[0260] 1.2.4.2 Unbiased genome prediction using new PCR method for genome features

[0261] The three traits of age at 100 kg, corrected backfat thickness at 100 kg and total litter size were predicted using a new genomic feature PCR (GF_PCR) method. When the genomic prediction accuracy reached the highest, the GEBV of each individual was obtained and its unbiasedness was calculated. The unbiased standard deviation of the genomic prediction of age at 100 kg was the smallest, and the unbiased standard deviation of the genomic prediction of total litter size was the largest ( Fig.16 , Fig.17 , Fig.18 , Appendix Tables B-10 to B-12).

[0262] Table B-10 Unbiasedness of genomic prediction of different variation ratios of traits of Large White pigs reaching 100 kg body weight at age (number of principal components)

[0263]

[0264]

[0265] Table B-11 Unbiasedness of genomic prediction of different variation ratios of the corrected backfat thickness trait of Large White pigs reaching 100 kg body weight (number of principal components)

[0266]

[0267] Table B-12 Unbiasedness of genomic prediction of different variation ratios of total litter size in Large White pigs (number of principal components)

[0268]

[0269]

[0270] For the age trait of Large White pigs reaching 100 kg, based on the variation explanation ratio obtained from the GWAS analysis results, SNPs with a variation explanation ratio of 1%, 2% or 3% were selected as genomic features. The genomic prediction bias of the new GF_PCR method was slightly greater than that of the PCR method. When the genomic features were used as fixed effects or covariates, the difference in the unbiasedness of the genomic prediction was small ( Fig.16 ).

[0271] For the trait of backfat thickness corrected for Large White pigs up to 100 kg body weight, SNPs with a variation explanation ratio of 1%, 2% or 3% were selected as genomic features based on the variation explanation ratio obtained from the GWAS analysis results. When genomic features were used as fixed effects, as the number of SNPs explaining the overall variation ratio increased, the regression coefficient of the new GF_PCR method became closer and closer to 1. SNPs explaining 3% of the overall variation were selected, and the unbiasedness was optimal. When genomic features were used as covariates, the skewness of the genomic prediction of the new GF_PCR method was slightly higher than the skewness of the genomic features as fixed effects ( Fig.17 ).

[0272] For the total litter size trait of Large White pigs, according to the variation explanation ratio obtained from the GWAS analysis results, SNPs with a variation explanation ratio of 1%, 2% or 3% were selected as genomic features. The unbiasedness of the GF_PCR new method was better than that of the PCR method. When the genomic features were used as fixed effects or covariates, as the number of SNPs explaining the overall variation ratio increased, the unbiasedness of the GF_PCR new method changed, and its regression coefficient became closer and closer to 1. The genomic prediction accuracy trends of the two methods were the same ( Figure 5-17 ).

[0273] 1.2.5 Time-consuming statistics of genome prediction

[0274] 80.6% of the BGLR packages used by the Bayes method were written in R language, the PCR method was a self-written R script, and the BayesR and BLUP methods both used Forturn language, so only the Bayes method and PCR method were compared for genome prediction time consumption. The time consumed for calculating the three traits of Large White pigs at 100 kg body weight, 100 kg body weight corrected back fat thickness and total litter size based on the marker effect method was as follows: Fig.19 The results show that the Bayes method takes about the same amount of time, but takes a long time; the PCR method takes the least time among the three traits. The HAT method reduces the cross-validation process and takes the shortest time.

[0275] In summary, the Large White Pig genome prediction method based on improved genome feature principal component regression provided in this embodiment is as follows: Figure 1 As shown, the steps include:

[0276] S1, screening genome prediction sample individuals, and performing correction phenotype calculation on the screened sample individuals;

[0277] S2, constructing a genome-wide association analysis model combining genome and pedigree information, and using the genome-wide association analysis model to analyze and obtain genome features containing genome and pedigree information;

[0278] S3, using the improved genomic feature principal component regression method to build a genomic prediction model, and using the genomic features analyzed in step S2 as model input, predicting and outputting genomic prediction results.

[0279] The GF_PCR method provided in this embodiment divides the traditional PCR method into two parts, performs a simple regression on the genome feature part, and performs PCR analysis on the remaining SNPs. When the genome feature contains more causative mutations, the model can obtain higher genome prediction accuracy. However, the accuracy of GF_PCR is affected by the genome feature. If the proportion of causative mutations in the pre-selected genome feature is very small, or even not selected into the genome feature, the accuracy of GF_PCR will be reduced due to the false positive of the genome feature, which may also be the reason why the GF_PCR genome prediction accuracy of the total litter size trait is lower than that of the PCR method. The GF_PCR model genome prediction accuracy of the 100Kg body weight age trait is higher than that of the PCR model, because the proportion of causative mutations in the pre-selected genome feature is large, and GF_PCR further improves the genome prediction accuracy.

[0280] In this example, the heritability of the corrected backfat thickness of Large White pigs up to 100 kg body weight is the highest, and its genomic prediction accuracy is also the highest. The accuracy of the GBLUP genomic prediction of the three traits of Large White pigs is higher than that of the traditional method because the reference population in this example is larger.

[0281] In this embodiment, different GEBV estimation methods were used to conduct genomic prediction studies on individuals with genotypes and phenotypes, and the genomic prediction accuracy and unbiasedness of different estimation methods were calculated. The results show that the genomic prediction accuracy of PCR is higher than the prediction values ​​obtained using the GBLUP model and the Bayes model, and the genomic prediction unbiasedness is generally better than the GBLUP model and the Bayes model. The GF_PCR new method was used to estimate GEBV using the GWAS results as genomic features. As the proportion of genetic variation explained increases, the genomic prediction accuracy of the new method remains unchanged or improves to varying degrees, its skewness is smaller, and the standard deviation of the regression coefficient is reduced. When SNPs that explain 1% of the genetic variation are used as covariates, the GF_PCR new method has the highest genomic prediction accuracy.

[0282] Whether it is the cross-validation method or the HAT method, as the number of principal components increases, the accuracy of the three traits increases first, then stabilizes, and finally decreases. The optimal number of principal components selected for different traits is different. When the number of principal components is less than 1000, the SNPs with different proportions of genetic variation explained are used as fixed effects. In the two traits of age at 100 kg body weight and total litter size of Large White pigs, the number of principal components selected by the HAT method is basically the same as that of the 5-fold cross-validation method of the new GF_PCR method. In the trait of backfat thickness corrected for 100 kg body weight, the accuracy of the HAT method is lower than that of the 5-fold cross-validation method of the new GF_PCR method. When SNPs with the same proportion of genetic variation explained are used as covariates, in the three traits, the number of principal components selected by the HAT method is basically the same as that of the 5-fold cross-validation method of the new GF_PCR method. Therefore, it is possible to consider using the HAT method instead of the cross-validation method to select the appropriate number of principal components.

[0283] It should be noted that the above specific implementations are only preferred embodiments of the present invention and the technical principles used. Those skilled in the art should understand that various modifications, equivalent substitutions, changes, etc. can be made to the present invention. However, as long as these changes do not deviate from the spirit of the present invention, they should be within the scope of protection of the present invention. In addition, some terms used in the specification and claims of this application are not restrictive, but are only for the convenience of description.

Claims

1. A method for predicting the genome of Large White pigs based on improved principal component regression of genome features, characterized in that the steps include: S1, screening genome prediction sample individuals, and performing correction phenotype calculation on the screened sample individuals; S2, constructing a genome-wide association analysis model combining genome and pedigree information, and using the genome-wide association analysis model to analyze and obtain genome features containing genome and pedigree information; S3, using the improved genomic feature principal component regression method to build a genomic prediction model, and using the genomic features analyzed in step S2 as model input, predicting and outputting genomic prediction results.

2. The method for predicting the Large White pig genome based on improved genome feature principal component regression according to claim 1, characterized in that: The screening of sample individuals is based on the condition that the traits meet the standards, including the age at which the body weight reaches 100 kg, the back fat thickness corrected for body weight reaching 100 kg, and the total number of piglets born.

3. The method for predicting the Large White pig genome based on improved genome feature principal component regression according to claim 2, characterized in that: The trait reaching standard conditions include genotype data quality control conditions and / or pedigree quality control conditions, wherein the genotype data quality control conditions are: individual genotype detection rate>90%, and / or marker genotype detection rate>90%, and / or minimum allele frequency>0.05; The pedigree quality control condition is: whether the kinship of the individual pair is greater than 0.7, if so, retain the two individual samples in the individual pair, if not, delete the individual pair.

4. The method for predicting the Large White pig genome based on improved genome feature principal component regression according to claim 1, characterized in that: The method for calculating the corrected phenotype for the two traits of age at reaching 100 kg body weight and back fat thickness at reaching 100 kg body weight is expressed by the following expression (1): y ijkas =BL i +G j +TYS k +a a +e ijkas (1-1) In expression (1-1), y ijkas Corrected backfat thickness for age at 100 kg or body weight at 100 kg; BL i is the strain effect; G j is the gender effect; TYS k To determine the year-season effect; a a For the value of a sport, obey in is the additive genetic variance, A is the pedigree relationship matrix; e ijkas is a random residual, subject to distributed, is the remaining residual variance.

5. The method for predicting the Large White pig genome based on improved genome feature principal component regression according to claim 1, characterized in that: The method for calculating the corrected phenotype for the total litter size trait is expressed by the following expression (1-2): yes ijkas =BL i +G j +TYS k +a a +Pe a +e ijkas (1-2) In expression (1-2), y ijkas is the total litter size; BL i is the strain effect; P j Parity effect; TYS k To determine the year-season effect; a a For the value of a sport, obey in is the additive genetic variance, A is the pedigree relationship matrix; Pe a is the permanent environmental effect of individual a, subject to distributed; represents the permanent environmental effect variance; e ijks is a random residual, subject to distributed, represents the remaining residual variance.

6. The method for predicting the Large White pig genome based on improved genome feature principal component regression according to claim 1, characterized in that: In step S2, the constructed genome-wide association analysis model is expressed by the following expression: The genome-wide association analysis model for growth traits is (1-1), and the genome-wide association analysis model for reproductive traits is (1-2); the only difference is that the kinship matrix H used in ssGWAS combines genomic and pedigree information and assumed In expression (1-3), A 11 is the submatrix of the pedigree kinship matrix A that associates individuals with genotype information, A 22 is an individual without genotype in matrix A, 12 or A 21 is the submatrix of the A matrix associating the relationship between genotyped individuals and non-genotyped animals, and G is the genomic kinship matrix, which is calculated by the following formula: G=0.95G0+0.05A 22 (1-4) Among them, G0 is the correction matrix of G, making it consistent with A 22 Compatible; The inverse matrix H of H matrix -1 The calculation is as follows:

7. The method for predicting the Large White pig genome based on improved genome feature principal component regression according to claim 1, characterized in that: In step S3, the genome prediction model is constructed by the following expression (1-13): and c =1μ+Sd+Zg+e (1-13) In expression (1-13), y c is the phenotypic value after the age at which the pig reaches 100 kg body weight, or the corrected backfat thickness or total litter size at 100 kg body weight, μ is the mean of the analysis population, 1 is a vector with 1 elements, S is the genomic feature matrix, d is the corresponding effect vector of the genomic feature, Z is the design matrix constructed according to the genotype value, g is the individual genomic breeding value vector, and Where G is the genomic kinship matrix, is the additive genetic variance, e is the random residual, and distribution, where I is the identity matrix, is the residual variance.

8. The method for predicting the Large White pig genome based on improved genome feature principal component regression according to claim 7, characterized in that: In step S3, the process of performing genome prediction by the genome prediction model includes the following steps: S31, input the genome feature matrix S and the genome feature corresponding effect vector d into the genome prediction model expressed by expression (1-13), and the model predicts the output y c value; S32, add genomic features as covariates from y c After the value is split, the remaining polygenes and the remaining residuals constitute the corrected phenotype y * ; S33, for the correction phenotype y * Singular value decomposition was performed to obtain the remaining multi-genomic feature effects as the model prediction results.

9. The method for predicting the Large White pig genome based on improved genome feature principal component regression according to claim 8, characterized in that: In step S32, the genomic features are used as covariates from y c The method of dividing the value is: y*=y c -1μ-sd.

10. The method for predicting the Large White pig genome based on improved genome feature principal component regression according to claim 8, characterized in that: The corrected phenotype y constructed in step S32 * The following expression (1-14) expresses: y * =Zg+e=1μ+Xβ+e (1-14) In expression (1-14), μ is the population mean, 1 is a vector with 1 elements, β is the residual multi-genomic feature effect (p×1), and its variance is I represents the identity matrix, represents the variance of the remaining multi-genome characteristic effect, p represents the number of selected remaining multi-genotype group characteristic effects, X is (n×p) the remaining multi-genotype group characteristic effect matrix, n represents the number of individuals with chips, e is the random residual, and obeys distributed.

Images