Short tandem repeat variation detection method based on HMM-GMM

Abstract

The invention provides a short tandem repeat (STR) variation detection method based on HMM-GMM, which belongs to the technical field of biological information, integrates genomic sequences, known STR site information and three-generation sequencing comparison results, and realizes high-precision positioning by using the advantage of long read length. Hidden state transition of an HMM capture sequence and complex distribution of GMM modeling are combined to improve variation prediction accuracy. According to the method, the STR region is directly compared with a reference genome, errors are reduced by combining a CIGAR value, the number of repetitions is quantified, specific base variation is detected, and pathogenic STR variation is recognized. Compared with traditional software, the output image shows the variation condition in more detail, and a powerful tool is provided for genetic disease diagnosis and gene research.

Description

Technical Field

[0001] The present invention relates to the field of bioinformatics, in particular to a method for detecting genomic structural variation, and specifically to a method and system for detecting short tandem repeat sequence variation based on HMM-GMM. Background Art

[0002] Short tandem repeats (STRs), also known as microsatellites, are a type of DNA sequence that is widely present in the genome and consists of 1 to 6 nucleotides repeated multiple times as basic units. These repeating units are usually arranged in series in a head-to-tail manner to form a longer DNA sequence. These sequences account for 3% of the total human genome. Since the length of STRs is highly polymorphic between different individuals (i.e., the length of STR sequences may vary between different individuals), they have a wide range of application value in genetics, forensics, population genetics, and disease research.

[0003] Short tandem repeat sequence variation refers to changes in the length or structure of the repeating units of short tandem repeat sequences in the genome. STRs variation is usually manifested as changes in STR length caused by insertion and deletion of repeating units. The most common cause of this change is the "slippage" phenomenon during DNA replication, that is, during DNA replication, DNA polymerase mismatches when encountering the STR region, resulting in some repeating units not being correctly replicated, thereby causing changes in the number of repeats in the daughter DNA molecules. The high variability of STRs means that they occupy an important position in human genetic diversity and may contribute more new mutations than other types of genetic variation. STRs variation is closely related to the occurrence of a variety of genetic diseases, such as Huntington's disease, fragile X syndrome, etc. Therefore, accurate identification and analysis of STRs variation is crucial to revealing their biological functions and medical significance.

[0004] Traditional molecular biology methods usually use polymerase chain reaction (PCR) technology to amplify specific regions in the genome, and then use capillary electrophoresis, gel electrophoresis and other technologies to determine the number of STRs repeats. This method is easy to operate, but the detection throughput is limited; Sanger sequencing, as a representative of first-generation sequencing technology, has high sensitivity and specificity and is regarded as the international standard for genetic testing. However, its shortcomings in high-throughput applications limit its application in large-scale genome sequencing projects. Second-generation sequencing technology has significantly improved the coverage and efficiency of whole-genome STRs detection, but due to read length restrictions (usually 100 to 150 base pairs), it has limitations in detecting long STRs amplification. Tools that have complex modeling for errors generated by PCR reactions in sequencing (such as RepeatSeq, STRviper, HipSTR, GangSTR, Tredparse) can improve detection accuracy, but are still limited by the length of sequencing data.

[0005] With the rapid development of high-throughput sequencing technology, especially the emergence of third-generation sequencing technology, the detection ability of STRs has been significantly improved. For example, SMRT and Nanopore sequencing have broken through the limitations of second-generation sequencing technology with their ultra-long read lengths, effectively overcoming the noise problems in the short read length and PCR amplification process. A variety of detection tools (such as RepeatHMM, TandemGenotypes, RepLong, STRique, TRiCoLOR, DeepRepeat, warpSTR, TRGT) have improved the detection ability of STRs variation using third-generation sequencing. Although third-generation sequencing technology has advantages in detecting long STRs, it still faces challenges when dealing with STR variations in highly repetitive sequences, low sequencing depth or complex gene regions. In addition, the high base recognition error rate of third-generation sequencing technology may also affect the detection accuracy of STRs variation. Therefore, the development of more accurate and efficient STRs variation detection algorithms and tools is an important direction of current research. Summary of the invention

[0006] In view of the problems of false positives and high base recognition error rate in existing detection methods when dealing with STR variations in highly repetitive sequences, low sequencing depth or complex gene regions, how to improve the accuracy and range of STRs variation detection and ensure the accuracy and reliability of STRs variation analysis in genetic disease research, the present invention aims to provide a short tandem repeat sequence variation detection method and system based on HMM-GMM (hidden Markov model-Gaussian mixture model).

[0007] In order to achieve the above object, the technical solution adopted by the present invention is: The present invention provides a method for detecting short tandem repeat sequence variation based on HMM-GMM, comprising: S1, file preprocessing, the files include genome annotation files, reference genome sequence files and sequence alignment files obtained by third-generation sequencing technology, and the formats of input files are bed, vcf, fasta and bam respectively; S2, identifying and locating the short tandem repeat sequence repeat region, including the start and end position information of the short tandem repeat sequence repeat unit; S3, inputting the repeat region information of the short tandem repeat sequence located in S2 and the sequence of the short tandem repeat sequence into the pre-trained HMM-GMM model to predict the type and length of the short tandem repeat sequence variation in the repeat region; S4, cluster the prediction results of S3, and divide the short tandem variants into homozygous and heterozygous according to the clustering results; S5, use visualization tools to display short tandem mutation information in the form of images.

[0008] The S1 includes: Preprocessing of genome annotation files, filtering out STR sites located in other repetitive regions (such as segmental duplications and transposons), and screening for short tandem sites with different repeat unit lengths; Standard variant vcf,file preprocessing, screening variant sites located in microsatellite regions, extracting short tandem variant information, and calculating the number of variants.

[0009] The identifying and locating of the short tandem repeat sequence repeat region includes: locating the short tandem repeat sequence repeat region in the genome, comparing the reference genome with the sequencing data, and the position of the short tandem repeat sequence repeat region in the sequencing data.

[0010] Furthermore, the positioning of the short tandem repeat sequence repeat region in the genome is combined with the TRFrepeatmasker and MSIsensor-pro methods; the reference genome and test data are aligned with the pbmm2, ngmlr and minimap2 tools to generate a SAM or BAM format alignment result file; the pysam module of the python library is used to achieve the positioning of the short tandem repeat sequence repeat region in the sequencing data.

[0011] The S3 includes: S31, using an EM (Expectation-Maximization) algorithm of a hidden Markov model to evaluate the hidden Markov model parameters and obtain optimal hidden Markov model parameters; S32, using the Viterbi algorithm for prediction, and correcting the prediction results according to the previously verified knowledge match (Match) > deletion mutation (Deletion) > insertion mutation (Insertion); S33, prediction of short tandem repeat lengths for deletion variant status using Gaussian mixture models; S34, combining the sequence length of the deletion variation state predicted by the hidden Markov model in S33 with the prediction result of the Viterbi algorithm in S32 to obtain the type and length prediction result of the short tandem repeat sequence variation in the repeat region.

[0012] Furthermore, the S33 constructs a data set based on the short tandem repeat sequence and the CIGAR value and trains a Gaussian mixture model to predict the deletion length.

[0013] The S4 calculates the support rate of each cluster. If the support rate is greater than 80%, the short tandem repeat sequence site is regarded as homozygous; if there is no cluster with a support rate greater than 80%, the short tandem repeat sequence site is regarded as heterozygous.

[0014] The present invention provides a system for implementing the above-mentioned HMM-GMM-based short tandem repeat sequence variation detection method, comprising: Data acquisition module, used to obtain genome annotation files, reference genome sequence files, and sequence alignment files obtained through third-generation sequencing technology; Data preprocessing module, used for file format conversion, data cleaning, screening, etc.; A short tandem repeat sequence recognition and positioning module, used to identify and locate the repetitive regions of short tandem repeat sequences in the genome; The model prediction module is used to input the sequence and location information of short tandem repeats into the pre-trained HMM-GMM model to predict the variant type and length; The result clustering and classification module is used to cluster the model prediction results and classify the short tandem variants into homozygous and heterozygous according to the clustering results; The result visualization module is used to display the short tandem mutation information in the form of images; Model evaluation module, used to verify and evaluate the model.

[0015] The present invention provides a computer device, comprising a memory and a processor, wherein the memory stores a computer program, and the processor implements the steps of the above-mentioned HMM-GMM-based short tandem repeat sequence variation detection method when executing the computer program.

[0016] The present invention provides a computer-readable storage medium, wherein the computer-readable storage medium stores a computer program, and when the computer program is executed by a processor, the steps of the above-mentioned method for detecting short tandem repeat sequence variation based on HMM-GMM are implemented.

[0017] Compared with the prior art, the beneficial effects of the present invention are: The present invention provides a method for detecting short tandem repeat sequence variation based on HMM-GMM, which shows remarkable technical effects in detection accuracy, comprehensiveness and meticulousness of variation prediction, automation and systematization, as well as clinical application and genetic research support. The method combines genome annotation, reference genome sequence and sequencing data, uses HMM-GMM model for variation prediction, and improves the accuracy and efficiency of STRs variation detection; through file preprocessing, short tandem repeat region identification, variation prediction, result clustering and visualization display, the complete process from data input to result output is realized; compared with traditional methods, the method emphasizes the diversity and integrity of input files, which not only includes genome sequence files and known STRs site information, but also includes sequence comparison results obtained by third-generation sequencing technology. The integration of such multi-source data provides a richer and more accurate information basis for subsequent variation detection, and utilizes the long read length advantage of third-generation sequencing technology (such as PacBio, Nanopore) and combines known STR site information to achieve high-precision positioning of short tandem repeat regions. Compared with short-read sequencing technology, third-generation sequencing can span the complex repetitive structure of the STR region, reducing splicing errors and missed detections; the combined model of HMM and GMM is introduced for variation prediction. HMM can capture the hidden state transition rules in the sequence, while GMM is good at modeling complex distributions. The combination of the two enables the model to more accurately describe the repetitive patterns and internal variations of STRs sequences. It can not only quantify the number of repetitions of STRs variations, but also detect the variation type at the specific base position. Considering the length factor of the missing variation, a large number of amplified pathogenic STRs variations can be identified; the accuracy and reliability of variation prediction are improved.

[0018] Furthermore, non-STR sites and STR sites located in other repetitive regions were filtered out to improve the pertinence and accuracy of subsequent analysis, short tandem sites with different repeat unit lengths were screened to increase the comprehensiveness and meticulousness of the analysis, variant sites located in microsatellite regions were screened, STR variant information was extracted, and high-quality data was provided for subsequent model prediction; the MSIsensor-pro method was used to locate STRs in the genome, and its advantages in STR detection were utilized to improve the accuracy of positioning; the minimap2 tool was used to compare the reference genome and sequencing data, and high-quality SAM or BAM files were generated, and the STR repeat region was realized through the pysam module. The precise positioning of the domain in the sequencing data provided reliable data support for subsequent analysis; the EM algorithm of HMM was used to evaluate the model parameters to ensure that the model could accurately reflect the inherent laws of the data; the Viterbi algorithm was used to predict and obtain preliminary results of STR variation; the correction of prediction results and the use of GMM model further improved the accuracy and reliability of variation prediction; the output results of HMM and GMM were combined to obtain more comprehensive and accurate predictions of STR variation type and length; the data set was constructed according to STR sequences and CIGAR values ​​to ensure the accuracy and representativeness of the training data; the GMM model was used to predict the missing length to improve the accuracy and applicability of the prediction results.

[0019] The system provided by the present invention realizes the automation and systematization of STRs variation detection through modules such as data collection, data processing and model evaluation, improves the detection efficiency and accuracy, reduces the errors caused by human intervention, and provides strong support for genetic research and clinical applications. Compared with the commonly used software IGV, the images output by the present invention can show the variation in more detail, which is convenient for researchers to intuitively understand and analyze STRs variation; by accurately detecting STR variation, it can help researchers better understand the pathogenesis of genetic diseases and provide an important basis for early diagnosis, prevention and treatment of diseases. BRIEF DESCRIPTION OF THE DRAWINGS

[0020] Figure 1 is a flow chart of the method of the present invention; Figure 2 It is the visualization output of the variation results detected by the HMM-GMM model of the present invention. DETAILED DESCRIPTION

[0021] In order to enable those skilled in the art to better understand the scheme of the present invention, the technical scheme in the embodiments of the present invention will be clearly and completely described below in conjunction with the drawings in the embodiments of the present invention. Obviously, the described embodiments are only part of the embodiments of the present invention, not all of the embodiments. Based on the embodiments of the present invention, all other embodiments obtained by ordinary technicians in this field without creative work should fall within the scope of protection of the present invention.

[0022] The present invention is described in further detail below: Example 1 Step 1: File preprocessing Step 1: Preprocess the genome annotation file, filter out STR loci in other repetitive regions, and screen STR loci with different repeat unit lengths. When filtering STR loci in other repetitive regions, mainly consider simple repeats (SR), segmental duplications (SD), and variable number of tandem repeats (VNTR). The positions of STR loci in the above three types of repetitive regions can be divided into three categories: The SRT portion spans other repeat regions; STRs completely span other repetitive regions; The STR is completely covered by other repetitive sequences.

[0023] Filter the bed file containing STR locus information, calculate the maximum length leni of each chromosome chri, and create two arrays of length leni for each chromosome to record the number of three types of repeat regions str_doct_i[] and the position information of each data point in the genome str_pos_i[]. Traverse the bed file, record the sites where SR, VNTR, and SD repeats appear, assign 1 to the corresponding data point value in the str_doct_i[] array, and record the position information of the site in the array str_pos_i[]. Secondly, sum all sites. If sum is greater than 0, it is considered that the STR site overlaps with other repeat regions, indicating that the STR site is located in other repeat regions and the site is filtered out. Otherwise, the site is retained. During the processing, when the bed file contains multiple STR repeat sites of different lengths, only the repeat units of the set length are considered.

[0024] Step 2: Preprocess the standard variant VCF file, filter out variant sites located in microsatellite regions, and extract STR variant information. Filtering variant sites is similar to filtering STR sites in other repetitive regions in the first step. The processed bed file is a reference file containing STR repeats. Traverse the variant sites in the standard VCF. Similar to the first step, if the sum is greater than 0, it means the variant site is located in the microsatellite, and this site is retained; otherwise, it is filtered out to obtain the pre-screened vcf. First, compare whether the REF part overlaps with the STR region. If there is an overlap, extract the STR variant information, including insertion, substitution, and deletion information, as well as the lengths of REF and ALT at each site: If ALT > REF, an insertion mutation has occurred, and it is considered a deletion mutation with the exchange of ALT and REF; If ALT = REF, a substitution mutation has occurred, which does not belong to STR mutation and is directly filtered; If ALT < REF, a deletion mutation has occurred. Extract the bases that are different between ALT and REF. If the position of a single base overlaps with the SRT site where the variant site occurs, and the genotype of this base is in the mofit genotype of the repeat unit of the corresponding STR site, the number of STR site mutations is incremented by 1. The mutation situation is statistically analyzed according to the following calculation method:

[0025] Among them, Tvar refers to the number of mutated repeats, Nvar refers to the number of statistically analyzed mutations, and Lmotif refers to the length of the STR repeat unit.

[0026] Step 2: Locate the short tandem repeat region and determine the position of the STR repeat region in the sequencing data The localization of short tandem repeat regions mainly includes the localization of STR repeat regions in the genome, the comparison of reference genome and sequencing data, and the position of STR repeat regions in sequencing data. The MSIsensor-pro method is used to realize the localization of STR repeat regions in the genome and generate genome annotation bed files. Minimap2 is used to realize the comparison of reference genome and test data, and then converted into bam files using samtools, and indexes are created. The pysam module of python is used to realize the localization of STR repeat regions in sequencing data, extract the base sequence read_str of sequencing data and the alignment value CIGAR of the sequence ("=" indicates a complete match, "M" indicates a match, "X" indicates a sequence mismatch, "I" indicates an insertion, "D" indicates a deletion, and "N" indicates skipping the reference sequence), and process it into a string sequence read_list that corresponds to the reference genome one by one. For example, read_str = "GCATATAATACG", CIGAR = "6MI3M1D2M", converted into read_list = [G, C, A, T, A, TA, A, T, A, '', C, G].

[0027] Step 3: Quantification and prediction of short tandem repeat sequence variation Use the Markov model to predict short tandem repeat mutations. The STR mutation detection problem is regarded as a state sequence prediction problem in a first-order hidden Markov model. The predicted state sequence contains specific mutation information for each base position. In addition, the hidden state consists of deletion (D), insertion (I), and match (M). The state of any position only depends on the state of the previous position. The observed sequence corresponds to the STR repeat sequence in the sequencing data (consisting of A, T, C, and G). The base at any position only depends on the state of the corresponding position in the state sequence. The element λ of the HMM model is (A, B, π), where the state transition probability matrix A represents the transition probability between different hidden states, and B represents the probability of a state generating different bases. Here we express the problem in mathematical form: For a STR repetitive sequence, it contains N hidden states Q={q1,q2,...,qn}, M observation symbols V={v1,v2,...,vm}, the observation sequence S={s1,s2,...,sT} represented by the observation symbols, the state sequence I={i1,i2,...,in} corresponding to the observation sequence, the state transition probability matrix B=[P(st=vk | it=qj)]N×M, and the initial state probability matrix π=(P(i1=qi)). Hidden state Q={q1=M,q2=I,q3=D}, The EM algorithm of HMM unsupervised learning method is used to evaluate the model parameter λ, determine the number of hidden states (the present invention recommends 3 or 12), set the random seed to a fixed value, call the fit algorithm to start training, and obtain the trained model.

[0028] The Viterbi algorithm is used to implement the prediction of the HMM model and predict the variation status corresponding to each STR repeat sequence in the sequencing data.

[0029] The prediction results are corrected, mainly for the matching and insertion states. According to the prior knowledge of the probability of the hidden state M>D>I, the number of occurrences of each number in the prediction results is counted. According to the assumptions and the number from most to least, the different results are matched to the three situations of matching, deletion, and insertion in turn, so as to obtain the corrected state sequence. The variation of each base position can be predicted, and the predicted STR repeat sequence length can be calculated.

[0030] The GMM model is used to accurately predict the deletion state. The STR repeat sequence and the corresponding CIGAR value are used to cluster the lengths of the deletion variants in different situations. The diversity and inaccuracy of the length of the deletion variant are considered to make the prediction results closer to the actual situation. The input of the GMM model consists of multiple Gaussian distributions, and eight situations are considered, namely AM, TM, CM, GM, AI, TI, CI, and GI. AM indicates the deletion that occurs after the STR repeat sequence base is A and the CIGAR value is M. If two deletions occur after the base A, that is, MDD, 2 is added to the AM category. Finally, the above eight situations are clustered into two categories, and the one close to the probability of 1 is taken as the length of the missing variant in this case.

[0031] Step 4: Cluster the predicted short tandem repeat sequence results and cluster the scattered prediction results into homozygous and heterozygous The same STR repeat site may have multiple sequencing data, and multiple variant prediction results will be obtained. The probability model GMM model is also used to cluster the prediction results of each site into two categories: homozygous and heterozygous. If the support rate of a certain category is greater than 80%, the STR site is considered to be homozygous, and the sequence greater than or equal to 80% in the clustering result is determined as the final prediction result. If there is no sequence greater than or equal to 80%, the STR site is considered to be heterozygous, and the two clustering results are taken as the prediction result.

[0032] Step 5: Visualize the output results. The output VCF file is displayed as an image to show the mutation type and mutation length. The variation results predicted by the STR repeat region are stored in VCF format, including chromosome number, site information and detected variation information. The commonly used IGV tool cannot display many variation conditions in the sequencing data due to the limited window size. The output visualization image of the present invention can more clearly and accurately display the details of the STR repeat variation, directly displaying the bases of the STR repeat unit, the repeat unit bases of the STR site, the chromosome name and the start and end positions, etc.

[0033] Example 2 The effect of the present invention can be further demonstrated by the following simulation: Simulation experiment conditions The experimental platform of the present invention is Linux operating system, Python 3.10 and Rstudio integrated development software.

[0034] The required number of simulations comes from the sequence alignment bam file and standard VCF file in the standard data set provided by GIAB. The 7 standard sample data HG001~HG007, the VCF file of the standard benchmark comes from ftp: / / ftp-trace.ncbi.nlm.nih.gov / giab / ftp / release / , and the benchmark files of the 7 samples contain 3893341, 4048342, 4000097, 4031346, 3856856, 3839315, and 3859704 variant data respectively. The number of STR variants obtained after preprocessing are 20180, 23374, 22434, 22967, 15261, 16466, and 16623 respectively.

[0035] Experimental result evaluation index: The performance of the present invention is measured by comparing the experimental results with the gold standard, using four indicators: precision, recall, the proportion of variants that can be detected in VCF, and the proportion of variants that can be detected but are not large enough to be detected by VCF.

[0036]

[0037]

[0038]

[0039]

[0040] in, It indicates the number of results detected by the method and the detection is accurate; Indicates the number of results detected by the method and the number of errors; Indicates the number of results detected by the standard VCF file and not detected by the method; indicates the number of results detected by both the method and the standard VCF file; It indicates the number of results detected by the method and not detected by the standard VCF file. With these four evaluation indicators, we can conduct real data experiments and analyze the experimental results.

[0041] Simulation content and result analysis (1) Model training process. The present invention uses an unsupervised EM algorithm to evaluate the λ parameter of the hidden Markov model, λ = (A, B, Π), uses the sequencing data of the Sequek sequencing platform of the HG002 sample for preprocessing, and obtains a sequencing data length of 15kb. Then, the sequence is aligned by pbmm2 to obtain a bam file, and the 4056 STRs of chromosome 1 are repeatedly trained. The initial state parameters obtained are as follows:

[0042]

[0043]

[0044] (2) All Hifi sequencing data of samples HG001~HG007 were tested, among which HG002, HG003 and HG004 were from the same family trio, HG002 was a child, HG003 and HG004 were his biological father and biological mother respectively, while HG005, HG006 and HG007 were from the same family trio, HG005 was a child, HG006 and HG007 were his biological father and biological mother respectively. At the same time, for the bed file to be tested, the STR loci located in other repeat regions were not filtered, and every 10 loci were processed once, and 187322 STR loci were obtained. This method was used to predict the STR repeats of these 7 different human individual samples on a genome-wide scale. The four indicator results of each individual and the total number of STR variations detected are shown in the following table. As can be seen from Table 1, this method can be reliable on different human individual samples and can be widely applied to different human individuals.

[0045] Table 1: Genome-wide prediction of STR repeats from 7 different human individual samples

[0046] Secondly, we used the variation detection method based on the HMM-GMM model to detect the variation in the sequencing data of the first 10Mbp length of chromosome 1 in the third-generation data HG002 sample HiFi data. Compared with the TRGT detection tool, the number of STs detected and the four evaluation indicators are shown in Table 2.

[0047] Table 2: Comparison of the number of STs detected by the present invention and the TRGT detection tool and four evaluation indicators

[0048] From the data in Table 2, we can see that the number of STs (short tandem repeats) detected by HMM-GMM is 2758, while the number detected by TRGT is 2106, which indicates that the HMM-GMM model can identify more ST sites, which is attributed to its more refined model structure and stronger data fitting ability; The accuracy of HMM-GMM is 1, that is, all detected ST sites are correctly identified, and there are no false positives. The accuracy of TRGT is 0.9652, which is also very high, but slightly lower than HMM-GMM. This further proves the superiority of the HMM-GMM model in identifying true ST sites; The recall rate of HMM-GMM is also 1, which means that all true ST sites are successfully detected without omission. The recall rate of TRGT is 0.9916, but there are still a small number of true ST sites that are not detected. This once again shows the advantage of the HMM-GMM model in comprehensive detection; The internal consistency index measures the accuracy of the detection tool in identifying multiple ST sites in the same region. The ms_In of HMM-GMM is 1, indicating that it performs perfectly in this regard. The ms_In of TRGT is 0.9919, which is slightly inferior to HMM-GMM; In summary, the variation detection method based on the HMM-GMM model of the present invention shows higher performance in terms of ST quantity, precision, recall rate and internal consistency when detecting ST variation in third-generation sequencing data compared with the TRGT detection tool. These results not only verify the effectiveness of the HMM-GMM model in complex genomic data analysis, but also provide strong support for its wide application in genetic research and clinical diagnosis.

[0049] In order to evaluate the ability of the method to predict longer STR variations, in the prediction results of the whole genome of the HG002 sample, 13 STR sites with variation lengths greater than 500bp were detected, of which the maximum insertion variation length that could be detected was 3085bp. For the chr10_8806613 site of human chromosome 10 in the HG002 sample, the STR repeats were 11 [AAGG] when the site was not mutated, indicating that the number of repeats was 11 and the repeat unit was AAGG. The method was used to detect the variation of the chr10_8806613 site, and the prediction results obtained were 676.5 [AAGG] and 782.25 [AAGG], that is, it was predicted that a pair of alleles had insertion variations of 2662bp and 3085bp, respectively. In order to verify the accuracy of the predicted variation results, the IGV tool was used to visualize the distribution of sequencing data at this site, as shown in the attached figure. Figure 2 shown.

[0050] By the attached Figure 2 The data show that the prediction results obtained by this method are close to the actual situation. Therefore, the HMM-based detection method not only has the ability to detect highly amplified STR variants, but also its prediction results show high accuracy and reliability for STR variants at different positions, which provides a new and powerful tool for the detection of STR variants in genetic research and clinical diagnosis, and helps to more comprehensively understand the diversity and complexity of the human genome.

[0051] The above contents are only for explaining the technical idea of ​​the present invention and cannot be used to limit the protection scope of the present invention. Any changes made on the basis of the technical solution in accordance with the technical idea proposed by the present invention shall fall within the protection scope of the claims of the present invention.

Claims

1. A method for detecting short tandem repeat sequence variation based on HMM-GMM, characterized in that: include: S1, file preprocessing, the files include genome annotation files, reference genome sequence files and sequence alignment files obtained by third-generation sequencing technology, and the formats of input files are bed, vcf, fasta and bam respectively; S2, identifying and locating the short tandem repeat sequence repeat region, including the start and end position information of the short tandem repeat sequence repeat unit; S3, inputting the repeat region information of the short tandem repeat sequence located in S2 and the sequence of the short tandem repeat sequence into the pre-trained HMM-GMM model to predict the type and length of the short tandem repeat sequence variation in the repeat region; S4, cluster the prediction results of S3, and divide the short tandem variants into homozygous and heterozygous according to the clustering results; S5, use visualization tools to display short tandem mutation information in the form of images.

2. The method for detecting short tandem repeat sequence variation based on HMM-GMM according to claim 1, characterized in that: The S1 includes: Preprocessing and filtering of genome annotation files to screen for short tandem sites with different repeat unit lengths; Standard variant vcf files were preprocessed to screen variant sites located in microsatellite regions, extract short tandem variant information and calculate the number of variants.

3. The method for detecting short tandem repeat sequence variation based on HMM-GMM according to claim 1, characterized in that: The identifying and locating of the short tandem repeat sequence repeat region includes: locating the short tandem repeat sequence repeat region in the genome, comparing the reference genome with the sequencing data, and the position of the short tandem repeat sequence repeat region in the sequencing data.

4. The method for detecting short tandem repeat sequence variation based on HMM-GMM according to claim 3, characterized in that: The positioning of the short tandem repeat sequence repeat region in the genome is combined with the TRF repeatmasker and MSIsensor-pro methods; the reference genome and test data are compared with the pbmm2, ngmlr and minimap2 tools to generate a SAM or BAM format comparison result file; Use the pysam module of the Python library to locate the repeat regions of short tandem repeats in sequencing data.

5. The method for detecting short tandem repeat sequence variation based on HMM-GMM according to claim 1, characterized in that: The S3 includes: S31, using the EM algorithm of the hidden Markov model to evaluate the hidden Markov model parameters and obtain the optimal hidden Markov model parameters; S32, using the Viterbi algorithm for prediction, and correcting the prediction results according to the first verified knowledge matching > deletion variation > insertion variation; S33, prediction of short tandem repeat lengths for deletion variant status using Gaussian mixture models; S34, combining the sequence length of the deletion variation state predicted by the hidden Markov model in S33 with the prediction result of the Viterbi algorithm in S32 to obtain the type and length prediction result of the short tandem repeat sequence variation in the repeat region.

6. The method for detecting short tandem repeat sequence variation based on HMM-GMM according to claim 5, characterized in that: The S33 constructs a data set based on short tandem repeat sequences and CIGAR values ​​and trains a Gaussian mixture model to predict the deletion length.

7. The method for detecting short tandem repeat sequence variation based on HMM-GMM according to claim 1, characterized in that: The S4 calculates the support rate of each cluster. If the support rate is greater than 80%, the short tandem repeat sequence site is regarded as homozygous; if there is no cluster with a support rate greater than 80%, the short tandem repeat sequence site is regarded as heterozygous.

8. A system for implementing the method for detecting short tandem repeat sequence variation based on HMM-GMM according to any one of claims 1 to 7, characterized in that: include: Data acquisition module, used to obtain genome annotation files, reference genome sequence files, and sequence alignment files obtained through third-generation sequencing technology; Data preprocessing module, used for file format conversion, data cleaning, screening, etc.; A short tandem repeat sequence recognition and positioning module, used to identify and locate the repetitive regions of short tandem repeat sequences in the genome; The model prediction module is used to input the sequence and location information of short tandem repeats into the pre-trained HMM-GMM model to predict the variant type and length; The result clustering and classification module is used to cluster the model prediction results and classify the short tandem variants into homozygous and heterozygous according to the clustering results; The result visualization module is used to display the short tandem mutation information in the form of images; Model evaluation module, used to verify and evaluate the model.

9. A computer device comprising a memory and a processor, wherein the memory stores a computer program, wherein: When the processor executes the computer program, the steps of the method for detecting short tandem repeat sequence variation based on HMM-GMM according to any one of claims 1 to 8 are implemented.

10. A computer-readable storage medium storing a computer program, characterized in that: When the computer program is executed by a processor, the steps of the method for detecting short tandem repeat sequence variation based on HMM-GMM according to any one of claims 1 to 8 are implemented.

Images