Method for screening and breeding economic character related molecular markers of litopenaeus vannamei by using apparent annotation information

Abstract

The invention discloses a method for screening molecular markers related to economic traits of litopenaeus vannamei and breeding by using apparent annotation information, which comprises the following steps of: A, performing phenotype determination and genome sequencing on growth traits of a litopenaeus vannamei breeding population to obtain phenotype data and SNP (Single Nucleotide Polymorphism) typing data; step B, performing GWAS analysis by using the phenotype data and the SNP typing data to obtain a P value of each SNP; step C, using the P value of each SNP to carry out Bonferroni multiple correction to obtain a P value lt; 1.9 E-08 is used as a significant site for judgment; sNPs significantly related to growth traits are screened out; the key histone modification region is utilized to accurately locate the SNPs significantly related to the growth traits and the expression regulation element of the corresponding gene thereof. By applying the apparent annotation information to GWAS analysis, a more accurate molecular marker is provided for genetic analysis of growth traits of litopenaeus vannamei, a foundation is laid for genetic breeding research of other economic traits, and the breeding process is accelerated.

Description

Technical Field

[0001] The invention belongs to the technical field of molecular breeding of Penaeus vannamei, and specifically relates to a method for screening molecular markers related to economic traits of Penaeus vannamei and breeding using epigenetic annotation information. Background Art

[0002] Selective breeding is an important means of genetic improvement of farmed animals. In recent years, with the development of molecular biotechnology and the continuous reduction of high-throughput sequencing costs, single nucleotide polymorphism (SNP) molecular markers widely present in the genome have become a trend in breeding of farmed animals. Whole-genome selective breeding uses molecular markers, mainly single nucleotide polymorphism markers, in the whole genome to construct statistical models for genetic breeding value estimation. At present, high-throughput typing technologies mainly include whole-genome resequencing and SNP chip technology. Whole-genome resequencing can obtain SNP site typing information in the whole genome when the sequencing depth is sufficient, but the cost is high. In actual breeding applications, in order to reduce the cost of genotyping and model calculation costs, SNP genotyping chips are usually developed by screening SNP sites related to target traits for genotyping of the test population. Although the universal high-density SNP chip has a wide genome coverage and a wide range of applications, the cost of use is still too high for production applications due to the large number of SNPs used, which seriously hinders the practical application of whole-genome selection methods. In addition, due to the changes in SNP genotype frequencies in different populations, the whole genome selection model constructed through the training population can often only obtain good breeding value estimation results in the internal population, and the breeding value estimation in cross-populations is often difficult to achieve satisfactory accuracy.

[0003] Litopenaeus vannamei is one of the most widely farmed aquatic species in the world and has important economic value. Growth, as a ubiquitous biological process, is one of the most concerned traits in the breeding of aquatic animals, including most other economic animals. At present, the research on genetic analysis of growth traits and genetic breeding mainly uses genome-wide association studies (GWAS) to screen loci and genes related to growth traits, and genome-wide selection analysis (GS) to select breeding by using high-density SNPs covering the whole genome.

[0004] GWAS analysis uses the genotype data of a large number of individuals and the corresponding phenotypic data to find molecular marker sites (such as SNPs) and genotypes that are significantly associated with a specific phenotype, and finds the strong linkage region of significant SNPs based on linkage disequilibrium (LD). The genes located in the region are used as candidate genes related to the trait. However, due to the existence of LD, when a locus with functional association with a trait is closely linked to another locus that is not related to the trait, it is often difficult to define which locus is the real cause of the trait difference, which may lead to false positive and false negative results in GWAS gene positioning studies. GS analysis uses high-density SNPs covering the entire genome for selective breeding. It is based on the LD state of each SNP and its adjacent quantitative trait locus (Quantitative Trait Locus, QTL) to trace the effects of all potential QTLs, and to achieve the prediction of individuals with unknown phenotypes while avoiding QTL positioning. Therefore, a sufficiently high SNP density and genetic communication between the test population and the training population are required to ensure that the prediction performance of the model breeding value is optimal. This not only leads to higher genotyping costs in subsequent applications, but also due to the randomness of the screened SNP sites (no functional association), some sites with no actual functionality are retained in the breeding model, while other functional sites are eliminated.

[0005] Therefore, it is necessary to establish a method for screening and breeding relevant molecular markers using epigenetic annotation information, to provide more accurate molecular markers for the genetic analysis of growth traits of Penaeus vannamei, significantly improve the accuracy of genetic breeding value estimation, and the effect of cross-population breeding value estimation, and lay the foundation for genetic breeding research on other economic traits. Summary of the invention

[0006] The purpose of the present invention is to apply histone modification information and epigenetic annotation information to molecular marker screening and breeding, to provide more accurate molecular markers for the genetic analysis of growth traits of cultured Penaeus vannamei, to significantly improve the accuracy of genetic breeding value estimation, as well as the effect of cross-population breeding value estimation, and to lay a foundation for genetic breeding research of cultured animals and other economic traits.

[0007] To solve the above technical problems: The present invention provides a method for screening molecular markers related to economic traits of Penaeus vannamei, comprising the following steps: step A: performing phenotypic measurement of growth traits and genome sequencing on a breeding population of Penaeus vannamei to obtain phenotypic data and SNP typing data; step B: performing GWAS analysis using the phenotypic data and the SNP typing data to obtain a P value for each SNP; step C: using the P value of each SNP to perform a Bonferroni multiple correction method, taking a P value <1.9E-08 as a significant site judgment, and screening out SNPs significantly related to growth traits; and using a key histone modification region to accurately locate SNPs significantly related to growth traits and expression regulatory elements of their corresponding genes.

[0008] Furthermore, step A comprises: step A-1: ​​performing phenotypic measurement of growth traits on the breeding population of Penaeus vannamei to obtain the phenotypic data; the growth trait is body length or weight; step A-2: extracting total DNA from the breeding population of Penaeus vannamei to obtain total DNA; constructing a DNA library using the total DNA; sequencing the DNA library to obtain the SNP typing data; wherein the construction process of the DNA library comprises: breaking the total DNA into small fragments by ultrasonic treatment, and performing end modification, linker connection, PCR amplification and magnetic bead purification.

[0009] Furthermore, step C includes: step C-1: sorting the P value of each SNP in ascending order, and judging the significant site with P value <1.9E-08 through Bonferroni multiple correction method; screening out SNPs significantly related to growth traits; step C-2: using CUT&Tag to screen the key histone modification regions, and using the key histone modification regions to accurately locate the SNPs significantly related to growth traits and the expression regulatory elements of their corresponding genes; the key histone modification regions are H3K4me1, H3K4me3, H3K27me3 and H3K27ac.

[0010] On the other hand, the present invention provides a method for breeding Penaeus vannamei based on whole genome data, comprising the following steps: step A: performing phenotypic measurement of growth traits and genome sequencing on a breeding population of Penaeus vannamei to obtain phenotypic data and SNP typing data; step B: performing GWAS analysis using the phenotypic data and the SNP typing data to obtain a P value of each SNP in the GWAS analysis; step C: using the P value of each SNP to perform a Bonferroni multiple correction method, taking a P value <1.9E-08 as a significant site judgment, and screening out SNPs significantly associated with growth traits; using a key histone modification region to accurately locate SNPs significantly associated with growth traits and expression regulatory elements of their corresponding genes; step D: determining a preferred genotype in the breeding population of Penaeus vannamei according to the SNP site subset, and breeding by retaining parents with the preferred genotype for pairing and breeding.

[0011] On the other hand, the present invention provides a computer device for a method for breeding Penaeus vannamei based on whole genome data, comprising a memory, a processor, and a computer program stored in the memory and executable on the processor, wherein when the processor executes the computer program, the method for breeding Penaeus vannamei based on whole genome data is implemented.

[0012] On the other hand, the present invention provides an application of a SNP molecular marker in assisted breeding of Penaeus vannamei, wherein the sequence of the SNP molecular marker is shown in SEQ ID NO.1; the SNP site is located at the 233rd position from the 5' end of the sequence shown in SEQ ID NO.1; the polymorphism of the SNP molecular marker is A / G type; and the application in assisted breeding of Penaeus vannamei is for screening the body length or weight of Penaeus vannamei.

[0013] The present invention innovatively uses different histone modifications for epigenetic annotation, obtains cis-regulatory element information and applies it to GWAS analysis, which can obtain expression regulatory element information of SNPs significantly associated with growth traits, which is conducive to a deeper study of its molecular mechanism. It provides more accurate molecular markers for the genetic analysis of growth traits of Penaeus vannamei, lays a foundation for genetic breeding research of other economic traits, and accelerates the breeding process. BRIEF DESCRIPTION OF THE DRAWINGS

[0014] The above content and the following specific embodiments of the present invention will be better understood when read in conjunction with the accompanying drawings. It should be noted that the accompanying drawings are only examples of the technical solutions claimed for protection.

[0015] Figure 1Manhattan plot of genome-wide association analysis of body length traits of Penaeus vannamei (the black line indicates the significance threshold, and after Bonferroni correction for multiple tests, sites with P < 1.9E-08 are considered to be significantly associated); Figure 2 Manhattan plot of genome-wide association analysis of body length traits of Penaeus vannamei (the black line indicates the significance threshold, and after Bonferroni correction for multiple tests, sites with P < 1.9E-08 are considered to be significantly associated); Figure 3 This is the Venn diagram of significant loci in the genome-wide association analysis of body length / weight traits of Penaeus vannamei (the blue area length_SNPs is the number of SNP loci for body length traits; the pink area weigth_SNPs is the number of SNP loci for body weight traits; the intersection of the blue area and the pink area is the number of SNP loci for body length / weight traits); Figure 4 The following are the growth trait-related significant loci and key genes that are modified by H3K4me3 in each key period of embryonic development of Litopenaeus vannamei (the scale of pfdn6-like gene is 2 kb; the black frame marks the SNP site 44_1839577; the brackets mark the H3K4me3 histone modification at the promoter; blas is the blastocyst stage; gast is the gastrula stage; lbe is the limb bud larvae stage; lim is the intramembranous larvae stage; N1 is the nauplii stage I; N3 is the nauplii stage III; N6 is the nauplii stage VI; Z1 is the zoa stage I; Z3 is the zoa stage III; M1 is the mysid larvae stage I; M3 is the mysid larvae stage III; the control is IgG); Figure 5 This is the conservative analysis of the significant loci of Penaeus vannamei in this population and other populations (where A is the body length trait and genotyping statistics of the SNP locus (44_1839577) in this population; B is the body weight trait and genotyping statistics of the SNP locus (44_1839577) in this population; C is the body length trait and genotyping statistics of the SNP locus (44_1839577) in other populations; D is the SNP locus (44_1839577) in other populations. Statistical chart of body weight trait and genotyping; body length unit is mm; body weight unit is g; group1 is the first batch population; group2 is the second batch population; group3 is the third batch population; "*" represents p value < 0.05, there is a difference; "**" represents p value < 0.01, there is a significant difference; "****" represents p value < 0.001, there is an extremely significant difference; AA genotype is 0 / 0; AG genotype is 0 / 1; TAB represents other groups); Figure 6Box plot of the prediction accuracy of the breeding value of body length phenotype for different marker densities and different models of Penaeus vannamei (where the horizontal axis is the SNP marker density and the vertical axis is the prediction accuracy; purple is the GBLUP model; green is the BayesA model; orange is the BayesB model; blue is the BayesC model); Figure 7 This is a chromatin state analysis diagram for the critical period of embryonic development of Penaeus vannamei (where A is the IGV map of chromatin state trajectories during embryonic development. ChromHMM learns the chromatin state definition based on these trajectories and assigns each genomic position to the corresponding state; B is the emission parameter heat map of 8 chromatin state extension models based on H3K4me1, H3K4me3, H3K27me3 and H3K27ac. The red depth of the heat map reflects the probability of histone marks appearing, and darker red indicates a higher probability of marks; C is an overlapping enrichment heat map of chromatin states and genomic annotated regions, showing the enrichment of each chromatin state in specific genomic features; D is a spatial enrichment heat map of chromatin states in the vicinity of transcription start sites (TSS), showing the enrichment pattern of chromatin states around transcription start sites; E is the annotation of each state and its abbreviation map; the annotations are: Flanking bivalent TSS / Enh 1, BivFInk1), Quiescent / low (Quies), Weak repressed Polycomb (ReprPCWk), Repressed Polycomb (ReprPC), Weak Enhancer (EnhWk), Flanking bivalent TSS / Enh 2 (BivFInk2), Weak transcription (TxWk), and Flanking TSS downstream (TssFlnkD); Figure 8This is a box plot of the prediction accuracy of the SNP marker density in different chromosomal states during the critical period of embryonic development of Penaeus vannamei for the breeding value of the phenotypic trait of body length (wherein the horizontal axis is the SNP marker density and the vertical axis is the prediction accuracy; embryo_E1 is: SNPs located in the BivFlnk1 state; embryo_E2 is: SNPs located in the Quies state; embryo_E3 is: SNPs located in the ReprPCWk state; embryo_E4 is: SNPs located in the ReprPC state; embryo_E5 is: SNPs located in the EnhWk state; embryo_E6 is: SNPs located in the BivFlnk2 state; embryo_E7 is: SNPs located in the TxWk state; embryo_E8 is: SNPs located in the TssFlnkD state; random is a randomly selected SNP site). DETAILED DESCRIPTION

[0016] The detailed features and advantages of the present invention are described in detail in the specific implementation modes below, and the contents are sufficient to enable any person skilled in the art to understand the technical contents of the present invention and implement them accordingly. Moreover, according to the description, claims and drawings disclosed in this specification, those skilled in the art can easily understand the relevant objects and advantages of the present invention.

[0017] It should be noted that in this specification, similar reference numerals and letters denote similar items in the following drawings, and therefore, once an item is defined in one drawing, it does not need to be further defined and explained in the subsequent drawings.

[0018] In order to make the purpose, technical scheme and advantages of the present invention clearer, the embodiments of the present invention will be further described in detail below with reference to the accompanying drawings. The experimental methods described in the embodiments of the present invention are conventional methods unless otherwise specified, and the materials, reagents, etc. used in the following embodiments can be obtained from commercial channels unless otherwise specified.

[0019] Example 1 A method for screening molecular markers related to economic traits of Litopenaeus vannamei by using key histone modification region information and GWAS joint analysis, comprising the following steps: S1. The phenotypic data collection of individual Litopenaeus vannamei and the extraction of total DNA from muscle tissue are as follows: (1) A total of 972 three-month-old Penaeus vannamei (purchased from Hebei Huanghua Hatchery) were randomly selected and divided into three batches (the first batch, the second batch, and the third batch). Among them, there were 313 Penaeus vannamei in the first batch, 331 in the second batch, and 328 in the third batch.

[0020] (2) The three batches of Penaeus vannamei sampled above were individually numbered, and their body length and weight were measured using an electronic vernier caliper with an accuracy of 0.02 mm and an electronic balance with an accuracy of 0.01 g (Quintix). The number of each individual Penaeus vannamei was repeated three times and the average value was taken.

[0021] Among them, the body length of the first batch of Penaeus vannamei individuals was between 160~280 mm, and the weight was between 36~108 g; the body length of the second batch of Penaeus vannamei individuals was between 160~240 mm, and the weight was between 24~84 g; the body length of the third batch of Penaeus vannamei individuals was between 120~200 mm, and the weight was between 12~48 g.

[0022] (3) DNA was extracted from the muscle of each individual Penaeus vannamei using a marine animal genomic DNA extraction kit (Tian Gen, DP324-03, the following reagents are included in the kit). The DNA extraction process was as follows: about 30 mg of tissue sample was added with GA buffer and Proteinase K (Protease K, used to decompose RNA and protein in tissue samples and protect DNA from damage), then treated and ground in a 56°C water bath until uniform, and RNase A was added to degrade RNA molecules in the sample. After adjusting the lysis time to optimize the DNA extraction effect, GB buffer was added for further treatment, and the marine animal genomic DNA extraction kit was used to extract DNA from the muscle tissue of Penaeus vannamei. The specific process included tissue grinding, DNA precipitation, column elution, and DNA dissolution and elution. The total DNA of the muscle tissue of Penaeus vannamei was obtained and stored at -20°C for future use. The concentration and purity of the DNA were tested using a nucleic acid analyzer. The OD value of the total DNA of the muscle tissue of 972 Penaeus vannamei individuals was 0.0347 W / m. 260 / OD 280 All of them were between 1.8 and 2.0, indicating that total DNA with good purity was obtained. The extraction of total DNA was complete as detected by 1.0% agarose gel electrophoresis and could be used for subsequent experiments.

[0023] S2. Construction of DNA library and sequencing The total DNA was broken into small fragments by ultrasonic treatment, and a DNA library was constructed by end modification, adapter ligation, PCR amplification and magnetic bead purification; the DNA library was used to perform PE150 double-end sequencing on the Illumina NovaSeq 6000 platform, that is, the original data file was quality controlled using Trimmomatic v3.29 software; the reads were aligned to the genome using BWA v0.7.17 software to obtain a SAM file; the SAM format file was converted into a BAM format file using Samtools v1.17 software; site variation detection was performed using GATK v4.4.0.0 software; SNP site variation was extracted using GATK v4.4.0.0 software and VCFtools v0.1.16 software, and missing genotypes were supplemented using Beagle v5.5 software, and qualified SNP sites with MAF>0.05 and Harwin equilibrium test P value<0.0001 were screened to obtain the SNP typing data. The specific reference steps are as follows: (1) The qualified total DNA is randomly broken into small fragments by the Covaris ultrasonic processing system, and these small fragments are used for end modification. An A base is added to the 3' end of the small fragment by Klenow enzyme (T4 DNA polymerase & DNA polymerase I). After the addition of the A base, the original small fragment changes from a flat end to a sticky end, which is easier to connect to the subsequent primers and adapters; the end of the small fragment after end modification has a protruding A tail, and the adapter has a protruding T tail. T4 DNA ligase is used to connect the adapter to the two ends of the DNA fragment. After the adapter is successfully connected, the low-cycle technology is used to modify the adapter.

[0024] (2) The system after adding the adapter contains polymerase and ligase, etc., and there may be large fragments. Therefore, double screening is performed by magnetic beads to remove large fragments and impurities, so as to obtain library fragments with successful addition of adapters. The DNA fragments with adapters are amplified with primers complementary to the adapters. After PCR, magnetic bead purification is required again to separate the PCR products from impurities. The obtained PCR products are the constructed DNA library. The PCR products are quantified by Qubit DNA HS ASSAY KIT; 2100 High Sensitivity DNA Chip electrophoresis is performed to determine whether the fragment size meets the requirements of subsequent sequencing (the fragment size is generally about 400bp); the molar concentration is calculated based on the Qubit quantitative results and the fragment size detected by the 2100 chip.

[0025] (3) The PCR products obtained above were used to perform PE150 double-end sequencing on the Illumina NovaSeq platform with a sequencing depth of 2G / sample to generate raw sequencing data. The sequencing data was aligned with the unpublished reference genome of the Penaeus vannamei chromosome in our laboratory for SNP typing, and the haplotype library of our laboratory was used to fill in (fill in missing genotypes) using Beagle v5.5 software, filtering SNPs with a minimum allele frequency (MAF) greater than 0.05 and a Hawes-Wenzhou equilibrium test (HWE) P value less than 0.0001, and a total of 2,609,549 qualified SNP sites were obtained.

[0026] S3. GWAS analysis, the specific reference steps are as follows: (1) GWAS analysis was performed on 2,609,549 qualified SNP loci screened by sequencing from the above 972 different batches of Penaeus vannamei. The genotype and phenotype data of Penaeus vannamei were imported separately, and GWAS analysis was performed using the key parameters gemma -bfile DX -gk 2 -p phe.txt -o length, gemma -bfile DX -k length.sXX.txt-lmm 1 -p length.txt -c group.txt -o length_GWAS.

[0027] (2) According to the P value of each SNP obtained by GWAS analysis, in each SNP analysis data, the P value of each SNP was sorted in ascending order. After Bonferroni correction for multiple tests in traditional GWAS analysis, sites with a P value < 1.9E-08 were considered to be significantly associated. Among the 2,609,549 qualified SNP sites, a total of 883 significantly associated SNP sites were screened with body weight as the phenotype (such as Figure 1 As shown); Taking body length as the phenotype, a total of 493 significantly associated SNP loci were screened (such as Figure 2 The significant loci of body length phenotype are all located within the significant loci of body weight phenotype, that is, the two are overlapped (such as Figure 3 Therefore, 883 SNP loci significantly associated with body weight were subsequently analyzed.

[0028] S4. Histone modification information is applied to GWAS analysis. The specific reference steps are as follows: (1) Through CUT&Tag experiment and data processing (using fastp software to perform quality control on raw data files; using BWA v0.7.17 software to align reads to the genome to obtain SAM files; using awk to filter the fragment length of the generated SAM files, retaining the alignment results with insert fragment lengths between 10 and 1000 bp; using Samtools v1.17 software to convert SAM format files into BAM format files; using MACS2 v2.1.4 software to perform PeakCalling analysis on BAM format alignment files, using IgG samples as controls and setting the significance threshold to q < 0.05; using IGV v2.19.1 software to visualize specific genes), peak files of different histone modifications were obtained, including peak chromosomes, start sites, and end sites. Genotype bim files included chromosome and SNP site information. The data.table R package was used first to use match <- peak[site, .(chr, peak.start, peak.end,peak.id), on = .(chr = chr,peak.start <= pos, peak.end >= pos)] command was used to screen SNPs located in different histone modification peak intervals.

[0029] The annotation results of SNP sites with key histone modification regions are shown in Table 1. Among the 883 SNP sites, there are 43 significant sites located in H3K4me1 peaks, accounting for 4.87% of the total number of significant sites, 47 significant sites located in H3K4me3peaks, accounting for 4.87% of the total number of significant sites, 25 significant sites located in H3K27me3 peaks, accounting for 2.83% of the total number of significant sites, and 39 significant sites located in H3K27ac peaks, accounting for 4.42% of the total number of significant sites. There are 104 SNP sites located in histone modification, accounting for 11.78% of the total number of significant sites. These SNPs located in histone modification regions mean that they are located in the genomic regulatory element region, such as active or silent enhancer regions and promoter regions.

[0030] Table 1 Significant sites located in histone modification regions

[0031] (2) Further screening of SNPs and key genes significantly associated with growth traits. The results are as follows: Figure 4As shown in the figure, taking the SNP site 44_1839577 (SEQ ID NO.1, 233rd position from the 5' end, polymorphism is A / G type) with significant correlation with growth traits in the H3K4me3 histone modification region as an example, this site is located in the pfdn6-like promoter region modified by H3K4me3, and the promoter region (including this site) is modified by H3K4me3 from the blastocyst stage to the mysid larvae stage, proving that H3K4me3 is a sign of gene transcription activation, and PFDN6 (Prefoldin Subunit 6) has an important function in the cytoplasmic folding of actin and tubulin monomers during the cytoskeleton assembly process, which is closely related to growth. Therefore, the use of histone modification information can effectively identify SNPs located in the regulatory region of the genome and reduce the probability of false positives and false negatives.

[0032] (3) To further verify the conservation of SNP site 44_1839577 in this population and other populations (TAB), statistical analysis was performed on the growth trait (body length, weight) phenotypic data and SNP typing data of this SNP site (44_1839577). Since the number of individuals with the GG genotype screened this time was too small and the GG genotype was not the preferred genotype, it was not included in the significance comparison.

[0033] The statistical analysis results of this group are shown in Table 2, Table 3 and Figure 5 As shown in A and B, the AA genotype is 0 / 0, and the AG genotype is 0 / 1. It can be seen that in the three batches of Penaeus vannamei in this group, the number of Penaeus vannamei with the AA genotype is greater than the number of Penaeus vannamei with the AG genotype.

[0034] In the first batch: the body length of the AA genotype was 224.35±1.14 mm, and the body length of the AG genotype was 214.26±1.83 mm, that is, there was a very significant difference between the two; the body weight of the AA genotype was 69.13±0.84 g, and the body weight of the AG genotype was 61.54±1.33 g, and there was a very significant difference between the two. In the second batch: the body length of the AA genotype was 201.04±0.83 mm, and the body length of the AG genotype was 197.87 plus or minus 1.18 mm, and there was a significant difference between the two; the body weight of the AA genotype was 49.95±0.51 g, and the body weight of the AG genotype was 47.24±0.80 g, and there was a significant difference between the two. In the third batch: the body length of the AA genotype was 169.95±0.71 mm, and the body length of the AG genotype was 169.79±1.48 mm, and there was a difference between the two; the body weight of the AA genotype was 31.11±0.34 g, and the body weight of the AG genotype was 31.11±0.70 g, and there was a difference between the two.

[0035] Table 2 Statistics of SNP genotyping results for growth traits of this population

[0036] Table 3 Statistics of phenotypic data corresponding to SNP genotyping of growth traits in this population

[0037] The statistical analysis results of other groups (TAB) are shown in Table 4 and Figure 5 As shown in C and D, the AA genotype is 0 / 0 and the AG genotype is 0 / 1. It can be seen that among the other groups of Penaeus vannamei, the number of Penaeus vannamei with the AA genotype is 17, and the number of Penaeus vannamei with the AG genotype is 14. The body length of the AA genotype is 113.22±0.15 mm, and the body length of the AG genotype is 104.66±0.22 mm, and there is a significant difference between the two; the body weight of the AA genotype is 10.55±0.37 g, and the body weight of the AG genotype is 8.67±10.42 g, and there is a significant difference between the two.

[0038] Table 4 Phenotypic data statistics corresponding to SNP genotyping of growth traits of other populations

[0039] In this population and other populations (TAB), it was found that individuals with the AA genotype at the SNP site 44_1839577 were older than those with the AG genotype, indicating that the combined method of histone modification and GWAS can effectively identify SNPs located in the regulatory regions of the genome and reduce the probability of false positives and false negatives.

[0040] Example 2 A method for breeding Litopenaeus vannamei by using epigenetic annotation information and GS joint analysis comprises the following steps: S1. The phenotypic data of individual Penaeus vannamei and the total DNA of muscle tissue were collected as in step S1 of Example 1. The total DNA of muscle tissue of Penaeus vannamei was obtained and stored at -20°C for later use. The concentration and purity of DNA were tested by a nucleic acid analyzer. The OD value of the total DNA of muscle tissue of 972 Penaeus vannamei individuals was 260 / OD 280 All of them were between 1.8 and 2.0, indicating that total DNA with good purity was obtained. The total DNA extraction was complete as detected by 1.0% agarose gel electrophoresis and could be used for subsequent experiments.

[0041] S2. Construction of DNA library and sequencing, as in step S2 of Example 1, the sequencing data was compared with the unpublished reference genome of the chromosome of Penaeus vannamei in this laboratory for SNP typing, and the haplotype library of this laboratory was used for filling (filling missing genotypes) using Beagle v5.5 software, and SNPs with a minimum allele frequency MAF less than 0.05 and a Hawes-Wenzhou equilibrium test HWE P value less than 0.0001 were filtered, and a total of 2,609,549 qualified SNP sites were obtained.

[0042] S2. GS analysis, the specific reference steps are as follows: (1) The heritability estimation results based on the GCTA algorithm and the 2,609,549 high-quality SNPs screened by sequencing showed that the body length and weight traits of Penaeus vannamei had medium to high heritability (body length 0.44, body weight 0.29), among which the genetic contribution of body length was significantly higher than that of body weight. The two are highly genetically correlated, and the heritability of body length is higher. Therefore, the body length phenotype was selected as the research indicator of growth traits for GS analysis (Plink v1.9 software).

[0043] The genotype and body length phenotype data of Penaeus vannamei were imported into the GBLUP model, and the prediction accuracy was calculated by using the sommer R package in the R language in the GBLUP model (Genomic Best Linear Unbiased Prediction, a best linear unbiased prediction method based on genomic information, which directly estimates the breeding value of individuals by constructing a genomic relationship matrix (G matrix) to replace the traditional kinship matrix (A matrix); or the prediction accuracy was calculated by using the BGLR R package in the BayesA, BayesB or BayesC model. The higher the correlation coefficient, the better the prediction ability of the method. The above prediction process randomly and uniformly extracts SNP markers of different densities, and each density is randomly extracted 5 times. Different models are used for five-fold cross validation and repeated 5 times to systematically evaluate the prediction performance of each model under different marker densities.

[0044] (2) The prediction results are as follows: Figure 6 As shown in the figure, with body length as the phenotype, GS analysis was performed on 2,609,549 qualified SNP sites of 972 different groups of vannamei shrimp. With the gradual increase of marker density, the prediction accuracy of the mainstream breeding models such as GBLUP, BayesA, BayesB and BayesC for the breeding value of traits with body length as the phenotype is also gradually increasing. When the marker density reaches 10k, the prediction performance of the model breeding value tends to be flat, among which BayesA has the highest prediction accuracy, and the prediction accuracy is 0.34±0.01 (mean±standard deviation) at 10k.

[0045] S3. ChromHMM analysis, the specific reference steps are as follows: (1) Well-developed Litopenaeus vannamei embryos at critical embryonic developmental stages were randomly selected for CUT&Tag library construction, covering a total of seven developmental stages (blastocyst stage, gastrula stage, limb bud larvae, intramembranous larvae, nauplii stage I, nauplii stage III, and nauplii stage VI). The CUT&Tag library was constructed using the Illumina HyperactiveUniversal CUT&Tag Assay Kit (Vazyme, China, #TD903). The sample was prepared using the nuclear extraction method. The basic principles and operation steps of the experiment were referred to CUT&Tag (Kaya-Okur HS et al. CUT&Tag for efficient epigenomic profiling of small samples and single cells. Nat Communications, 2019, full text), and the Litopenaeus vannamei embryo samples for library construction were obtained.

[0046] (2) The vannamei embryo samples obtained for library construction were sequenced on the Illumina NovaSeq 6000 platform (Novogene, Tianjin). Sequencing data were obtained using Nova-PE150 paired-end sequencing (1G raw data / library). The sequencing data were aligned with the unpublished vannamei chromosome reference genome of our laboratory using the BWA-MEM algorithm (a new alignment algorithm used to align sequencing reads or assembled contigs to a large reference genome), and peak calling (a computational method used to identify regions in the genome where alignment reads obtained by sequencing are enriched) was performed using MACS2. The key parameters include: macs2 callpeak -t A.bam -c IgG.bam -f BAMPE--outdir ... / ... -B -n A -g 1.66e9 --keep-dup auto -q 0.05, the obtained bam file was subjected to chromHMM analysis (a computational tool based on Hidden Markov Model (HMM), which is used to systematically identify and characterize chromatin states in the genome, annotate their distribution locations across the genome, and assist the biological interpretation of these states by automatically calculating the enrichment of external annotations). The key parameters include java -mx4000M-jar ChromHMM.jar BinarizeBam -gzip -b 200 -f 0 -g 0 -p 0.0001 DUIXIA_chrom_sizes.txt bam sheet.txt binarization, java -mx4000M -jar ChromHMM.jarLearnModel -gzip -d 0.001 -color 129,0,0 -p 20 -i chrhmm binarizationlearnmodel 10 dx.

[0047] (3) The classification of chromatin states requires information about histone modifications, and there are many histone modifications with known functions. ChromHMM was used to analyze histone modifications in the key early embryonic development stage of Litopenaeus vannamei. Using the signal tracks marked by four histone modifications (H3K4me1: a mark of enhancer regions, H3K4me3: commonly found in promoter regions, H3K27me3: usually associated with gene expression repression, and H3K27ac: a marker of chromatin regions in an open state), ChromHMM first learned and defined a series of chromatin states, and then assigned each segment in the genome to one of these states. The LearnModel function was used to finally identify eight chromatin states, which were named according to the ChromHMM protocol (reference: Ernst J et al. Chromatin-state discovery and genome annotation with ChromHMM. Nature protocols, 2017, full text).

[0048] The results are as follows Figure 7 As shown in the figure, using the signal tracks of the four key histone modification regions H3K4me1, H3K4me3, H3K27me3, and H3K27ac in the genome during the critical period of embryonic development, a total of eight chromatin state information, including active transcription, enhancer, bivalent, repression, and quiescence, were annotated at the position 3506745-3731315 of chromosome 1 of Penaeus vannamei (chr1: 3506745-3731315), namely: bivalent TSS / enhancer region 1 (BivFlnk1), quiescent / low activity (Quies), weakly repressed by Polycomb (ReprPCWk), repressed by Polycomb (ReprPC), weak enhancer (EnhWk), bivalent TSS / enhancer region 2 (BivFlnk2), weak transcription (TxWk), and downstream TSS region (TssFlnkD), which can be visualized using IGV software.

[0049] S4, ChromHMM and GS joint analysis, the specific reference steps are as follows: (1) The match parameter in the R script was used to screen SNPs located in the peak interval of the ChromHMM bed file, and the --extract parameter of plink was used to retain the SNP sites in the genotype file. Then, the --thin-count parameter of plink was used to randomly extract SNP sites of different densities, and the --recodeA parameter of plink was used to convert the genotype file into a 012 format file. The prediction accuracy of the GBLUP model was calculated using the sommer R package, and the prediction accuracy of the BayesA / B / C model was calculated using the BGLR R package.

[0050] (2) The accuracy of the prediction was determined by using the sommer R package in the R language in the GBLUP model, or by using the BGLR R package in the BayesA, BayesB, or BayesC model. The higher the correlation coefficient, the better the predictive ability of the method. The above prediction process was performed by randomly and uniformly extracting SNP markers of different densities, with each density randomly extracted 5 times, and five-fold cross-validation was performed using different models and repeated 5 times to systematically evaluate the predictive performance of each model under different marker densities.

[0051] The results are as follows Figure 8 As shown in the figure, the epigenetic annotation information is used to carefully distinguish the chromatin state of the genome to screen SNPs. With the gradual increase in the marker density of SNPs located in different chromosome states during the key embryonic development period, the prediction accuracy of the BayesA model for the breeding value of the phenotype of body length is also gradually increasing. When the marker density is in the range of 1k~2.5k, the prediction accuracy of SNPs located on the E1, E5, E6 and E8 state Peaks is generally 6.25%~11.11% higher than that of randomly selected SNPs, and the prediction ability is more stable.

[0052] It can be concluded that the GWAS method provided by the present invention for screening significant SNPs and gene positioning, and applying histone modification information to the whole genome association analysis of growth traits of Penaeus vannamei, is more conducive to screening SNPs and gene positioning related to traits, can effectively identify SNP sites located in the genome regulatory region, and obtain expression regulatory element information of key SNPs related to traits, reduce the probability of false positives and false negatives, and is more conducive to a deeper study of its molecular mechanism. It provides more accurate molecular markers for the genetic analysis of growth traits of Penaeus vannamei, lays a foundation for genetic breeding research of Penaeus vannamei and other economic traits, and accelerates the breeding process.

[0053] The terms and expressions used herein are for descriptive purposes only and the present invention should not be limited to these terms and expressions. The use of these terms and expressions does not mean to exclude any equivalent features of the illustrations and descriptions (or parts thereof), and it should be recognized that various modifications that may exist should also be included in the scope of the claims. Other modifications, changes and substitutions may also exist. Accordingly, the claims should be deemed to cover all such equivalents.

[0054] Similarly, it should be pointed out that although the present invention has been described with reference to the current specific embodiments, ordinary technicians in this technical field should realize that the above embodiments are only used to illustrate the present invention, and various equivalent changes or substitutions may be made without departing from the spirit of the present invention. Therefore, as long as the changes and modifications to the above embodiments are within the scope of the essential spirit of the present invention, they will fall within the scope of the claims of the present invention.

Claims

1. A method for screening molecular markers related to economic traits of Penaeus vannamei, characterized in that: The following steps are involved: Step A: Perform phenotypic measurement of growth traits and genome sequencing on the breeding population of Penaeus vannamei to obtain phenotypic data and SNP typing data; Step B: performing GWAS analysis using the phenotypic data and the SNP typing data to obtain a P value for each SNP; Step C: Using the P value of each SNP and the Bonferroni multiple correction method, with P value <1.9E-08 as the significant site judgment, SNPs significantly associated with growth traits were screened out; using key histone modification regions to accurately locate SNPs significantly associated with growth traits and the expression regulatory elements of their corresponding genes.

2. The method for screening molecular markers related to economic traits of Penaeus vannamei according to claim 1, characterized in that: The step A comprises: Step A-1: ​​performing phenotypic measurement of growth traits on the breeding population of Penaeus vannamei to obtain the phenotypic data; the growth trait is body length or weight; Step A-2: extracting total DNA from the Litopenaeus vannamei breeding population to obtain the total DNA; constructing a DNA library using the total DNA; sequencing the DNA library to obtain the SNP typing data; The construction process of the DNA library includes: breaking the total DNA into small fragments by ultrasonic treatment, and undergoing end modification, linker connection, PCR amplification and magnetic bead purification.

3. The method for screening molecular markers related to economic traits of Penaeus vannamei according to claim 1, characterized in that: The step C comprises: Step C-1: sort the P value of each SNP in ascending order, and use the Bonferroni multiple correction method to judge the P value <1.9E-08 as a significant site; screen out SNPs that are significantly associated with growth traits; Step C-2: Use CUT&Tag to screen the key histone modification regions, and use the key histone modification regions to accurately locate SNPs significantly associated with growth traits and the expression regulatory elements of their corresponding genes; the key histone modification regions are H3K4me1, H3K4me3, H3K27me3 and H3K27ac.

4. A method for breeding Penaeus vannamei based on whole genome data, characterized in that: The following steps are involved: Step A: Perform phenotypic measurement of growth traits and genome sequencing on the breeding population of Penaeus vannamei to obtain phenotypic data and SNP typing data; Step B: performing GWAS analysis using the phenotypic data and the SNP typing data to obtain a P value for each SNP; Step C: using the P value of each SNP to undergo a Bonferroni multiple correction method, with a P value <1.9E-08 as a significant site judgment, to screen out SNPs significantly associated with growth traits; using key histone modification regions to accurately locate SNPs significantly associated with growth traits and expression regulatory elements of their corresponding genes; Step D: determining the preferred genotype in the breeding population of Penaeus vannamei according to the SNP locus subset, and breeding by retaining the parents with the preferred genotype for pairing and breeding.

5. An application of SNP molecular markers in assisted breeding of Penaeus vannamei, characterized in that: The sequence of the SNP molecular marker is shown in SEQ ID NO.1; the SNP site is located at the 233rd position from the 5' end of the sequence shown in SEQ ID NO.1; the polymorphism of the SNP molecular marker is A / G type; and the application of the SNP molecular marker in assisted breeding of Penaeus vannamei is the screening of the body length or weight of Penaeus vannamei.

Images