Method and system for predicting pre-mRNA shearing based on deep learning algorithm

Abstract

The invention relates to the technical field of RNA (Ribonucleic Acid) cleavage and deep learning, in particular to a method and a system for predicting pre-mRNA cleavage based on a deep learning algorithm, which can solve the problems that the existing method is relatively weak in recessive non-classical cleavage site prediction capability, not ideal in accuracy on a real data set and relatively slow in prediction speed to a certain extent. The method comprises the following steps: acquiring pre-mRNA sequence data based on a public database and an RNA-seq sequencing data set, and calibrating a cleavage site and a non-cleavage site in the pre-mRNA sequence data at the same time; preprocessing the pre-mRNA sequence data to obtain a data set required by model generation, the data set comprising a training set, a verification set and a test set; building a neural network model combining a convolutional neural network and a self-attention mechanism, and training the neural network model by using the training set and the verification set to obtain an optimal neural network model; and performing three-classification prediction on the test set by adopting the optimal neural network model to obtain the probabilities of donor / receptor shear sites and non-shear sites.

Description

Technical Field

[0001] The present application relates to the field of RNA splicing and deep learning technology, and specifically, to a method and system for predicting pre-mRNA splicing based on a deep learning algorithm. Background Art

[0002] RNA splicing is a complex biological process, which mainly refers to the process of removing introns from pre-mRNA transcribed from a DNA template chain, splicing exons together, and finally generating mRNA. The transcribed product is translated into corresponding proteins to play its specific function. It can be seen that the correct splicing of RNA plays a key role in the regulation of gene expression. The precise splicing of RNA mainly depends on the regular base sequence at the junction of exons and introns, generally the double base GU at the splicing donor (donor, starting from the 5' end of the intron) and the double base AG at the splicing acceptor (acceptor, starting from the 3' end of the intron). In addition, other motifs located on the pre-mRNA, such as splicing branch points, polypyrimidine tracts, and splicing silencers and enhancers, will participate in the entire RNA splicing process, such as Fig. 9 shown.

[0003] RNA splicing is a key step in the process of gene transcription regulation. RNA splicing refers to the process in which discontinuous genes encoding proteins remove introns from pre-mRNA in sequence through RNA splicing and splice exons into mature mRNA. Selective splicing determines the expression of proteins in different tissues and the exercise of their specific functions. It has been found that different abnormal RNA splicing can be caused according to the different locations of gene mutations in the DNA sequence, which are divided into exon skipping, inclusion of an intron fragment, inclusion of a cryptic exon, and loss of an exon fragment, which further induces abnormal function of the translated protein product and ultimately leads to the occurrence of diseases. For example, the c.892A>G mutation of the gene MAPT reduces RNA splicing, causing exon10 skipping, leading to the occurrence of frontotemporal dementia with Parkinson's disease; the c.2204+6T>C mutation of the gene ELP1 will cause the skipping of exon 20 of chromosome 9, leading to familial autonomic dysfunction, etc.

[0004] At present, the main experimental methods used to study abnormal RNA splicing include transcriptome sequencing, RT-PCR, etc. These traditional wet experiments are mostly limited by the source of samples, time-consuming and costly. Traditional machine learning methods developed based on logistic regression, decision trees and support vector machines, such as GeneSplicer, SpliceMachine and MaxEntScan models, have certain shortcomings. For example, linear regression is very poor in dealing with nonlinear relationships, is not flexible enough in identifying complex patterns, and has low classification accuracy. Decision trees are prone to overfitting and ignore the correlation between attributes in the data set. Support vector machines are not suitable for processing a large number of samples and have no general solution to nonlinear problems. As a result, the features extracted by the model are limited by the existing knowledge of splicing sites, and the prediction ability of non-classical splicing sites is poor, and the classic double bases can only be predicted in a small range of sequences.

[0005] In recent years, with the development of deep learning technology, researchers have used deep learning algorithms to learn the characteristics of gene transcripts in public databases and directly predict splicing sites. For example, SpliceAI invented by Kishore Jaganathan et al. uses a deep convolutional residual network to combine a convolutional structure with a depth of 32 and a residual structure, and trains on sequences of lengths of 80nt, 400nt, 2000nt, and 10,000nt, respectively. It was found that the longer the sequence length used for training, the higher the accuracy of predicting splicing sites, with the accuracy rate rising from 57% for 80nt to 95% for 10,000nt. This is also the first time that predictions have been made within the range of tens of thousands of bases and a relatively ideal accuracy rate has been achieved.

[0006] However, since most existing splicing methods use known splicing sites as training sets, they have weak prediction capabilities for hidden non-classical splicing sites, their accuracy on real data sets is not ideal (SpliceAI is only 75%), and their prediction speed is slow.

[0007] In addition, in the field of RNA shearing, ASO, as a short-chain DNA or RNA molecule antisense oligonucleotide drug, has become a research hotspot in recent years due to its advantages in precise gene regulation, treatment of genetic diseases and rapid development cycle. In particular, ASO based on the steric blocking mechanism can be used to design interventions in the RNA shearing process, affect the selection of shearing sites, and thus change the RNA shearing pattern.

[0008] However, the design of ASOs based on steric hindrance mechanisms faces multiple challenges, including ensuring specificity and reducing off-target effects, improving stability and delivery efficiency, considering the effects of RNA secondary structure, and determining the optimal dose to avoid toxicity. Effectively solving these problems is the key to the successful application of ASOs.

[0009] Deep learning technology can provide support for ASO design. Through rapid screening and high-accuracy prediction, deep learning can significantly speed up the ASO screening process, save R&D costs, and increase the success rate. At the same time, deep learning technology can also be used to personalize ASO therapy and improve treatment effects. Although deep learning models have great potential in ASO design, there is still no suitable deep learning model that can effectively predict and screen the effects of ASO. Summary of the invention

[0010] In order to solve the problems that most of the existing splicing methods use known splicing sites as training sets, have weak prediction capabilities for hidden non-classical splicing sites, have unsatisfactory accuracy on real data sets (SpliceAI is only 75%), and have slow prediction speed, the present application provides a method and system for predicting pre-mRNA splicing based on a deep learning algorithm.

[0011] The embodiment of the present application is implemented as follows:

[0012] In a first aspect, the present application provides a method for predicting pre-mRNA splicing based on a deep learning algorithm, comprising:

[0013] Based on the public database and RNA-seq sequencing data set, pre-mRNA sequence data is obtained, and the splicing sites and non-splicing sites in the pre-mRNA sequence data are marked;

[0014] Preprocessing the pre-mRNA sequence data to obtain a data set required for generating a model, wherein the data set includes a training set, a validation set, and a test set;

[0015] Building a neural network model that combines a convolutional neural network with a self-attention mechanism, and using the training set and the validation set to train the neural network model to obtain an optimal neural network model;

[0016] The optimal neural network model is used to perform three-category prediction on the test set to obtain the probabilities of donor / acceptor splicing sites and non-splicing sites.

[0017] In a possible implementation, the obtaining of pre-mRNA sequence data based on the public database and the RNA-seq sequencing data set, and the marking of splicing sites and non-splicing sites in the pre-mRNA sequence data, include:

[0018] Obtain the longest transcript information and splicing site location information of the gene from public databases and RNA-seq sequencing data sets;

[0019] A label is created for the data of each bit point in the position information of the splicing site.

[0020] In a possible implementation, the pre-processing of the pre-mRNA sequence data to obtain a data set required for generating a model includes:

[0021] Encoding the pre-mRNA sequence data and the label to obtain an array sequence;

[0022] Performing word encoding on the array sequence to obtain a multi-dimensional feature vector that maps each base;

[0023] Performing position coding on the position of the base in the sequence, wherein the position coding is performed by sine-cosine transform;

[0024] Combining the array sequences after word encoding and position encoding so that the input vector contains both the position feature of the sequence and the content information of the sequence;

[0025] The combined array sequence is divided into a training set, a validation set and a test set according to the proportion.

[0026] In a possible implementation, labeling the data of each point in the position information of the shearing site includes:

[0027] The donor site is coded as [0,1,0], the acceptor site is coded as [1,0,0], and the non-cleavage site is coded as [0,0,0];

[0028] The combined array sequence is matched one by one with the labels and converted into an array format, and the array format matches the format of the neural network model during training and prediction.

[0029] In a possible implementation, the architecture of the convolutional neural network includes a plurality of residual network modules, each of which includes two convolutional layer units;

[0030] Each of the convolutional layer units includes a layer normalization layer, a ReLu activation layer, a convolutional layer and a plurality of convolutional kernels;

[0031] The results of feature extraction by the residual network module in every four layers are added to the features extracted by the first four layers, and different convolution expansion rates are set on different residual network modules.

[0032] In a possible implementation, the neural network model combining a convolutional neural network with a self-attention mechanism includes:

[0033] The convolutional neural network uses equal-length zero-padding operations during each convolution, so that the total length of the sequence remains unchanged and the extracted features are concentrated on the middle neurons;

[0034] Remove the leading and trailing zeros to prepare for feeding into the self-attention mechanism network;

[0035] Based on the features extracted by the convolutional neural network, the features extracted by the convolution are linearly transformed into low-dimensional key, value and query vectors;

[0036] Using the low-dimensional key, value, and query vector, independently calculating a self-attention weight, wherein the self-attention weight is calculated to calculate the correlation between each feature extracted on the same sequence;

[0037] Take the weighted sum of the output values ​​of each attention head.

[0038] In a possible implementation, the using the data set to train the neural network model and saving the optimal weights in the training process includes:

[0039] Training the neural network model according to the training set, and using the validation set to validate the prediction accuracy of the trained neural network model;

[0040] The optimal neural network model is obtained by back-propagation updating, optimizing model parameters and adjusting hyperparameters.

[0041] In a second aspect, the present application provides a system for predicting pre-mRNA splicing based on a deep learning algorithm, comprising:

[0042] A data selection module, used to obtain pre-mRNA sequence data based on the public database and RNA-seq sequencing data set, and mark the splicing sites and non-splicing sites in the pre-mRNA sequence data;

[0043] A data preprocessing module, used to preprocess the pre-mRNA sequence data to obtain a data set required for generating a model, wherein the data set includes a training set, a validation set, and a test set;

[0044] A model building and training module, used to build a neural network model combining a convolutional neural network and a self-attention mechanism, and use the training set and the validation set to train the neural network model to obtain an optimal neural network model;

[0045] A prediction module is used to use the optimal neural network model to perform three-category prediction on the test set to obtain donor / acceptor splicing sites and non-splicing sites.

[0046] In a third aspect, the present application provides a prediction method for identifying gene mutations affecting splicing, using any of the methods in claims 1-7, comprising:

[0047] Centering on the mutant base, extract the 5000bp sequences on both sides to construct a mutant sequence with a total length of 10001bp. Use the same method to construct a wild-type sequence without the mutant base;

[0048] Based on the optimal neural network model, a three-category prediction of shearing is performed for each position of the wild-type and mutant sequences. The difference between the two is calculated as the shearing score of the mutation affecting the site;

[0049] Based on the size of the shearing score at each position, it is estimated whether the variant will affect normal shearing.

[0050] In a fourth aspect, the present application provides an ASO prediction method based on a steric hindrance mechanism, applying the method described in any one of claims 1 to 7, comprising:

[0051] Set the length of the ASO to be designed and the location of the central region to be acted on. This region can be the exon that needs to be inhibited from skipping, or a recessive exon caused by a mutation, or other abnormal splicing conditions such as intron retention and exon truncation;

[0052] Take the active region and the 200 bp before and after it to generate the wild-type sequence;

[0053] Design the length of ASO according to the needs and generate sequences with the same length but all N;

[0054] The sequence consisting of all N is replaced with the wild-type sequence, and the replacement method is to replace the first N positions of the wild-type sequence one by one with a step length of 1 starting from the beginning of the wild-type sequence;

[0055] Based on the optimal neural network model, the wild-type and mutant sequences are predicted to obtain the scores of donor and acceptor in the wild-type sequence, and the scores of donor and acceptor at the same position in the mutant are subtracted to obtain the final ASO prediction score; the score of each predicted ASO is output, as well as its corresponding sequence and position information.

[0056] The technical solution provided by this application can at least achieve the following beneficial effects:

[0057] 1. The present invention adds human lymphocyte, central nervous system tissue and motor neuron RNA-seq data, supplements non-classical splicing site features as input, and the enhanced training set increases the number of splicing sites by 42.9% (150,403 splicing acceptor sites and 140,071 splicing donor sites) compared with the original training set;

[0058] 2. The deep learning neural network model proposed in this invention combines the dilated convolutional neural network (CNN) with the neural network based on self-attention, taking into account the extraction of adjacent features and distal features on the sequence, and improves the accuracy of RNA splicing site recognition and the speed of prediction compared with existing tools;

[0059] 3. The present invention uses the validation set to evaluate the accuracy of the model and obtains the optimal neural network model by repeatedly adjusting the parameters. The trained deep learning neural network can make three-class predictions for each position in the input sequence and evaluate the effect of gene variation on RNA shearing. The prediction time is short, the use is convenient, and the prediction efficiency is improved.

[0060] 4. The deep learning model trained by the present invention can be used for the prediction of ASO sequences, with fast prediction speed and high accuracy, which promotes the development of ASO drugs based on steric hindrance mechanism, improves the efficiency of ASO sequence screening in the early stage of drug development, and reduces the cost of drug development. BRIEF DESCRIPTION OF THE DRAWINGS

[0061] In order to more clearly illustrate the embodiments of the present application or the technical solutions in the prior art, a brief introduction is given below to the drawings required for use in the embodiments or the description of the prior art. Obviously, the drawings described below are some embodiments of the present application. For ordinary technicians in this field, other drawings can be obtained based on these drawings without paying any creative labor.

[0062] Figure 1 It is a schematic diagram of a process for predicting pre-mRNA splicing based on a deep learning algorithm according to an exemplary embodiment of the present application;

[0063] Figure 2 is a schematic diagram of a data acquisition process shown in an exemplary embodiment of the present application;

[0064] Figure 3 is a schematic diagram of a data preprocessing process shown in an exemplary embodiment of the present application;

[0065] Figure 4 is a schematic diagram of a process of making a label according to an exemplary embodiment of the present application;

[0066] Figure 5is a schematic diagram of a flow chart of building a model according to an exemplary embodiment of the present application;

[0067] Figure 6 is a schematic structural diagram of a device for predicting pre-mRNA shearing based on a deep learning algorithm, shown in an exemplary embodiment of the present application;

[0068] Figure 7 is a structural schematic diagram of a computer device shown in an exemplary embodiment of the present application;

[0069] Figure 8 It is a flowchart of a prediction method for identifying gene mutations affecting splicing, shown in an exemplary embodiment of the present application.

[0070] Fig. 9 is a schematic diagram of RNA shearing shown in an exemplary embodiment of the present application;

[0071] Fig.10 is a structural schematic diagram of a neural network model shown in an exemplary embodiment of the present application;

[0072] Fig.11 is a schematic diagram of the input and output file formats of a neural network model software shown in an exemplary embodiment of the present application;

[0073] Fig.12 It is a schematic diagram of a training set and a test set of a neural network model and its performance on the test set shown in an exemplary embodiment of the present application;

[0074] Fig.13 It is a schematic diagram of speed evaluation of a neural network model on a GPU and a CPU shown in an exemplary embodiment of the present application;

[0075] Fig.14 is a schematic diagram showing training and testing of a neural network model for different tissues according to an exemplary embodiment of the present application;

[0076] Fig.15 is a schematic diagram showing the performance of a neural network model in distinguishing between pathogenic and benign variants annotated in the ClinVar database, as shown in an exemplary embodiment of the present application;

[0077] Fig.16 It is a schematic diagram of the prediction of 36 splicing variations in autistic patients by a neural network model shown in an exemplary embodiment of the present application;

[0078] Fig.17 is a schematic structural diagram of a phantom model shown in an exemplary embodiment of the present application;

[0079] Fig.18It is a schematic diagram of a process of training a motif model and analyzing a shearing motif shown in an exemplary embodiment of the present application;

[0080] Fig.19 is a schematic diagram of analyzing a shearing motif of a motif model shown in an exemplary embodiment of the present application;

[0081] Fig. 20 It is a schematic diagram of the process and performance of a neural network model predicting ASO based on steric hindrance mechanism shown in an exemplary embodiment of the present application. DETAILED DESCRIPTION

[0082] In order to make the purpose, implementation mode and advantages of the present application clearer, the exemplary implementation mode of the present application will be clearly and completely described below in conjunction with the drawings in the exemplary embodiments of the present application. Obviously, the described exemplary embodiments are only part of the embodiments of the present application, not all of the embodiments. It should be understood that the specific embodiments described here are only used to explain the present application and are not used to limit the present application.

[0083] It should be noted that the brief description of terms in this application is only for the convenience of understanding the embodiments described below, and is not intended to limit the embodiments of this application. Unless otherwise specified, these terms should be understood according to their ordinary and common meanings.

[0084] The terms "first", "second", "third", etc. in the specification and claims of this application and the above drawings are used to distinguish similar or similar objects or entities, and do not necessarily mean to limit a specific order or sequence, unless otherwise noted. It should be understood that the terms used in this way can be interchangeable under appropriate circumstances.

[0085] The terms "comprises," "comprising," and "having," and any variations thereof, are intended to cover but not exclude inclusion, for example, a product or device comprising a list of components is not necessarily limited to all the components expressly listed but may include other components not expressly listed or inherent to such product or device.

[0086] Before explaining the method for predicting pre-mRNA splicing based on a deep learning algorithm provided in an embodiment of the present application, the application scenario and implementation environment of the embodiment of the present application are first introduced.

[0087] RNA splicing is a complex biological process, which mainly refers to the process of removing introns from pre-mRNA transcribed from a DNA template chain, splicing exons together, and finally generating mRNA. The transcribed product is translated into corresponding proteins to play its specific function. It can be seen that the correct splicing of RNA plays a key role in the regulation of gene expression. The precise splicing of RNA mainly depends on the regular base sequence at the junction of exons and introns, generally the double base GU at the splicing donor (donor, starting from the 5' end of the intron) and the double base AG at the splicing acceptor (acceptor, starting from the 3' end of the intron). In addition, other motifs located on pre-mRNA, such as splicing branch points, polypyrimidine tracts, and splicing silencers and enhancers, will participate in the entire RNA splicing process.

[0088] Fig. 9 It is a schematic diagram of RNA shearing shown in an exemplary embodiment of the present application.

[0089] RNA splicing is a key step in gene transcription regulation. Fig. 9 As shown in the figure, RNA splicing refers to the process by which discontinuous genes encoding proteins remove introns from pre-mRNA in sequence through RNA splicing, and splice exons into mature mRNA. Selective splicing determines the expression of proteins in different tissues and the exercise of their specific functions. It has been found that different abnormal RNA splicing can be caused according to the different locations of gene mutations in the DNA sequence, which can be divided into exon skipping, inclusion of an intron fragment, inclusion of a cryptic exon, and loss of an exon fragment, which further induces dysfunction of the translated protein product and ultimately leads to the occurrence of the disease. For example, the c.892A>G mutation of the gene MAPT reduces RNA splicing, causing exon10 skipping, leading to the occurrence of frontotemporal dementia with Parkinson's disease; the c.2204+6T>C mutation of the gene ELP1 causes exon 20 skipping of chromosome 9, leading to familial autonomic dysfunction, etc.

[0090] At present, the main experimental methods used to study abnormal RNA splicing include transcriptome sequencing, RT-PCR, etc. These traditional wet experiments are mostly limited by the source of samples, time-consuming and costly. Traditional machine learning methods developed based on logistic regression, decision trees and support vector machines, such as GeneSplicer, SpliceMachine and MaxEntScan models, have certain shortcomings. For example, linear regression is very poor in dealing with nonlinear relationships, is not flexible enough in identifying complex patterns, and has low classification accuracy. Decision trees are prone to overfitting and ignore the correlation between attributes in the data set. Support vector machines are not suitable for processing a large number of samples and have no general solution to nonlinear problems. As a result, the features extracted by the model are limited by the existing knowledge of splicing sites, and the prediction ability of non-classical splicing sites is poor, and the classic double bases can only be predicted in a small range of sequences.

[0091] In recent years, with the development of deep learning technology, researchers have used deep learning algorithms to learn the characteristics of gene transcripts in public databases and directly predict splicing sites. For example, SpliceAI invented by Kishore Jaganathan et al. uses a deep convolutional residual network, combines a convolutional structure with a depth of 32 and a residual structure, and trains on sequences of lengths of 80nt, 400nt, 2000nt, and 10,000nt, respectively. It was found that the longer the sequence length used for training, the higher the accuracy of predicting splicing sites, and the accuracy rate increased from 57% for 80nt to 95% for 10,000nt. This is also the first time that predictions have been made within the range of tens of thousands of bases and a relatively ideal accuracy rate has been achieved.

[0092] Based on this, the present application provides a method and system for predicting pre-mRNA splicing based on a deep learning algorithm. By constructing a neural network model that combines convolution and attention mechanisms, the convolution layer is used to extract local features of the sequence, and the self-attention layer is used to learn the relationship between various positions in the global range. It can effectively extract splicing site features, and combine public databases with real data to improve the predictive ability of the model, thereby achieving accurate prediction of RNA splicing sites.

[0093] Next, the technical solution of the present application and how the technical solution of the present application solves the above-mentioned technical problems will be described in detail through embodiments and in combination with the accompanying drawings. The embodiments may be combined with each other, and the same or similar concepts or processes may not be repeated in some embodiments. Obviously, the described embodiments are part of the embodiments of the present application, not all of them.

[0094] Figure 1It is a schematic diagram of a process for predicting a pre-mRNA splicing method based on a deep learning algorithm, shown in an exemplary embodiment of the present application.

[0095] In an exemplary embodiment, Figure 1 As shown, a method for predicting pre-mRNA splicing based on a deep learning algorithm is provided. In this embodiment, the method may include the following steps:

[0096] Step 100: Based on the public database and RNA-seq sequencing data set, pre-mRNA sequence data is obtained, and splicing sites and non-splicing sites in the pre-mRNA sequence data are marked.

[0097] Step 200: Preprocess the pre-mRNA sequence data to obtain a data set required for generating a model, wherein the data set includes a training set, a validation set, and a test set.

[0098] Step 300: Build a neural network model that combines a convolutional neural network with a self-attention mechanism, and use the training set and the validation set to train the neural network model to obtain an optimal neural network model.

[0099] Step 400: Use the optimal neural network model to perform three-category prediction on the test set to obtain the probabilities of donor / acceptor splicing sites and non-splicing sites.

[0100] Figure 2 It is a flowchart of data acquisition shown in an exemplary embodiment of the present application.

[0101] In one possible implementation, Figure 2 As shown, the method of obtaining pre-mRNA sequence data based on the public database and RNA-seq sequencing data set, and simultaneously marking the splicing sites and non-splicing sites in the pre-mRNA sequence data, includes:

[0102] Step 110: Obtain the longest transcript information of the gene and the location information of the splicing site from the public database and RNA-seq sequencing data set.

[0103] This part can be divided into data processing of public databases, blood cell RNA-seq of 30 healthy people and 69 ALS patients, spinal cord and motor cortex RNA-seq of 115 ALS patients, and motor neuron RNA-seq of 240 ALS patients. For the processing of public databases, first download the human genome annotation file (GRCh37 / hg19) and human genome sequence file (hg19.fa) from GENCODE and UCSC official websites, respectively, and screen the transcripts of protein-encoding genes from the downloaded human genome annotation files. If the gene has multiple transcripts, select the longest one, and save each gene transcript and the position of the exon in the transcript.

[0104] For the data in the GTEx database, the downloaded file is the aligned RNA-seq read count file, which includes the location of the exon and the number of times the exon appears in each sample. The exon positions that appear in more than 5 people are retained and do not overlap with the exons already in GENCODE.

[0105] For RNA-seq data, we first use FastQC and Trimmomatic software to quality check the original sequencing files, and then use STAR to align the quality-checked data to obtain the SJ tab file (Splicing Junction file, containing the specific positions of the aligned exons and the number of read counts). Similarly, the exon positions that appear in more than 5 people are retained and do not overlap with the existing exons in GENCODE and GTEx.

[0106] The exon loci obtained from GTEx and RNA-seq data were integrated with the files deposited by GENCODE.

[0107] The bioMart tool on the Ensembl website was used to annotate the paralogous information of genes.

[0108] Step 120: creating a label for each bit point in the position information of the shearing site.

[0109] Figure 3 It is a flowchart of data preprocessing shown in an exemplary embodiment of the present application.

[0110] In one possible implementation, Figure 3 As shown, the pre-mRNA sequence data is preprocessed to obtain the data set required for generating the model, including:

[0111] Step 210: Encode the pre-mRNA sequence data and the label to obtain an array sequence;

[0112] Step 220: performing word encoding on the array sequence to obtain a multi-dimensional feature vector mapped to each base;

[0113] Step 230: performing position coding on the position of the base in the sequence, wherein the position coding is performed by sine-cosine transform;

[0114] Step 240: combining the array sequences after word coding and position coding so that the input vector contains both the position feature of the sequence and the content information of the sequence;

[0115] Among them, it is necessary to character encode the encoded sequence data, that is, map each base into a 64-dimensional feature vector; then the position of the base in the sequence is encoded through sine and cosine transforms, where one base is mapped to a 64-dimensional feature vector, and the corresponding position code is also 64-dimensional; in this way, the input vector contains both the position characteristics of the sequence and the content information of the sequence.

[0116] The calculation formula of position encoding is as follows:

[0117]

[0118] In the formula, pos represents the position of each base, dmodel represents the dimension of the mapping feature vector, which is 64 dimensions, and i represents the dimension index, i∈[0,2,4,...,60,62].

[0119] Step 250: Divide the combined array sequence into a training set, a validation set, and a test set according to proportion.

[0120] The encoded sequences were mapped one by one to the labels, divided into training sets (from chromosomes 2, 4, 6, 8, 10, 11-22, X and Y) and test sets (from chromosomes 1, 3, 5, 7, 9, excluding sequences of paralogous genes), converted into arrays, and saved in npz format.

[0121] Figure 4 It is a schematic diagram of a process of making a label according to an exemplary embodiment of the present application.

[0122] In one possible implementation, Figure 4 As shown, the step of labeling the data of each point in the position information of the shearing site includes:

[0123] Step 121: using [0,1,0] encoding for the donor site, [1,0,0] encoding for the acceptor site, and [0,0,0] encoding for the non-cleavage site;

[0124] Step 122: The combined array sequence is matched one by one with the labels, and converted into an array format, wherein the array format matches the format of the neural network model during training and prediction.

[0125] The specific steps of this part are:

[0126] According to the saved transcript location, the bedtools tool was used to extract the corresponding transcript sequence from the human genome sequence file;

[0127] Label the splicing sites according to the positions of the exons. In a transcript sequence, the splicing donor site is [1,0,0], the splicing acceptor site is [0,1,0], and the non-splicing site is [0,0,1];

[0128] Set a 15000bp window (the target sequence is the middle 5000bp, the flanking sequences are 5000bp on both sides, and the total is 10000bp), 5000bp step length, and extract a 15000bp sequence on the entire transcript and the corresponding tag. To ensure that every position on the transcript is extracted, add a 5000bp sequence with N bases at the beginning of the transcript and a m bp sequence with N bases at the end of the transcript. Length m and total number of sliding steps n step The calculation formula is as follows:

[0129] m = 10000 - len(seq) % 5000

[0130]

[0131] In the formula, len(seq) represents the length of each transcript, and % represents the remainder. Indicates rounding up.

[0132] The extracted transcript sequence contains four bases: N (representing unknown), A, T / U, C and G. These five bases are encoded as 0, 1, 2, 3, and 4 respectively.

[0133] Fig.10 It is a structural diagram of a neural network model shown in an exemplary embodiment of the present application.

[0134] In one possible implementation, Fig.10 As shown, the architecture of the convolutional neural network includes a plurality of residual network modules, each of which includes two convolutional layer units;

[0135] Each of the convolutional layer units includes a layer normalization layer, a ReLu activation layer, a convolutional layer and a plurality of convolutional kernels;

[0136] The results of feature extraction by the residual network module in every four layers are added to the features extracted by the first four layers, and different convolution expansion rates are set on different residual network modules.

[0137] Among them, the convolutional neural network uses a total of 16 residual network modules. Each residual module contains two convolutional layer units. A convolutional layer unit includes a layer normalization layer, a ReLu activation layer and a convolution layer.

[0138] The results of feature extraction by every four residual modules are added to the features extracted by the previous four layers. Each convolution layer uses 64 convolution kernels to extract features. In terms of the window size setting of the convolution kernel, the first 8 residual modules use a convolution kernel with a window size of 11; the 9th to 12th residual modules use a convolution kernel with a window size of 21; and the 13th to 16th residual modules use a convolution kernel with a window size of 41.

[0139] At the same time, in order to increase the receptive field of the convolution kernel and extract features in a wider range, different convolution dilation rates are set on different residual modules. The 1st to 4th residual modules do not perform dilated convolution; the 5th to 8th residual modules set the dilated convolution rate to 4; the 9th to 12th residual modules set the dilated convolution rate to 10; the 13th to 16th residual modules set the dilated convolution rate to 25.

[0140] Figure 5 It is a flowchart of building a model shown in an exemplary embodiment of the present application.

[0141] In one possible implementation, Figure 5 As shown, the neural network model combining the convolutional neural network and the self-attention mechanism includes:

[0142] Step 310: The convolutional neural network uses an equal-length zero-padding operation during each convolution to keep the total length of the sequence unchanged and to concentrate the extracted features on the middle neuron;

[0143] Step 320: Remove the leading and trailing zero-filled parts to prepare for feeding into the self-attention mechanism network;

[0144] In this case, the total length of the sequence is kept unchanged by using equal-length zero-padding operations for each convolution, and the extracted features are concentrated on the middle 5,000 neurons. Finally, the leading and trailing zero-padding parts are removed to prepare for the third part to be sent to the self-attention mechanism network.

[0145] Step 330: Based on the features extracted by the convolutional neural network, linearly transform the features extracted by the convolution into low-dimensional key, value and query vectors;

[0146] Step 340: Using the low-dimensional key, value and query vector, independently calculate the self-attention weight, where the self-attention weight is calculated to calculate the correlation between each feature extracted on the same sequence;

[0147] Step 350: Perform weighted summation on the values ​​output by each attention head.

[0148] Among them, the self-attention mechanism layer in some embodiments of the present application uses a multi-head attention mechanism method and combines local attention and sparse attention mechanisms.

[0149] For the multi-head attention mechanism, first, based on the convolutional neural network obtained in the second part, the features extracted by convolution (each feature is represented by a feature vector of size 64) are linearly transformed into low-dimensional (dmodel / nnhead) key, value and query vectors. Here, 8 attention heads are used, that is, 8 linear transformations are performed; then, the self-attention weights are independently calculated using the low-dimensional key, value and query vectors. The attention weight here refers to the calculation of the correlation between each feature extracted on the same sequence; finally, the values ​​output by each attention head are weighted and summed.

[0150] The multi-head attention mechanism allows different heads to learn different hidden features, thereby achieving better performance. The specific calculation process of a single attention head is as follows:

[0151] ① Multiply the eigenvector with the three matrices WQ, WK, and WV respectively, perform three linear transformations, and derive three new vector matrices: query vector Q, key vector K, and value vector V.

[0152] ② Multiply the query vector Q by the key vector K to obtain the attention weight of the feature vector at each position. At the same time, in order to stabilize the gradient of the network, normalization is used, that is, the obtained attention weight is divided by dk is the dimension of the key vector. After the softmax operation, the correlation between each position on the feature vector is obtained.

[0153] ③The obtained weight is multiplied by the corresponding value vector V, and the weighted value vector is summed to obtain the output of the feature vector at all positions. The calculation formula of multi-head self-attention is as follows:

[0154]

[0155] MultiHead(QK,V)=Concat(head1,…,head h )W O

[0156]

[0157] In order to speed up the training of the model, reduce the memory occupied by the model and improve the reasoning speed in the prediction stage, in the process of matrix multiplication of the multi-head attention mechanism, an attention calculation strategy combining Local attention and Sparse attention is adopted, that is, only the correlation of the 50 features before and after a certain feature is calculated; the correlation is calculated every 50 features in the same sequence.

[0158] After the multi-head sparse attention matrix operation, the features are added to the output of the second part of the convolutional neural network and then sent to the feed forward network layer through layer normalization, which contains two linear transformation layers and one ReLu activation layer.

[0159] In a possible implementation, the neural network model also includes a fully connected layer. After the multi-head sparse attention matrix operation, the features are added to the output of the second part of the convolutional neural network and then sent to the feed forward network layer (Feed forward network) through layer normalization. After the output is obtained, the output is sent to the fully connected layer. The fully connected layer transforms the high-dimensional hidden layer feature vector into three dimensions to perform three-classification of the shearing site.

[0160] In one possible implementation, Figure 5 As shown, the neural network model is trained using the data set and the optimal weights in the training process are saved, including:

[0161] Step 360: Train the neural network model according to the training set, and use the verification set to verify the prediction accuracy of the trained neural network model.

[0162] Step 370: Obtain the optimal neural network model by back-propagation updating, optimizing model parameters and adjusting hyperparameters.

[0163] In the training phase, the output results are used to calculate the loss rate through the cross entropy loss function, and back propagation is performed to continuously update the model parameters until the training loss rate and the loss rate verified by the validation set during the training process are at a low value and remain stable, and the model reaches the optimal solution. Save the optimal model parameters and weights.

[0164] During the training process, the Ranger optimizer is used to update and calculate the network parameters that affect model training and model output; the crossover function is used to calculate the loss rate; and the annealing algorithm is used to adjust the learning rate.

[0165] In this embodiment, the hyperparameters include the number of convolutional neural network layers, the number of convolution kernels, the optimizer, the data partition ratio, the number of training cycles (epochs), and the number of data items in each batch (batch_size).

[0166] It should be understood that, although the various steps in the flowcharts involved in the above-mentioned embodiments are displayed in sequence according to the instructions, these steps are not necessarily executed in sequence according to the order of the instructions. Unless there is a clear explanation in this article, the execution of these steps is not strictly limited in order, and these steps can be executed in other orders. Moreover, at least a portion of the steps in the flowcharts involved in the above-mentioned embodiments may include multiple steps or multiple stages, and these steps or stages are not necessarily executed at the same time, but can be executed at different times, and the execution order of these steps or stages is not necessarily to be carried out in sequence, but can be executed in turn or alternately with other steps or at least a portion of the steps or stages in other steps.

[0167] Corresponding to the aforementioned embodiment of the method for predicting pre-mRNA splicing based on a deep learning algorithm, adopting the same technical concept, the present application also provides an embodiment of a device for predicting pre-mRNA splicing based on a deep learning algorithm.

[0168] Figure 6 It is a schematic diagram of the structure of a device for predicting pre-mRNA based on a deep learning algorithm, shown as an exemplary embodiment of the present application.

[0169] In an exemplary embodiment, Figure 6 As shown, the device for predicting pre-mRNA splicing based on a deep learning algorithm includes:

[0170] Data selection module 1, used to obtain pre-mRNA sequence data based on the public database and RNA-seq sequencing data set, and mark the splicing sites and non-splicing sites in the pre-mRNA sequence data;

[0171] Data preprocessing module 2, used to preprocess the pre-mRNA sequence data to obtain a data set required for generating a model, wherein the data set includes a training set, a validation set, and a test set;

[0172] Model building and training module 3, used to build a neural network model combining a convolutional neural network and a self-attention mechanism, and use the training set and the verification set to train the neural network model to obtain an optimal neural network model;

[0173] Prediction module 4 is used to use the optimal neural network model to perform three-category prediction on the test set to obtain the probabilities of donor / acceptor splicing sites and non-splicing sites.

[0174] For the specific definition of the device for predicting pre-mRNA shearing based on a deep learning algorithm, please refer to the definition of the method for predicting pre-mRNA shearing based on a deep learning algorithm above, which will not be repeated here. Each module in the above-mentioned device for predicting pre-mRNA shearing based on a deep learning algorithm can be implemented in whole or in part by software, hardware and a combination thereof. The above-mentioned modules can be embedded in or independent of the processor in the computer device in the form of hardware, or can be stored in the memory in the computer device in the form of software, so that the processor can call and execute the operations corresponding to the above modules.

[0175] Figure 7 It is a structural schematic diagram of a computer device shown in an exemplary embodiment of the present application.

[0176] In an exemplary embodiment, the above-mentioned method for predicting pre-mRNA splicing based on deep learning algorithm can be applied to Figure 7 In the computer device 10 shown, at this time, the present application can utilize the computer device to construct a neural network model that combines convolution and attention mechanisms. The convolution layer is used to extract local features of the sequence, and the self-attention layer is used to learn the relationship between various positions in the global range. It can effectively extract splicing site features, and combine public databases with real data to improve the predictive ability of the model and achieve accurate prediction of RNA splicing sites.

[0177] In a possible implementation, the structure of the computer device is as follows: Figure 7 As shown, the computer device 10 includes at least a processor 11 , a memory 12 , a communication bus 13 , and a communication interface 14 .

[0178] The processor 11 may be a general-purpose central processing unit (CPU), a network processor (NP), a microprocessor, or may be one or more integrated circuits for implementing the solution of the present application, such as an application-specific integrated circuit (ASIC), a programmable logic device (PLD) or a combination thereof. The above-mentioned PLD may be a complex programmable logic device (CPLD), a field programmable gate array (FPGA), a generic array logic (GAL) or any combination thereof.

[0179] Optionally, the processor 11 may include one or more CPUs. The computer device 10 may include multiple processors 11. Each of these processors 11 may be a single-core processor (single-CPU) or a multi-core processor (multi-CPU).

[0180] It should be noted that the processor 11 here may refer to one or more devices, circuits and / or processing cores for processing data (such as computer program instructions).

[0181] The memory 12 can be a read-only memory (ROM) or other types of static storage devices that can store static information and instructions, or a random access memory (RAM) or other types of dynamic storage devices that can store information and instructions, or an electrically erasable programmable read-only memory (EEPROM), a compact disc (CD-ROM) or other optical disc storage, optical disc storage (including compressed optical disc, laser disc, optical disc, digital versatile disc, Blu-ray disc, etc.), a magnetic disk storage medium or other magnetic storage device, or any other medium that can be used to carry or store the desired program code in the form of instructions or data structures and can be accessed by a computer, but is not limited to this.

[0182] Optionally, the memory 12 may exist independently and be connected to the processor 11 via the communication bus 13 ; the memory 12 may also be integrated with the processor 11 .

[0183] The communication bus 13 is used to transmit information between components (such as between the processor and the memory). The communication bus 12 can be divided into an address bus, a data bus, a control bus, etc. For ease of representation, Figure 7 Only one communication bus is used for illustration, but this does not mean that there is only one bus or one type of bus.

[0184] The communication interface 14 is used for the computer device 10 to communicate with other devices or communication networks. The communication interface 14 includes a wired communication interface or a wireless communication interface. Among them, the wired communication interface can be, for example, an Ethernet interface. The Ethernet interface can be an optical interface, an electrical interface or a combination thereof. The wireless communication interface can be a wireless local area network (WLAN) interface, a cellular network communication interface or a combination thereof.

[0185] In some embodiments, the computer device 10 may also include an output device 15 and an input device 16 ( Figure 7 The output device 15 communicates with the processor 11 and can display information in a variety of ways. For example, the output device 15 can be a liquid crystal display (LCD), a light emitting diode (LED) display device, a cathode ray tube (CRT) display device or a projector. The input device 16 communicates with the processor 11 and can receive user input in a variety of ways. For example, the input device 16 can be a mouse, a keyboard, a touch screen device or a sensor device.

[0186] In some embodiments, the memory 12 is used to store a computer program for executing the solution of the present application, and the processor 11 can execute the computer program stored in the memory 12. For example, the computer device 10 can call and execute the computer program stored in the memory 12 through the processor 11 to implement the steps of predicting the pre-mRNA splicing method based on the deep learning algorithm provided in the embodiment of the present application.

[0187] It should be understood that the method for predicting pre-mRNA splicing based on a deep learning algorithm provided in the present application can be applied to a device for predicting pre-mRNA splicing based on a deep learning algorithm, which can be implemented as part or all of the processor 11 through software, hardware, or a combination of software and hardware, and integrated into the computer device 10.

[0188] Figure 8 It is a flowchart of a prediction method for identifying gene mutations affecting splicing, shown in an exemplary embodiment of the present application.

[0189] In an exemplary embodiment, Figure 8 As shown, a prediction method for identifying gene mutations affecting splicing is also provided, applying any method described in claims 1-7. In this embodiment, the method may include the following steps:

[0190] Step 1000: Taking the mutant base as the center, extract the 5000 bp sequences on both sides to construct a mutant sequence with a total length of 10001 bp. Use the same method to construct a wild-type sequence without the mutant base;

[0191] Step 2000: Based on the optimal neural network model, a three-category prediction of shearing is performed for each position of the wild-type and mutant sequences. The difference between the two is calculated as the shearing score of the mutation affecting the site;

[0192] Step 3000: Estimate whether the variation will affect normal shearing based on the shearing score of each position.

[0193] Fig.11 It is a schematic diagram of the input and output file formats of the neural network model software shown in an exemplary embodiment of the present application.

[0194] In an exemplary embodiment, a computer-readable storage medium is also provided, on which a computer program is stored, wherein the execution of the program can achieve the prediction of RNA splicing sites as described above and the evaluation of the effect of gene variation on RNA splicing. The evaluation of the effect of gene variation on RNA splicing mainly includes the following steps:

[0195] Step 1) The user needs to prepare a vcf file containing mutation information. The input format is shown in Fig.11 ;

[0196] Step 2) The model automatically extracts the sequence before mutation and the sequence containing the mutation site;

[0197] Step 3) The model predicts each position of the sequence before mutation (wild type) and after mutation (mutant type);

[0198] Step 4) The predicted splicing donor score at each position of the mutant sequence is subtracted from the corresponding wild-type splicing donor score to obtain the acceptorΔscore; the donorΔscore of the splicing acceptor is calculated similarly;

[0199] Step 5) For acceptor / donor Δscore greater than 0 (indicating the emergence of a new splicing site), sort from large to small, and for Δscore less than 0 (indicating the disappearance of the original splicing site), sort from small to large; according to the user's choice, output the largest Δscore (minimum if negative) and its corresponding position within the range of n bp before and after the mutation site. The maximum value of N is 5000bp.

[0200] Output format see Fig.11 The Output section starts with Fig.11 The output result is explained by taking the first mutation in as an example, and the input mutation (chr16:31202282A>G) is predicted. This mutation may cause a new acceptor site to appear at a position of +36bp away from the mutation site (Δscore=0.62), and the original +2bp acceptor site to disappear (Δscore=0.99). At the same time, it will also cause the -173bp donor site to disappear (Δscore=0.06).

[0201] The application scenarios of some embodiments of the present application are evaluated and tested below.

[0202] Fig.12 It is a schematic diagram of the training set and test set of a neural network model and its performance on the test set shown in an exemplary embodiment of the present application.

[0203] like Fig.12 As shown in a, the training set contains 20,270 sequences and 673,670 splicing sites. The test set contains 7,748 sequences and 67,303 splicing sites.

[0204] Fig.12 bf shows the performance of the model on GENCODE and RNA-seq datasets, and compares it with several existing models.

[0205] Among them, SpliceAI and Pangolin use the deep learning method of convolutional neural network. MaxEntScan uses the machine learning algorithm of maximum entropy model. It can be seen that Spliformer's overall performance in the combined dataset of GENCODE and RNA-seq of different tissues (n = 67303 splicing sites) is slightly better than other models, with an AUPRC of 0.72. In the subgroups, Spliformer predicted GENCODE splicing sites (n = 27061 splicing sites) with an AUPRC of 0.93, predicted lymphocyte splicing sites (n = 44684 splicing sites) with an AUPRC of 0.69, predicted central nervous system tissue splicing sites (n = 57686 splicing sites) with an AUPRC of 0.71, and predicted motor neuron splicing sites (n = 28534 splicing sites) with an AUPRC of 0.62.

[0206] Some embodiments of the present application propose speed evaluation of the model on GPU and CPU.

[0207] Fig.13 FIG. 1 is a schematic diagram of speed evaluation of a neural network model on a GPU and a CPU, shown in an exemplary embodiment of the present application. Fig.13 As shown, the splicing prediction speed of Spliformer, SpliceAI and Pangolin for the same number of SNVs on CPU and GPU is demonstrated.

[0208] Fig.13 a is the speed of predicting the effect of 11 million SNVs on RNA splicing on an A100-PCIE GPU (40GB) graphics card. The prediction speed of Spliformer is 6.5 days, SpliceAI is 31 days, and Pangolin is 9 days; Fig.13b is the speed of predicting the effect of 10,000 SNVs on RNA shearing on an A100-PCIE GPU (40GB) graphics card. The prediction speed of Spliformer is 9 minutes, SpliceAI is 78 minutes, and Pangolin is 18 minutes; Fig.13 c is the speed of predicting the effect of 10,000 SNVs on RNA splicing on an Intel(R) Xeon(R) Gold5218 CPU; the prediction speed of Spliformer is 142 minutes, SpliceAI is 201 minutes, and Pangolin is 177 minutes.

[0209] Fig.14 It is a schematic diagram showing the training and testing of a neural network model for different tissues according to an exemplary embodiment of the present application.

[0210] like Fig.14 As shown, Spliformer was trained and tested using RNA-seq data from different human tissues (frontal cortex, substantia nigra, skeletal muscle, cerebellum, hypothalamus and spinal cord) to evaluate the accuracy of different models on different tissues.

[0211] Fig.14 a is a flowchart for training and testing models using RNA-seq data from GTEx. The aligned RNA-seq read count files of six tissues (frontal cortex, substantia nigra, skeletal muscle, cerebellum, hypothalamus, and spinal cord) are from GTEx. The frontal cortex includes 31 samples, the substantia nigra includes 17 samples, the skeletal muscle includes 121 samples, the cerebellum includes 18 samples, the hypothalamus includes 21 samples, and the spinal cord includes 15 samples.

[0212] For different tissues, if the exon fragments that were aligned appeared in 5 or more people, the position of the exon was retained as the splicing site. The steps of extracting sequences and making labels, the framework of the model is the same as in Example 1, and will not be described in detail here. Sequences from chromosomes 2, 4, 6, 8, 10-22, X, and Y were used for model training, and sequences from chromosomes 1, 3, 5, 7, and 9 (sequences without any paralogous homologs) were used for model testing.

[0213] Fig.14 b shows the performance of the model based on GTEx multi-tissue RNA-seq data. Spliformer tissue-specific splicing prediction models were constructed in frontal cortex, substantia nigra, skeletal muscle, cerebellum, hypothalamus and spinal cord tissues. The radar chart shows the prediction performance of Spliformer, Pangolin and SpliceAI in different tissues.

[0214] The models were comprehensively evaluated using AUPRC, Top0.5, Top1, precision, and recall. The AUPRC of Spliformer, SpliceAI, and Pangolin were similar, but Spliformer was more accurate in different tissues (frontal cortex: 0.89, substantia nigra: 0.83, skeletal muscle: 0.9, cerebellum: 0.9, hypothalamus: 0.89, spinal cord: 0.9).

[0215] Fig.15 It is a schematic diagram of the performance of a neural network model in distinguishing between pathogenic and benign variants marked in the ClinVar database, as shown in an exemplary embodiment of the present application.

[0216] The models proposed in some embodiments of the present application are shown to be effective in distinguishing between variants annotated as pathogenic and benign in the ClinVar database.

[0217] like Fig.15 As shown, the model is used to predict pathogenic and benign variants in the ClinVar database and analyze the relationship between splicing variation and pathogenicity.

[0218] The ClinVar dataset was processed using the following criteria.

[0219] 1) The variant is a SNV or INDELS (the REF or ALT field is a single base in the case of insertion / deletion; the insertion and deletion bases are less than 5000bp).

[0220] 2) The variant is contained in the longest transcript of the gene.

[0221] 3) Each variant can only exist in one gene.

[0222] 4) Mutations included those near exons (≤50 bp from exon boundaries), deep intronic regions (≥50 bp from exon boundaries), and synonymous mutations.

[0223] 5) Only benign and pathogenic variants that can be scored by Spliformer, SpliceAI, and Pangolin were considered.

[0224] 6) In order to observe the effect of the variant on splicing at the longest distance, the parameter -D4999 (maximum distance between the variant and the gain / loss splicing site) was used when making predictions using the three software.

[0225] 7) The larger value of gainΔscore and lossΔscore in the shearing donor or acceptor is selected as the final result.

[0226] Fig.15a shows the proportion of pathogenic variants near exons (≤50bp from exon boundaries, n=9184) and deep introns (≥50bp from exon boundaries, n=195) predicted to be splicing variants at different thresholds. Pathogenic splicing variants are significantly enriched in regions near exons (Fisher's exact test, P<2.2E-16).

[0227] Fig.15 b shows the ratio of splicing variants to non-splicing variants in high pLI (pLI ≥ 0.9, the probability of being loss-of-function intolerant) genes at different thresholds. Splicing variants are significantly enriched in high pLI genes (Fisher’s exact test, P < 2.2E-16).

[0228] Fig.15 c shows the proportion of pathogenic / benign variants predicted as splicing variants at different thresholds (intronic variants and synonymous mutations located in non-splicing regions, including 100,069 pathogenic variants and 165,479 benign variants). Splicing variants are significantly enriched in pathogenic variants (Fisher's exact test, P < 2.2E-16).

[0229] Fig.16 It is a schematic diagram of the prediction of 36 splicing variations in autistic patients by a neural network model shown in an exemplary embodiment of the present application.

[0230] like Fig.16 Shown is the prediction of 36 splicing variants in autism patients by Spliformer.

[0231] Fig.16 a The bar graph shows the number and ratio of successful and failed splicing mutations predicted by Spliformer, SpliceAI and Pangolin. The validation rate of Spliformer is 0.79, the validation rate of SpliceAI is 0.75, and the validation rate of Pangolin is 0.71.

[0232] Fig.16 b shows two examples of Spliformer successfully predicting splicing variations.

[0233] The RNA-seq samples of the corresponding samples were visualized using ggsashimi. The visualization results verified the predictions of Spliformer, and the gene mutation (FH:c.378+11T>C) led to the generation of a new splicing donor on FHexon3 (Spliformer predicted donor gainΔscore was 0.98); the gene mutation (ELOVL1:c.619-3T>G) led to the generation of a new splicing acceptor site on ELOVL1 exon8 (Spliformer predicted acceptor gainΔscore was 0.40).

[0234] Fig.17 is a schematic structural diagram of a phantom model shown in an exemplary embodiment of the present application, Fig.18 is a schematic diagram of a process of training a motif model and analyzing a shearing motif, shown in an exemplary embodiment of the present application. Fig.19 It is a schematic diagram of the analysis of the shear regulation motif of the motif model shown in an exemplary embodiment of the present application.

[0235] In an exemplary embodiment, the present application also proposes a motif model based on the cutting method proposed in some embodiments of the present application, which can be based on the visualization of the self-attention layer and the analysis of the motif.

[0236] like Fig.17 The model framework of the Spliformer-motif shown in the figure, as well as the attention layer that extracts the predictions of the Spliformer-motif model, and the analysis and visualization results of the cutting motif.

[0237] The Spliformer-motif model takes a pre-mRNA sequence of 81 nt as input. After the sequence is hot-encoded, three convolutional layers are used to extract local sequence features. Each layer uses 64 convolution kernels with window sizes of 3, 5, and 7. Then, six self-attention layers are used to calculate the weight of the interrelationship of the features. The architecture of the self-attention layer is the same as that of the vanilla transformer encoder layer. Finally, the probability of any position on the sequence being a splice acceptor, splice donor, and non-splicing site is predicted. The self-attention weight score (AWS) of the motif is calculated as follows: After matrix multiplication of matrices Q and K in the ScaledDotProductAttention function layer, the AWS of the motif of the six self-attention layers is directly extracted. Matrices Q and K are low-dimensional matrices generated by linear projection of the sequence. Then, the AWS of the six attention layers are summed according to the position of each motif in the input sequence. Finally, the summed AWS matrix (81*81) is normalized.

[0238] like Fig.18 As shown, the process of training the Spliformer-motif model and analyzing the shearing motif.

[0239] We trained the model using sequences from chromosomes 2, 4, 6, 8, 10-22, X, and Y from GENCODE, GTEx, and 30 healthy human lymphocyte RNA-seq data. We extracted the self-attention weights when the model predicted sequences from chromosomes 1, 3, 5, 7, and 9 (excluding paralogous genes), visualized the self-attention weight values, and analyzed the shear motifs.

[0240] Fig.18 b is based on the analysis of attention weights of 21,480 sequences containing splice donors.

[0241] The bar graph shows the top 10 target donor motifs, all of which contain the splicing donor canonical sequence N / GT in the center. The pie chart shows that the target donor motif N / GTRAG accounts for 43.15% of all motifs, N / GT (excluding GTRAG) accounts for 55.12%, and N / GC accounts for 1.58%.

[0242] Fig.18 c is the visualization of the proportion of unknown target donor motifs (type of others in the pie chart) at each position using the R package (ggseqlog). The results show that the unknown target donor motif is characterized by N / AUAUCC.

[0243] Fig.18 d is an example of target donor motif heat map visualization.

[0244] From the heat map, we can see that the target donor motif AUUACAGGUGUGA located on the gene DNAJC13 (chr3:132140515-132140595, hg19) has a higher self-attention weight value than other motifs (left column).

[0245] Fig.18 e Based on analysis of attention weights for 22,037 splice receptor-containing sequences.

[0246] The bar graph shows the top 10 target acceptor motifs. Seven target acceptor motifs contain the classic splicing receptor sequences CAG / N, TAG / N, and AAG / N. The pie chart shows that the target acceptor motif CAG / N accounts for 58.10% of all motifs, TAG / N accounts for 30.47%, AAG / N accounts for 8.77%, and GAG / N accounts for 2.53%.

[0247] Fig.18 f is the visualization of the proportion of unknown target acceptor motifs (type of others in the pie chart) appearing at each position using the R package (ggseqlog). The results show that the unknown target acceptor motif is characterized by U / AC.

[0248] Fig.18 g is an example of a heat map visualization of the target acceptor motif. From the heat map, we can see that the target acceptor motif UUUGAGAUGGAGU located on the gene CCDC30 (chr1:42980300-42980380, hg19) has a higher self-attention weight value than other motifs (left column).

[0249] Fig.19 a is the visualization of the proportion of donor regulatory motif sequences appearing at each position using the R package (ggseqlog). The analysis included a total of 81,940 donor regulatory motifs (AWS>0.5) and the results showed that donor regulatory motifs were characterized by N / AG.

[0250] Fig.19 b is a heat map visualization of donor regulatory motifs. In this example, the donor-related motif UUCCAGGCCCACG (the darkest motif in the heat map row) located on the gene PLEKHN1 (chr1:905935-906014, hg19) is most correlated with sequences containing donor sites.

[0251] Fig.19 c Bar graph shows the composition of acceptor regulatory motifs (AWS>0.5), a total of 51,037 acceptor regulatory motifs (AWS>0.5). 44.43% of the splicing regulatory motifs do not contain known splicing receptor regulatory sequences.

[0252] Fig.19d Heatmap visualization of acceptor regulatory motifs containing splice branch point sequences. In this example, the acceptor regulatory motif UACAAGGCUCACC (the darkest motif in the heatmap) located on the gene HEMK1 (chr3:506088316-50608395, hg19) contains the canonical splice branch point sequence yUNAy.

[0253] Fig.19 e is the visualization of the proportion of sequences that do not contain known acceptor regulatory motifs at each position using the R package (ggseqlog). The results show that the acceptor regulatory motif is characterized by N / AG.

[0254] Fig.19 f Heatmap visualization of the acceptor regulatory motifs that do not contain the splice branch point sequence. In this example, the acceptor regulatory motif CCCCAGUGUACCC (the darkest motif in the heatmap row) located in the gene OGG1 (chr3:9798699-9798791, hg19) contains the N / AG sequence feature.

[0255] Fig.19 g is the frequency distribution diagram of acceptor regulatory motifs. The acceptor-related motifs contain the classic splice branch site sequence yUNAy. Such motifs are mainly distributed about 30bp upstream of the splice acceptor site. The dotted line is the intron-exon junction.

[0256] Fig.19 h is the frequency distribution diagram of donor regulatory motif, which is mainly distributed around 40bp upstream of the splicing donor site. The dotted line is the intron-exon junction.

[0257] Fig. 20 To use Spliformer to predict the process and performance of ASO based on steric hindrance mechanism.

[0258] Based on the three currently known ASOs for the treatment of spinal muscular atrophy (SMA), ataxia telangiectasia (AT), and Duchenne muscular dystrophy (DMD), we evaluated whether Spliformer has the ability to predict ASO drugs and their corresponding performance.

[0259] The process of evaluating the model is as follows Fig. 20As shown in a, the specific process is as follows: (1) Determine the sequence information of the known ASO and the exon it acts on; (2) Extract the exon and the 200bp sequence on both sides of the exon, which is the wild-type sequence; (3) Generate a sequence of the corresponding length N according to the length of the ASO sequence. The sequence with N is replaced with the extracted wild-type sequence from beginning to end in steps of 1, and the replaced sequence is the mutant type; (4) All sequences are input into the model prediction to obtain the prediction scores of the exons donor and acceptor; (5) Subtract the prediction scores of the wild-type sequence donor and acceptor from the prediction scores of the mutant sequence donor and acceptor, and convert the scores using the following formula to obtain the final ASO prediction score:

[0260] score=max(l(WT a )-l(VT a )),l(WT d )-l(VT d ))

[0261] l(n)=log((n+10 -10 ) / (1-n+10 -10 )

[0262] Where WT a Represents the prediction score of the wild-type sequence acceptor; WT d Represents the prediction score of the wild-type sequence donor; VT a Represents the prediction score of the mutant sequence acceptor; VT d Represents the prediction score of the mutant sequence donor.

[0263] The ASO information used in the evaluation is as follows: Nusinersen is used to treat SMA, the sequence is: 5'-UCACUUUCAUAAUGCUGG-3', the sequence position is chr5:69372411-69372428 (hg19), and its function is to prevent SMN2exon7 skipping. AT043 and AT046 are used to treat AT. The sequences are 5'-GCUGAAUAAGAUCCUGAAGA-3' and 5'-AAAAACGGGUAAAGACAUGC-3'. The sequence positions are chr11:108179773-108179792 (hg19) and chr11:108179825-108179844 (hg19), respectively. Their function is to induce skipping of the pseudoexon (chr11:108179696-18179832, hg19) produced by ATM c.5763-1050A>G. Casimersen is used to treat DMD. The sequence is 5'-CAAUGCCAUCCUGGAGUUCCUG-3' and the sequence position is chrX:31986613-31986634 (hg19). Its function is to promote DMD exon 45 skipping.

[0264] Fig. 20 b shows the results of using Spliformer to predict ASOs that inhibit SMN2 exon7skipping. Among 437 ASOs, the ASO predicted to be Nusinersen ranks 16th, which is in the top 3.66% of the total ranking.

[0265] Fig. 20 c is the ASO result predicted to cause skipping of the pseudoexon generated by ATMc.5763-1050A>G. Among the 518 ASOs, the ASO predicted to be AT043 ranked 47th and the ASO predicted to be AT046 ranked 13th, which were the top 9.07% and top 2.51% of the total ranking respectively.

[0266] Fig. 20 d is the predicted ASO result that caused DMDexon 45skipping. Among the 555 ASOs, the ASO predicted to be Casimersen ranked 7, which is in the top 1.26% of the total ranking.

[0267] These results show that Spliformer has a high accuracy in predicting ASOs based on steric hindrance mechanisms. This high-efficiency and low-cost method has important application value in the initial screening stage of ASO drugs.

[0268] The technical features of the above embodiments may be combined arbitrarily. To make the description concise, not all possible combinations of the technical features in the above embodiments are described. However, as long as there is no contradiction in the combination of these technical features, they should be considered to be within the scope of this specification.

[0269] The above-described embodiments only express several implementation methods of the present application, and the descriptions thereof are relatively specific and detailed, but they cannot be understood as limiting the scope of the invention patent. It should be pointed out that, for a person of ordinary skill in the art, several variations and improvements can be made without departing from the concept of the present application, and these all belong to the protection scope of the present application. Therefore, the protection scope of the patent of the present application shall be subject to the attached claims.

Claims

1. A method for predicting pre-mRNA splicing based on a deep learning algorithm, characterized in that: include: Based on the public database and RNA-seq sequencing data set, pre-mRNA sequence data is obtained, and the splicing sites and non-splicing sites in the pre-mRNA sequence data are marked; Preprocessing the pre-mRNA sequence data to obtain a data set required for generating a model, wherein the data set includes a training set, a validation set, and a test set; Building a neural network model that combines a convolutional neural network with a self-attention mechanism, and using the training set and the validation set to train the neural network model to obtain an optimal neural network model; The optimal neural network model is used to perform three-category prediction on the test set to obtain the probabilities of donor / acceptor splicing sites and non-splicing sites.

2. The method for predicting pre-mRNA shearing based on a deep learning algorithm according to claim 1, characterized in that The method of obtaining pre-mRNA sequence data based on the public database and RNA-seq sequencing data set, and simultaneously marking the splicing sites and non-splicing sites in the pre-mRNA sequence data, includes: Obtain the longest transcript information and splicing site location information of the gene from public databases and RNA-seq sequencing data sets; A label is created for the data of each bit point in the position information of the splicing site.

3. The method for predicting pre-mRNA shearing based on a deep learning algorithm as claimed in claim 2, characterized in that The pre-processing of the pre-mRNA sequence data to obtain a data set required for generating a model includes: Encoding the pre-mRNA sequence data and the label to obtain an array sequence; Performing word encoding on the array sequence to obtain a multi-dimensional feature vector that maps each base; Performing position coding on the position of the base in the sequence, wherein the position coding is performed by sine-cosine transform; Combining the array sequences after word encoding and position encoding so that the input vector contains both the position feature of the sequence and the content information of the sequence; The combined array sequence is divided into a training set, a validation set and a test set according to the proportion.

4. The method for predicting pre-mRNA shearing based on a deep learning algorithm as claimed in claim 3, characterized in that The step of labeling the data of each point in the position information of the shearing site comprises: The donor site is coded as [0,1,0], the acceptor site is coded as [1,0,0], and the non-cleavage site is coded as [0,0,0]; The combined array sequence is matched one by one with the labels and converted into an array format, and the array format matches the format of the neural network model during training and prediction.

5. The method for predicting pre-mRNA shearing based on a deep learning algorithm according to claim 1, characterized in that: The architecture of the convolutional neural network includes a plurality of residual network modules, each of which includes two convolutional layer units; Each of the convolutional layer units includes a layer normalization layer, a ReLu activation layer, a convolutional layer and a plurality of convolutional kernels; The results of feature extraction by the residual network module in every four layers are added to the features extracted by the first four layers, and different convolution expansion rates are set on different residual network modules.

6. The method for predicting pre-mRNA shearing based on a deep learning algorithm according to claim 1, characterized in that: The neural network model combining a convolutional neural network and a self-attention mechanism includes: The convolutional neural network uses equal-length zero-padding operations during each convolution, so that the total length of the sequence remains unchanged and the extracted features are concentrated on the middle neurons; Remove the leading and trailing zeros to prepare for feeding into the self-attention mechanism network; Based on the features extracted by the convolutional neural network, the features extracted by the convolution are linearly transformed into low-dimensional key, value and query vectors; Using the low-dimensional key, value, and query vector, independently calculating a self-attention weight, wherein the self-attention weight is calculated to calculate the correlation between each feature extracted on the same sequence; Take the weighted sum of the output values ​​of each attention head.

7. The method for predicting pre-mRNA shearing based on a deep learning algorithm according to claim 1, characterized in that: The using the data set to train the neural network model and saving the optimal weights during the training process includes: Training the neural network model according to the training set, and using the validation set to validate the prediction accuracy of the trained neural network model; The optimal neural network model is obtained by back-propagation updating, optimizing model parameters and adjusting hyperparameters.

8. A system for predicting pre-mRNA splicing based on a deep learning algorithm, characterized in that: include: A data selection module, used to obtain pre-mRNA sequence data based on the public database and RNA-seq sequencing data set, and mark the splicing sites and non-splicing sites in the pre-mRNA sequence data; A data preprocessing module, used to preprocess the pre-mRNA sequence data to obtain a data set required for generating a model, wherein the data set includes a training set, a validation set, and a test set; A model building and training module, used to build a neural network model combining a convolutional neural network and a self-attention mechanism, and use the training set and the validation set to train the neural network model to obtain an optimal neural network model; A prediction module is used to use the optimal neural network model to perform three-category prediction on the test set to obtain donor / acceptor splicing sites and non-splicing sites.

9. A prediction method for identifying gene mutations affecting splicing, using any of the methods of claims 1-7, characterized in that: include: Centering on the variant base, 5000 bp sequences on both sides were extracted to construct a mutant sequence with a total length of 10001 bp. The wild-type sequence without the mutated base was constructed in the same way; Based on the optimal neural network model, a three-category prediction of shearing is performed for each position of the wild-type and mutant sequences. The difference between the two is calculated as the shearing score of the mutation affecting the site; Based on the size of the shearing score at each position, it is estimated whether the variant will affect normal shearing.

10. A method for screening ASO based on steric hindrance mechanism, using the method described in any one of claims 1 to 7, characterized in that: include: Set the length of the ASO to be designed and the location of the central region to be acted on. This region can be the exon that needs to be inhibited from skipping, or a recessive exon caused by a mutation, or other abnormal splicing conditions such as intron retention and exon truncation; Take the active region and the 200 bp before and after it to generate the wild-type sequence; Design the length of ASO according to the needs and generate sequences with the same length but all N; The sequence consisting of all N is replaced with the wild-type sequence, and the replacement method is to replace the first N positions of the wild-type sequence one by one with a step length of 1 starting from the beginning of the wild-type sequence; Based on the optimal neural network model, the wild-type and mutant sequences are predicted to obtain the scores of donor and acceptor in the wild-type sequence, and the scores of donor and acceptor at the same position in the mutant are subtracted to obtain the final ASO prediction score; the score of each predicted ASO is output, as well as its corresponding sequence and position information.

Images