Genomic selection method for breeding disease-resistant varieties of Tosei-boku

Abstract

The present invention relates to a genome selection method for breeding a disease-resistant variety of Toseisho. In the method, a reference group for disease resistance is generated using juveniles of East Star spotted from non-family groups originating from different geographies, and the East Star spotted juveniles for generating the disease resistance reference group are obtained after screening for disease resistance. A validation group is generated of individuals having a disease resistance phenotype but not included in the reference group to verify the effectiveness of calculating a genomic estimated breeding value using the SNP sites provided by the present invention. These SNP sites include SNP sites related to the disease resistance properties of East Star spotted. The genomic estimated breeding value of the candidate group is calculated using the SNP sites provided by the present invention, and individuals with high genomic estimated breeding value are recognized as having strong disease resistance and are used to breed a disease-resistant variety of East Star spotted. The method for breeding a disease-resistant variety of East Star spotted based on the genomic selection technology provided by the present invention can be used to screen parents with strong disease resistance. The selected parents can be used for selective breeding of East Star spotted. Thus, the present invention provides an efficient technical means for improving the survival rate of a breeding group of East Star spotted.

Description

[Technical field]

[0001] The present invention relates to the technical field of breeding genetic varieties of aquatic products, and in particular to a genome selection method for breeding a disease-resistant variety of Azuma star spot. [Background technology]

[0002] The scientific name of the Eastern starfish is Plectropomus leopardus, which belongs to the order Perciformes, family Sparratidae, subfamily Serraninae, and genus Plectropomus. Eastern starfish has dense and tasty flesh, a bright red color all over its body, and is a valuable marine farmed fish with high economic and ornamental value. Due to high market demand, wild populations have been greatly reduced due to overfishing, making Eastern starfish a rare fish in the market. In recent years, with the advancement and development of artificial farming technology for Eastern starfish, the range and scale of farming have gradually expanded. However, due to the popularity of high-density farming models, as well as the deterioration of the farming environment and the deterioration of germplasm, diseases have frequently occurred during the farming process, causing huge economic losses. Among them, "rot disease" caused by Vibrio harveyi has become a major problem for farmers. For example, in 2011, a bacterial disease caused by Vibrio harveyi occurred in Hainan Province, China, causing a large number of deaths of diseased fish and symptoms such as tail ulcers and whole-body whitening. Usually, aquaculture farmers focus on prevention when dealing with diseases, and treat the disease by spraying antibiotics such as oxytetracycline, enrofloxacin, and norfloxacin in the entire pool after the disease occurs. However, the large-scale use of antibiotics not only affects the quality of aquatic products, but also causes food safety problems, affects the value of commercially available fish, and brings potential damages to consumers. Therefore, it is very necessary to develop a disease-resistant variety of Dong Xing Bin and improve the disease resistance of the Dong Xing Bin aquaculture group itself.

[0003] How to accurately select reproductive individuals is a core issue for breeding superior varieties. In fish disease resistance breeding, disease resistance phenotypes are usually obtained using large-scale established families and artificial infection experiments, and then selective breeding is performed based on the estimated breeding value (EBV) estimated from the best linear unbiased prediction (PBLUP). However, to avoid vertical infection of pathogens, it is common to select and breed uninfected healthy individuals from families with high EBV, rather than using individuals that survived the artificial infection experiment as breeding parents. Therefore, the selection efficiency is not high, and the genetic progress obtained each year is limited. If selectively bred fish cannot establish a family line, selection will be based on phenotypes and experience, which greatly reduces the efficiency of selective breeding. Currently, artificial breeding of Dong Xing mainly relies on natural spawning and natural fertilization, and family lines cannot be established. Therefore, there is a lack of effective selective breeding methods for fish that are difficult to establish family lines, such as Dong Xing. The present invention aims to describe in detail the method for cultivating excellent disease-resistant varieties of Tosei spotted fish by using non-family groups, provide a technical method for cultivating excellent varieties of Tosei spotted fish, and provide new ideas for the selective breeding of fish species that are difficult to establish family lines. Summary of the Invention

[0004] The objective of the present invention is to establish a genomic selection method for breeding excellent disease-resistant varieties of Tohoshinki based on non-family groups, overcome the problems of the need for pedigree establishment and low selection accuracy in traditional selective breeding, provide a molecular breeding method for breeding excellent disease-resistant varieties of Tohoshinki, improve the disease resistance of Tohoshinki cultured groups, and promote the breeding process of excellent disease-resistant varieties of Tohoshinki.

[0005] The genomic selection method for breeding a disease-resistant variety of Tosei mosaic provided by the present invention comprises the steps of: 1) To generate disease resistance reference groups for Tosei mackerel using Tosei mackerel juveniles from different geographies; Furthermore, in order to verify the selection effect by the above method of the present invention, a part of individuals having a disease-resistant phenotype but not included in the reference group is generated as a verification group, wherein the juveniles of East Star Spot for generating the disease-resistant reference group are obtained after screening for disease resistance; 2) performing whole genome resequencing on the reference group and the validation group to capture high-quality single nucleotide polymorphism (SNP) sites within the genome; 3) Comparing two methods, Genomic Best Linear Unbiased Prediction (G-BLUP) and Bayes Cπ, with respect to the tendency of the accuracy of predicting disease resistance trait genomic estimated breeding value (GEBV) to change depending on the number of SNPs, thereby selecting the number of SNPs and the method for calculating the GEBV of the candidate group; 4) selecting a SNP marker containing a SNP site associated with the disease resistance properties of Tonoshiro spot for calculating the GEBV of the Tonoshiro spot candidate group based on the results obtained in step 3); 5) using the validation group to validate the predictive effect of the genome selection method selected in step 3) and the SNP markers selected in step 4), and if the average GEBV of the surviving validation individuals is greater than the average GEBV of the dead individuals, it indicates that there is a positive correlation between GEBV and disease resistance of the individuals, and since individuals with high GEBV also have strong disease resistance, it is possible to breed offspring with strong disease resistance by selecting candidate parents using GEBV; 6) calculating the GEBV of the candidate group using the genome selection method selected in step 3) and the SNP markers selected in step 4), and recognizing individuals with high GEBV as having strong disease resistance and using their offspring as parents for breeding.

[0006] The SNP sites, which are 20,000 (20k) in total, are located at positions 36 and 107 of any of the sequences of SED ID NO:1 to SED ID NO:10000, and the sequences within 35 bp before and after the sites are flanking sequences.

[0007] The SNP sites include 30 SNP sites associated with the disease resistance of Eastern star spot, which are located at positions 36 and 107 of any of the sequences of SED ID NO:1 to SED ID NO:15, and the sequences within 35 bp before and after the sites are flanking sequences.

[0008] The method for breeding disease-resistant varieties of Eastern Star Spot based on the genome selection technology and 20,000 SNP sites provided by the present invention can be used to screen parents with strong disease resistance. The selected parents can be directly used for selective breeding of Eastern Star Spot, so the present invention provides an efficient technical means to improve the survival rate of Eastern Star Spot breeding groups. [Brief description of the drawings]

[0009] [Figure 1] Box plot of common site deletion rates for the Eastern Star reference and validation groups. [Diagram 2] 1 is a graph showing the change in accuracy of predicting the disease resistance property GEBV of eastern star spot by G-BLUP and Bayes Cπ. [Diagram 3] FIG. 11 is a density distribution diagram of 20k SNPs for estimating validation group GEBV. [Figure 4] FIG. 11 is a comparison diagram of the average disease resistance GEBV values ​​of surviving and dead individuals in the validation group. [Diagram 5] FIG. 11 is a comparison of the mean disease-resistant GEBV values ​​of highly disease-resistant individuals and susceptible individuals in the candidate group. DETAILED DESCRIPTION OF THE PREFERRED EMBODIMENTS

[0010] The method for cultivating the anti-Vibrio harveyi Tosei cultured group of the present invention is described in detail below. EXAMPLES

[0011] Example 1: Generating an Eastern Star Spot Resistant Reference Group We collected juveniles of Vibrio harveyi from many large-scale fish farms in Hainan and Shandong provinces of China, transported the collected juveniles to a site where artificial infection experiments can be performed, and artificially infected these juveniles by intraperitoneal injection of Vibrio harveyi. The time of injection of the bacterial solution was set as the starting point (0:00), and the survival status of the juveniles after one injection was observed every 8 hours. Dead individuals were removed at appropriate times, individual information (total length, weight, time of death, etc.) was recorded, and the caudal fin rays were cut off and preserved in absolute ethanol. The experiment was performed for a total of 5 days, and after completion, the caudal fin rays of surviving individuals were collected and preserved in absolute ethanol, and data such as individual origin, total length, weight, and time of death were recorded. Some individuals were selected from them and subjected to whole genome resequencing to generate an anti-Vibrio harveyi reference group of Vibrio harveyi of Tohoshinki (Table 1). The survival status of the Tohoshinki juveniles tested in the artificial infection experiment was used as a disease-resistant phenotype for subsequent analysis. The disease resistance phenotype is defined as a binary trait, where 0 represents juveniles that died during the infection experiment and 1 represents juveniles that survived the experimental period.

[0012] In the experiment, the collected juvenile fish were divided into six groups, and after injection, each was placed into a fiber-reinforced plastic (FRP) tank of the same specifications. The average weight and average total length of the juveniles of A. moniliforme from Shandong Province, China were both higher than those from Hainan Province, China. One-way variance test showed that there was no significant difference in infection survival rate among the six experimental groups (p = 0.05). The survival rate range after infection was 31.85% - 43.41%, and a total of 534 A. moniliforme juveniles were selected to generate the reference group (Table 1).

[0013] [Table 1] EXAMPLES

[0014] Example 2: Genome resequencing of the Eastern Star Spot Resistance Reference Group and screening of high-quality SNPs 1. Whole genome resequencing and mutation detection in the reference group The standard procedure provided by the DNA extraction kit (TIANGEN, Beijing) was used to extract and purify the genomic DNA of the East Star Spot disease-resistant reference group. In addition to the genomic DNA of the reference group, the genomic DNA of 298 East Star Spots with disease-resistant phenotypes was also extracted to verify the effect of the above method of the present invention. The purified genomic DNA was used to generate an Illumina paired-end library, which was then sequenced using the Illumina HiSeq2000 sequencing platform. After filtering out low-quality sequencing reads, these reads were compared with the East Star Spot reference genome using the software BWA. The software GATK was used to detect mutation information (SNPs and INDELs). Finally, in the reference group, 498 individuals were successfully sequenced, the amount of individual sequencing data was 10G, a total of 34,026,719 mutations were obtained, the number of chromosomal mutations was 28,789,387, the number of non-chromosomal mutations was 5,237,332, the average individual sequencing depth was 7.66X, and the mutation results were saved in the file "variances.Ref.vcf.gz". In the validation group, 298 individuals were successfully sequenced, the amount of individual sequencing data was 5G, a total of 12,452,483 mutations were obtained, the number of chromosomal mutations was 11,355,373, the number of non-chromosomal mutations was 1,097,110, the average individual sequencing depth was 3.66X, and the mutation results were saved in the file "variances.Val.vcf.gz".

[0015] 2. Screening for high-quality SNP mutation sites Follow the steps below to filter the above mutations, retain high-quality SNP sites, and then fill in the missing site information. Run the following command in a Linux environment:

[0016] #Use VCFtools to output SNP sites on chromosomes vcftools --gzvcf variances.Ref.vcf.gz --remove-indels --chr chr1 --chr chr2 --chr chr3 --chr chr4 --chr chr5 --chr chr6 --chr chr7 --chr chr8 --chr chr9 --chr chr10 --chr chr11 --chr chr12 --chr chr13 --chr chr14 --chr chr15 --chr chr16 --chr chr17--chr chr18 --chr chr19 --chr chr20 --chr chr21 --chr chr22 --chr chr23 --chr chr24 --recode-INFO-all --recode --stdout | bgzip -c -f > raw.snps.Ref.vcf.gz vcftools --gzvcf variances.Val.vcf.gz --remove-indels --chr chr1 --chr chr2 --chr chr3 --chr chr4 --chr chr5 --chr chr6 --chr chr7 --chr chr8 --chr chr9 --chr chr10 --chr chr11 --chr chr12 --chr chr13 --chr chr14 --chr chr15 --chr chr16 --chr chr17--chr chr18 --chr chr19 --chr chr20 --chr chr21 --chr chr22 --chr chr23 --chr chr24 --recode-INFO-all --recode --stdout | bgzip -c -f > raw.snps.Val.vcf.gz

[0017] After execution, the SNP sites on the chromosomes of the reference and validation groups are saved in the files “raw.snps.Ref.vcf.gz” and “raw.snps.Val.vcf.gz”, respectively.

[0018] #Use PLINK2 to read the above VCF file and reserve biallelic sites plink2 --vcf raw.snps.Ref.vcf.gz --autosome-num 24 --max-alleles 2 --min-alleles 2 --make-bpgen --out raw.Ref plink2 --vcf raw.snps.Val.vcf.gz --autosome-num 24 --max-alleles 2 --min-alleles 2 --make-bpgen --out raw.Val

[0019] #Use R to modify the "raw.Ref.bim" and "raw.Val.bim" files library(data.table) bim <- fread(“raw.Ref.bim”) bim$V2 <- paste(“rs”, bim$V1, ':', bim$V4, sep = “”) fwrite(bim, “raw.Ref.bim”, sep = “\t”, col.names = F, quote = F) bim <- fread(“raw.Val.bim”) bim$V2 <- paste(“rs”, bim$V1, ':', bim$V4, sep = “”) fwrite(bim, “raw.Can.bim”, sep = “\t”, col.names = F, quote = F)

[0020] #Count the deletion rate at each site using PLINK2 plink2 --bpfile raw.Ref --autosome-num 24 --missing --out Ref plink2 --bpfile raw.Val --autosome-num 24 --missing --out Val

[0021] #Use R to count the deletion rate of common sites in two files and output sites with deletion rate less than 0.05 library(data.table) vmiss.ref <- fread(“Ref.vmiss”) vmiss.val <- fread(“Val.vmiss”) comm.loci <- intersect(vmiss.ref$ID, vmiss.val$ID) ref.comm <- vmiss.ref[ID %in% comm.loci] val.comm <- vmiss.val[ID %in% comm.loci] new.vmiss <- data.table(SNP = comm.loci, MISSING = (ref.comm$MISSING_CT + val.comm$MISSING_CT), OBS = (ref.comm$OBS_CT + val.comm$OBS_CT)) new.vmiss$F_MISS <- round(new.vmiss$MISSING / new.vmiss$OBS, 3) kept.loci <- new.vmiss[F_MISS < 0.05, SNP] write.table(kept.loci, “kept.loci.txt”, sep = “\t”, col.names = F, row.names = F, quote = F) After execution, the file "kept.loci.txt" will contain the sites with deletion rates of less than 0.05 in the two files (Figure 1).

[0022] #Merge two files ##Generate VCF files using PLINK2 plink2 --bpfile raw.Ref --autosome-num 24 --extract kept.loci.txt --export vcf-4.2 --out snps.ref plink2 --bpfile raw.Val --autosome-num 24 --extract kept.loci.txt --export vcf-4.2 --out snps.val Compress VCF files bgzip -c -f snps.ref.vcf > snps.ref.vcf.gz bgzip -c -f snps.val.vcf > snps.val.vcf.gz ##Create indexes using bcftools bcftools index -t snps.ref.vcf.gz bcftools index -t snps.val.vcf.gz Merge VCF files using bcftools bcftools merge -m snps -O snps.ref.vcf.gz snps.val.vcf.gz -o snps.merged.vcf.gz

[0023] #Use PLINK2 to perform quality control on the merged SNP sites ##Read the merged file plink2 --vcf snps.merged.vcf.gz --autosome-num 24 --max-alleles 2 --min-alleles 2 --make-bpgen --out snps.merged ##Quality control plink2 --bpfile snps.merged --autosome-num 24 --maf 0.05 --hwe 0.001 --make-bpgen --out snps.qc ##Converting quality-controlled data into VCF plink2 --bpfile snps.qc --autosome-num 24 --export vcf --out snps4impute Compress VCF files bgzip -c -f snps4impute.vcf > snps4impute.vcf.gz ##Filling in missing SNP sites java -Xmx20g -jar beagle.05Apr21.9b7.jar gt=snps4impute.vcf.gz ne=500 iterations=30 burnin=10 out=snps.imputed After the run, a total of 1,213,496 high-quality SNPs are left for subsequent analysis and the filled-in genotypes are saved in the file “snps.imputed.vcf.gz”. EXAMPLES

[0024] Example 3: Two genome selection methods predict the accuracy of disease-resistant GEBV in Eastern spot In order to examine the method and the appropriate number of SNP markers for predicting anti-Vibrio Harveyi GEBV in Eastern star spots, the present invention used the phenotype and genotype data of the generated reference group and compared the change trend of the accuracy of predicting anti-Vibrio Harveyi GEBV in Eastern star spots by G-BLUP and BayesCπ according to the decrease in the number of SNPs by 5-fold cross-validation. The prediction accuracy was evaluated using the area under the subject operating characteristic curve (AUC), which was calculated a total of 25 times, and the final accuracy is shown as the average of the 25 times. Before calculation, the R packages: data.table, ASReml, BGLR, Rcpp, pROC and Parallel must be installed.

[0025] 1. Extract different numbers of SNP sets (SNP subsets) from the Tosei disease resistance reference group Different numbers of SNP subsets are extracted from the high-quality SNPs obtained in Example 2. The individual numbers of the disease-resistant reference group are stored in the file "indiv.ref.txt". Run the following command in a Linux environment to extract the SNP subsets: #Read the VCF file and prepare related data ##Use PLINK2 to read the VCF file and run quality control again plink2 --vcf snps.imputed.vcf.gz --autosome-num 24 --maf 0.05 --make-bpgen --out snps.imputed ##Counting allele frequencies using PLINK2 plink2 --bpfile snps.imputed --autosome-num 24 --freq --out snps.imputed ##Execute the following code in R to prepare the files required for genotype output. library(data.table) freq <- fread(“snps.imputed.afreq”) alleles <- rbind(freq[ALT_FREQS <= 0.5, c(“ID”, “ALT”)], freq[ALT_FREQS > 0.5, c(“ID”, “REF”)], use.names = F) write.table(alleles, “export.A.txt”, sep = “\t”, col.names = F, quote = F)

[0026] After the run, 1,211,259 high-quality SNPs and 796 individuals were left for subsequent analysis, of which 498 individuals were in the reference group and 298 individuals were in the validation group.

[0027] #Using R, extract subsets containing different numbers of SNPs from the reference group and save them in binary format. The numbers of SNPs are 500, 700, 1k, 2k, 8k, 10k, 20k, 50k, 80k, and 100k, respectively. ##Generate genotype AA / Aa / aa files encoded in 0 / 1 / 2 mode using PLINK2 plink2 --bpfile snps.imputed --autosome-num 24 --export-allele export.A.txt --keep indiv.ref.txt --export A --out geno.ref ##Output the individual number awk '{print $2}' geno.ref.raw | sed '1d' > iid.ref.txt ##Run the following script in R to extract the SNP subset and save it in binary format library(data.table) Den.snp <- c(5e2, 7e2, 1e3, 2e3, 8e3, 1e4, 2e4, 5e4, 8e4, 1e5) geno <- fread(paste(“geno.ref.raw”)) geno[, c(1:6):= NULL] geno <- as.matrix(geno) for (i in 1: length(Den.snp)) { coor <- floor(seq(from = 1, to = ncol(geno), length.out = Den.snp[i])) M <- geno[, coor] save(M, file = paste(“ref.snp.s”, i, “.bin”, sep = “”)) } After execution, all SNP subsets are saved in a binary format for easy reading later.

[0028] 2. Evaluate the prediction accuracy of G-BLUP In ASRreml-R, the following models are used to fit the phenotypic data and evaluate the predictive accuracy of G-BLUPs: y ijk =Φ(Ori.Location i +b.weight j +a k ) In the above model, y ijk is the phenotype value, which indicates the phenotype when an individual k collected from location i is infected at body weight j. i is a fixed effect, indicating the location of the subject. j is a covariate and indicates the body weight of the subject at the time of the infection experiment. k is the random additive effect of individuals. Φ denotes the cumulative function of the standard normal distribution. Run the following command in a Linux environment. ##Run the following script in R to generate the inverse of the 3-column G matrix: library(data.table) library(Rcpp) sourceCpp(“matrix_multiplication.cpp”) ginv <- function(invG) { Ginv <- data.frame(row = rep(1 : nrow(invG), nrow(invG)), column = rep(1 : nrow(invG), each = nrow(invG)), value = as.numeric(invG), lower.mat = as.logical(lower.tri(invG, diag = T))) Ginv <- Ginv[Ginv$lower.mat == T, c("row", "column", "value")] Ginv <- Ginv[order(Ginv$row, Ginv$column), ] return(Ginv) } No.subsets <- 10 IID <- fread(“iid.ref.txt”, header = F)$V1 A22 <- diag(1, length(IID)) w <- 0.95 for (i in 1: No.subsets) { load(paste(“ref.snp.s”, i, “.bin”, sep = “”)) pi <- colSums(M) / (2 * nrow(M)) P <- matrix((2 * pi), nrow = nrow(M), ncol = ncol(M), byrow = T) Z <- as.matrix(M - P) sum2pq <- 2 * (sum(pi * (1- pi))) Zt <- t(Z) ZZt <- eigenMatMult(Z, Zt) G <- (ZZt / sum2pq) Ga <- (w * G) + ((1 - w) * A22) invG <- solve(Ga) Ginv <- ginv(invG) attr(Ginv, “rowNames”) <- IID attr(Ginv, “INVERSE”) <- TRUE save(Ginv, file = paste(“Ginv.ref.s”, i, “.bin”, sep = “”)) } ##Before execution, save the following script to the file "matrix_multiplication.cpp" and place it in the working path. / / [[Rcpp::depends(RcppArmadillo, RcppEigen)]] #include<RcppArmadillo.h> #include<RcppEigen.h> / / [[Rcpp::export]] SEXP armaMatMult(arma::mat A, arma::mat B){ arma::mat C = A * B; return Rcpp::wrap(C); } / / [[Rcpp::export]] SEXP eigenMatMult(Eigen::MatrixXd A, Eigen::MatrixXd B){ Eigen::MatrixXd C = A * B; return Rcpp::wrap(C); } / / [[Rcpp::export]] SEXP eigenMapMatMult(const Eigen::Map <eigen::matrixxd>A, Eigen::Map <eigen::matrixxd>B) { Eigen::MatrixXd C = A * B; return Rcpp::wrap(C); } ##Run the following script in R to generate the file "coor4cv.ref.csv" required for 5-fold cross-validation. library(data.table) fold <- 5 No.cv <- 5 No.individual <- 500 res <- matrix(nrow = (No.indiv / fold), ncol = (fold * No.cv)) for (i in 1: fold) { coor <- data.table(seq = c(1: No.indiv), or = rnorm(No.indiv)) setorder(coor, -or) res[, (i * 5 - 4): (i * 5)] <- coor$seq } for (i in 1: ncol(res)) res[, i] <- sort(res[, i]) res[which(res == 499)] <- NA res[which(res == 500)] <- NA write.table("coor4cv.ref.csv", sep = ",", col.names = F, row.names = F, quote = F) ##Run the following script in R to evaluate the prediction accuracy of G-BLUP. In it, the name of the disease resistance phenotype variable is sur.status, and Ori.Location and b.weight are fixed items, which are the collection location and weight, respectively. library(data.table) library(pROC) library(asreml)s No.subset <- 10 ref.pheno <- asreml.read.table("ref.pheno4cv.csv”, sep = ",”, header = T) cv.coor <- as.matrix(fread("coor4cv.ref.csv”)) cv.res <- matrix(NA, nrow = (ncol(cv.coor) + 1), ncol = No.subset) colnames(cv.res) <- paste(‘s’, 1: ncol(cv.coor), sep = "”) rownames(cv.res) <- c(paste("CV”, 1: ncol(cv.coor), sep = "”), "Average”) for (i in 1: No.subset) { load(paste("Ginv.ref.s”, i, ".bin”, sep = "”)) for (j in 1: ncol(cv.coor)) { pheno4cv <- ref.pheno[-na.omit(cv.coor[, j]), ] gblup <- asreml(sur.status ~ Ori.Location + b.weight, random = ~vm(IID, Ginv), dense = ~vm(IID, Ginv), family = asr_binomial("probit”), data = pheno4cv, maxiter = 100) y4val <- ref.pheno$sur.status[na.omit(cv.coor[, j])] gebv4val <- coef(gblup)$random[na.omit(cv.coor[, j])] cv.res[j, i] <- roc(y4val, gebv4val)$auc } } cv.res [nrow(cv.res), ] <- colMeans(cv.res, na.rm = T) write.table(cv.res, "Ref.cv.results.gblup.csv", sep = ",", row.names = F, quote = F) After execution, the information in the last line of the file “Ref.cv.results.gblup.csv” is the accuracy with which G-BLUP predicts disease-resistant GEBV according to the change in the number of SNPs.

[0029] 3. Evaluate the prediction accuracy of BayesCπ Using the same model as G-BLUP, we evaluate the predictive accuracy of BayesCπ in R. Since the Bayes method is time-consuming, we specify one command for each SNP subset. Below, we use the SNP subset "s1" as an example to explain how to estimate the predictive accuracy of BayesCπ. In a Linux environment, run the following command: #Run the following script in R to evaluate the predictive accuracy of BayesCπ when using SNP subset "s1". library(data.table) library(BGLR) library(pROC) set.seed(123) i <- 1 #Change the value of i to use a different SNP subset load(paste("ref.snp.s", i, ".bin", sep = "")) ref.pheno <- fread ("ref.pheno4cv.csv") cv.coor <- as.matrix(fread("coor4cv.ref.csv")) No.CV <- ncol(cv.coor) gebv <- matrix(nrow = nrow(ref.pheno), ncol = No.CV) cv.res <- matrix(nrow = No.CV + 1) colnames(cv.res) <- paste("s”, i, sep = "”) rownames(cv.res) <- c(paste("CV”, 1: ncol(dr.coor), sep = "”), "Average”) for (i in 1: No.CV) { pheno4trn <- ref.pheno[-na.omit(cv.coor[, i]), ] geno4trn <- M[-na.omit(cv.coor[, i]), ] ETA4trn <- list(FIX = list(~factor(Ori.Location) + b.weight, data = pheno4trn, model = 'FIXED'), MAR = list(X = geno4trn, model = 'BayesC')) bc <- BGLR(y = pheno4trn$sur.status, ETA = ETA4trn, response_type = "ordinal”, nIter = 20000, burnIn = 5000, thin = 10, saveAt = "bc_”, verbose = F) MarEff <- bc$ETA$MAR$b gebv[, i] <- M %*% MarEff y4val <- ref.pheno$sur.status[na.omit(cv.coor[, i])] gebv4val <- gebv[na.omit(cv.coor[, i]), i] cv.res[i, ] <- roc(y4val, gebv4val)$auc } cv.res[nrow(cv.res), ] <- colMeans(cv.res, na.rm = T) write.table(cv.res, paste("Ref.results.cv.BC.s", i, ".csv", sep = ""), sep = ",", row.names = F, quote = F) After execution, the last line information of the file "Ref.results.cv.BC.s1.csv" is the accuracy of disease resistance predicted by BayesCπ when using SNP subset "s1". The above script is executed 10 times (i ranges from 1 to 10) to obtain the accuracy of BayesCπ predicting the resistance of Eastern star spot to Vibrio harveyi GEBV according to the change in the number of SNPs (Figure 2).

[0030] 4. Compare the prediction accuracy of the two genome selection methods under different numbers of SNPs The change in the accuracy of predicting anti-Vibrio harveyi GEBVs in East Star Spot by G-BLUP and BayesCπ according to the change in the number of SNPs is shown in Figure 2. As can be seen, the prediction accuracy by G-BLUP and BayesCπ tends to decrease with the decrease in the number of SNPs. When the number of SNPs is 50k, the prediction accuracy by G-BLUP and BayesCπ is the highest, 0.670 and 0.681, respectively. When the number of SNP markers decreases from 20k to 3k, the prediction accuracy by G-BLUP and BayesCπ decreases within an acceptable range, and when the number of markers is less than 2k, the prediction accuracy by the two methods decreases rapidly to 0.606 and 0.611. In short, the prediction accuracy of BayesCπ is slightly better than G-BLUP, and BayesCπ can be used as a method to estimate anti-disease GEBVs. The more markers used to estimate GEBV, the higher the cost of genotyping. Therefore, by using 20k SNP sites to estimate disease-resistant GEBV, the present invention not only ensures the prediction accuracy of the genome selection method, but also appropriately reduces the cost of genotyping and saves calculation time. EXAMPLES

[0031] Example 4: Selecting SNP markers for predicting Eastern Star Candidate Group GEBV Based on the results of Example 3, the present invention selects 20,000 (20k) SNP sites to predict the disease-resistant GEBV of the candidate group, including 30 sites related to Vibrio harveyi resistance in the Tonoseropsis variegata and 19,970 sites evenly distributed on the Tonoseropsis genome. The sites related to disease resistance are located at positions 36 and 107 of any of the sequences of SED ID NO:1 to SED ID NO:15, and are obtained by genome-wide association analysis. The other SNP sites are evenly distributed on the Tonoseropsis genome, and the sites are located at positions 36 and 107 of any of the sequences of SED ID NO:16 to SED ID NO:10000.

[0032] In a Linux environment, run the following command: #Prepare genotype data using PLINK plink2 --bpfile snps.imputed --autosome-num 24 --keep indiv.ref.txt --make-bpgen --out geno4gwas #Generate the GRM matrix using GCTA gcta --bpfile geno4gwas --autosome-num 24 --make-grm --out grm.gcta #Generate GRM sparse matrix using GCTA gcta --grm grm.gcta --make-bK-sparse 0.05 --out sparse.grm.gcta # Organize the disease resistance phenotype and covariates into text files "pheno4gcta.phe" and "qcov4gcta.qcov", respectively. The covariates include body weight and the top 10 principal component vectors from the PCA analysis. #Perform genome-wide association studies using GCTA gcta --bpfile geno4gwas --grm-sparse sparse.grm.gcta --pheno pheno4gcta.phe-qcovar qcov4gcta.qcov --fastGWA-mlm --autosome-num 24 --out res.gwas #Use R to read the analysis results and output the 30 SNP information related to disease resistance properties library(data.table) res <- fread("res.gwas.fastGWA") setorder(res, P) sel.loci <- res[1: 30, ] setorder(sel.loci, CHR, POS) write.table(sel.loci, "selected.loci.info.txt", sep = "\t", quote = F) After execution, information such as the site numbers, chromosome numbers, and physical locations of the 30 SNPs related to disease resistance are saved in the file "selected.loci.info.txt". Then, 20k SNPs used for calculating GEBV are selected, and these 20k SNPs include the above 30 SNPs.

[0033] In a Linux environment, run the following command: #Select 20k SNPs that evenly cover the genome ##Using PLINK2, 19,970 SNPs that cover the genome evenly were extracted. Increasing the value of --bp-space resulted in fewer SNPs being retrieved, whereas decreasing the value of --bp-space resulted in more SNPs being retrieved. plink2 --bpfile snps.imputed --extract comm.loci.txt --bp-space 37267 --make-bpgen --out snps.20k.part1 ##Using R to clean up the files to extract the last 20k SNPs library(data.table) bim.p1 <- fread("snps.20k.part1.bim") loci.info <- fread("trait.realted.loci.txt") loci.20k <- c(bim.p1$V2, loci.info$SNP) write.table(loci.20k, "snps.20k.txt", sep = "\t", col.names = F, row.names = F, quote = F) After execution, the file "snps.20k.txt" contains 20k SNP numbers, which can be read by PLINK2 to extract the corresponding sites. The file contains 20,000 SNP sites in total, including 30 SNP sites related to the disease resistance of Toushiroshima and 19,970 SNP sites that uniformly cover the Toushiroshima reference genome. Figure 3 shows the distribution of these 20k SNPs on the Toushiroshima reference genome in every 1M range, and it can be seen that these sites basically cover the Toushiroshima genome evenly. EXAMPLES

[0034] Example 5: Verifying the predictive effect of Bayes Cπ and preferred 20k SNP markers 1. Evaluate prediction effectiveness using cross-validation The same model as described in Example 3 is used to verify the 20k SNPs obtained in Example 4. Compared with the results of Example 3, it can be seen that the prediction accuracy of BayesCπ using the evenly distributed 20k SNP sites of the Eastern Star Spotted genome (not including the 30 SNP sites related to disease resistance selected in Example 4) is 0.670 (Figure 2), and after using 20k SNP sites including these 30 disease resistance-related sites, the prediction accuracy of BayesCπ is 0.682, which is increased by 0.012 and is close to the prediction accuracy (0.681) of BayesCπ using 50k SNPs in Example 3. Therefore, the use of the 20k SNPs provided by the present invention not only can successfully complete genome selection, but also has good prediction effect.

[0035] 2. Evaluate the predictive effect using disease-resistant GEBV in the validation group The individuals with disease-resistant phenotype but not included in the reference group are used as validation individuals, and the disease-resistant GEBV of the validation group is calculated using BayesCπ and the 20k SNPs obtained in Example 4, and the average GEBV of the surviving individuals and the dead individuals is compared to validate the feasibility of the above method of the present invention. Run the following command in a Linux environment: #Extract the genotypes of candidate individuals ##Execute the following command on PLINK2 plink2 --bpfile snps.imputed --autosome-num 24 --extract snps.20k.txt --keep indiv.ref.txt --export-allele export.A.txt --export A --out geno.Ref.20k plink2 --bpfile snps.imputed --autosome-num 24 --extract snps.20k.txt --keep indiv.val.txt --export-allele export.A.txt --export A --out geno.Val.20k awk '{print $2}' geno.Val.20k.raw | sed '1d' > iid.val.20k.txt #The GEBV of the New York City library ( data . table ) library(BGLR) set.seed(123) pheno.ref <- fread("pheno.Ref.csv”) pheno.val <- fread("pheno.Val.csv”) iid . val < - fread ( " iid . val . 20k . txt ”, header = F ) $ V1 M.ref <- fread("gen.Ref.20k.raw”) M.val <- fread("geno.Val.20k.raw”) M.ref <- as.matrix(M.ref[, c(1:6):= NULL]) M.value <- as.matrix(M.value[,c(1:6):=NULL]) ETA4MarEst <- list(FIX = list(~factor(Ori.Location)+b.weight, data=pheno.ref, model='FIXED'), MAR=list(X=M.ref, model='BayesC')) BC <- BGLR(y = pheno$south.status, ETA = ETA4SeaEast, response_type = "ordinal", nIter = 20000, burnIn = 5000, thin = 10,000. saveAt = "BC4Val_", verbose = T) MarEff <- BC$ETA$MAR$b gebv <- data.frame(IID = iid.val, GEBV = M.val %*% MarEff) gebv <- merge(gebv, pheno.val, sort = F) mean(gebv$GEBV[gebv$sur.status == 0]) mean(gebv$GEBV[gebv$sur.status == 1]) length(gebv$GEBV[gebv$sur.status == 0]) length(gebv$GEBV[gebv$sur.status == 1]) After implementation, the disease resistance GEBV of the 298 test groups was successfully predicted, and the average GEBV of the surviving and dead individuals is shown in Figure 4. Among them, there were 110 surviving individuals, with an average GEBV of 0.413, which was higher than the overall average GEBV, and there were 188 dead individuals, with an average GEBV of 0.255, which was lower than the overall average GEBV. This shows that GEBV is positively correlated with the disease resistance of individuals, and individuals with high GEBV have strong disease resistance, so GEBV can be used to select individuals, and then those with strong disease resistance can be selected for breeding. EXAMPLES

[0036] Example 6: Calculating the disease-resistant GEBV of the Eastern Star Spot Candidate Group 100 juveniles of East Star Spot are selected as a candidate group, and whole genome resequencing is used to obtain the genotype of the candidate group, and the genotype information is saved in the file "candidates.geno.vcf.gz". After the sequencing is completed, the same 20k SNPs as described in Example 4 are extracted and combined with the reference group in Example 3, and the disease resistance GEBV of the candidate group is calculated using BayesCπ. To further illustrate the strong disease resistance of individuals with high disease resistance GEBV, the present invention uses infection experiments to measure the disease resistance of the candidate group, and surviving individuals are identified as individuals with high disease resistance, and dead individuals are identified as individuals susceptible to infection. These 100 East Star Spot candidate groups each contain 50 individuals with high disease resistance and 50 individuals susceptible to infection. Run the following command in a Linux environment: # Extract the 20k SNPs described in Example 4 from the candidate group genotype file. ##Run the following script in R library(data.table) snps.20k <- fread("snps.20k.txt") pos <- as.integer(matrix(unlist(strsplit(snps.20k$V1, ':')), ncol = 2, byrow = T)[, 2]) chr <- gsub("rs", "chr", matrix(unlist(strsplit(snps.20k$V1, ':')), ncol = 2, byrow = T)[, 1] info <- data.table(CHR = chr, POS = pos) fwrite(info, "info.20k.txt", sep = "\t", col.names = F, quote = F) ##Run the following command with VCFtools vcftools --gzvcf candidates.geno.vcf.gz --positions info.20k.txt --recode -- recode-INFO-all --stdout | bgzip -c -f -@ 10 > candidates.geno.20k.vcf.gz ##Execute the following command on PLINK2 plink2 --vcf candidates.geno.20k.vcf.gz --autosome-num 24 --make-bpgen --out Can.geno.20k ##Run the following script in R library(data.table) bim <- fread("Can.geno.20k.bim") bim$V2 <- paste("rs", bim$V1, ':', bim$V4, sep = "") fwrite(bim, "Can.geno.20k.bim", sep = "\t", col.names = F, quote = F) ##Execute the following command on PLINK2 plink2 --bpfile Can.geno.20k --autosome-num 24 --export-allele export.A.txt --export A --out geno.Can.20k awk '{print $2}' geno.Val.20k.raw | sed '1d' > iid.can.20k.txt #Calculate the candidate group GEBV. The first 50 individuals are susceptible individuals, and the last 50 individuals are highly disease resistant individuals. ##Run the following script in R library(data.table) library(BGLR) setseed(123) pheno.ref <- fread("pheno.Ref.csv”) iid.can <- fread("iid.can.20k.txt”, header = F)$V1 M.ref <- fread("geno.Ref.20k.raw”) M.can <- fread("geno.Can.20k.raw”) M.ref <- as.matrix(M.ref[, c(1: 6):= NULL]) M.can <- as.matrix(M.can[, c(1: 6):= NULL]) ETA4MarEst <- list(FIX = list(~factor(Ori.Location) + b.weight, data = pheno.ref, model = ‘FIXED’), MAR = list(X = M.ref, model = ‘BayesC’)) BC <- BGLR(y = pheno$sur.status, ETA = ETA4MarEst, response_type = "ordinal”, nIter = 20000, burnIn = 5000, thin = 10, saveAt = "BC_”, verbose = F) MarEff <- BC$ETA$MAR$b gebv <- data.frame(IID = iid.can, GEBV = M.val %*% MarEff, pheno = c(rep(0, 50), rep(1, 50))) fwrite(gebv, "GEBV.candidates.csv”, sep = ",”, quote = F) After implementation, the disease resistance GEBV of 100 East Star Spot candidate groups was successfully calculated, and the average GEBV of 50 highly disease-resistant individuals was 0.643, which was higher than the average GEBV of susceptible individuals (0.525) (Figure 5). In addition, among the top 20% of GEBV individuals, there were 15 highly disease-resistant individuals, which was much more than the susceptible individuals (5). This further indicates that individuals with high GEBV have strong disease resistance and can be used as parents to breed offspring to improve the disease resistance of offspring. At the same time, it also shows that the Bayes Cπ method selected in Example 3 of the present invention and the 20k SNPs preferably selected in Example 4 have good predictive effect and can be used to screen individuals with strong disease resistance.< / eigen::matrixxd> < / eigen::matrixxd>

Claims

1. A method for breeding a disease-resistant variety of Tosei based on genomic selection, comprising: 1) To generate a reference group for disease resistance of Tosei with Tosei juveniles from different geographical origins; and generating a validation group with individuals having a disease-resistant phenotype but not included in the reference group to validate the selection effect of the method; 2) performing whole genome resequencing on the reference and validation groups to capture high-quality single nucleotide polymorphism (SNP) sites within the genome; 3) comparing two methods, the best linear unbiased prediction G-BLUP and Bayes Cπ, for the tendency of the accuracy of predicting disease resistance trait genomic estimated breeding value GEBV to change depending on the number of SNPs, thereby selecting the number of SNPs and the method for calculating the GEBV of the candidate group; 4) selecting a SNP marker containing a SNP site associated with the disease resistance properties of the Eastern Star Spot for calculating the GEBV of the Eastern Star Spot candidate group based on the results obtained in step 3); 5) calculating the GEBV of the candidate group using the genome selection method selected in step 3) and the SNP markers selected in step 4), and recognizing individuals with high GEBV as having strong disease resistance and using their offspring as parents for breeding. A method for breeding a disease-resistant variety of Tosei blotch based on genomic selection, comprising:

2. In step 1), the juveniles of Tosei for generating the disease resistance reference group are obtained by screening for disease resistance; 2. The method of claim 1 .

3. In the method, a validation group is used to validate the prediction effect of the genome selection method selected in step 3) and the SNP marker selected in step 4); 2. The method of claim 1 .

4. The SNP markers are located at positions 36 and 107 of any of the sequences of SED ID NO: 1 to SED ID NO: 10000; 2. The method of claim 1 .

5. The SNP sites associated with the disease-resistant properties of eastern star spot are located at positions 36 and 107 of any of the sequences of SED ID NO: 1 to SED ID NO: 15; The genome selection method according to claim 1 .

Images