Base resolution prediction method based on full convolutional neural network

Abstract

The invention relates to a base resolution prediction method based on a full convolutional neural network, and the method comprises the following steps: obtaining a to-be-predicted DNA sequence, inputting the to-be-predicted DNA sequence into a pre-established and trained full convolutional neural network base resolution primitive mining model, and carrying out the prediction of a base resolution signal; the full convolutional neural network base resolution primitive mining model comprises an encoder structure used for down-sampling input data and extracting specific features of a sequence; the decoder structure is used for carrying out up-sampling on the features and carrying out modeling on the signal of each basic group; the skip structure is used for combining the features learned in the encoder structure with the discrimination information in the decoder structure; in the training process of the full convolutional neural network base resolution primitive mining model, the base resolution is combined with a signal value to serve as a label for supervised learning. Compared with the prior art, the method realizes the prediction and analysis of the specific binding loci of the primitives on the base resolution level.

Description

Technical Field

[0001] The present invention relates to the technical field of transcription factor binding site recognition, and in particular to a base resolution prediction method based on a fully convolutional neural network. Background Art

[0002] The identification of transcription factor binding sites is the basis for studying gene regulatory networks in biological systems, because transcription factors can activate or inhibit gene transcription by binding to specific DNA non-coding regions, thus playing an indispensable role in gene expression. In recent years, in order to study the specificity of transcription factor and DNA binding sites and the prediction of motifs, a lot of computational work has been done and many algorithms, computational tools and databases have been developed. However, a deep understanding of the binding mechanism of transcription factors and DNA is still lacking, especially the need for comprehensive and accurate computational methods for the whole genome to systematically reveal the binding mechanism of TFBSs is becoming more and more urgent.

[0003] The rapid development of high-throughput sequencing technology has brought a large amount of data related to transcription factors and DNA binding, which can be used to study TFBSs-related tasks. However, these methods often face defects such as poor generalization ability, high time complexity, and sensitivity to data noise, which makes them often unable to cope with the mining and analysis of massive data.

[0004] Therefore, there is an urgent need to improve such computational methods to meet the requirements of accurately identifying cis-regulatory motifs on massive data. In recent years, deep learning has achieved amazing performance in many fields, such as computer vision and natural language processing. At the same time, with the development of high-throughput sequencing technology and deep learning, researchers have been inspired to develop deep learning-based methods to predict transcription factor binding sites. These deep learning-based methods have been developed and used to systematically study transcription factor binding site recognition, and many of these algorithms have achieved very good performance. However, these deep learning-based methods mainly focus on identifying whether a DNA sequence contains a motif sequence, and cannot accurately locate the motif position.

[0005] To remedy such problems, DESSO uses a CNN model to mine motifs from a given ChIP-seq peak and uses a statistical model based on binomial distribution to optimize the recognition ability of motifs. FCNA predicts motifs on ChIP-seq data by using a fully convolutional neural network (FCN) and global average pooling (GAV). D-AEDNet uses an encoder and decoder architecture to identify the location of transcription factor DNA binding sites in DNA sequences. BPNet introduces dilated convolution to predict the base-resolution transcription factor DNA binding probability of TFs. Although these deep learning-based methods have achieved good accuracy in the tasks of classification, positioning or predicting motifs, there is still a lack of a high-precision integration method for multiple tasks of transcription factor binding sites. Summary of the invention

[0006] The purpose of the present invention is to overcome the defects of the above-mentioned prior art and provide a more accurate base resolution prediction method based on a fully convolutional neural network.

[0007] The purpose of the present invention can be achieved by the following technical solutions:

[0008] A base resolution prediction method based on a fully convolutional neural network comprises the following steps:

[0009] Obtain the DNA sequence to be predicted, input it into the pre-established and trained fully convolutional neural network base resolution primitive mining model, and predict the base resolution signal;

[0010] The fully convolutional neural network base resolution primitive mining model includes:

[0011] The encoder structure is used to downsample the input data and extract specific features of the sequence;

[0012] The decoder structure is used to upsample the features and model the signal of each base;

[0013] The skip structure is used to combine the features learned in the encoder structure with the discriminative information in the decoder structure;

[0014] During the training process of the fully convolutional neural network base resolution primitive mining model, base resolution combined with signal value is used as a label for supervised learning.

[0015] Furthermore, the encoder structure includes a convolutional module, a GRU layer and a global average aggregation layer connected in sequence, the number of the convolutional modules is multiple, and each convolutional module includes a convolutional layer, a ReLU layer, a maximum aggregation layer and a dropout layer connected in sequence, which are used to gradually reduce the dimension of the features and encode the specific features of the sequence;

[0016] The GRU layer includes a bidirectional GRU layer and a dropout layer connected to each other, and the bidirectional GRU layer is used to learn long-term dependencies in sequence-specific features;

[0017] The global average aggregation layer is used to capture the global background of DNA sequences.

[0018] Further, the decoder structure includes a plurality of groups of upsampling and mixing block combinations connected in sequence, each group of upsampling and mixing block combinations includes an upsampling layer and a mixing block connected in sequence, and the mixing block includes a batch normalization layer, a ReLU layer and a convolution layer connected in sequence;

[0019] The upsampling layer is used to restore the size of the downsampled features;

[0020] The mixing block is used to rescale the values ​​of the upsampled features.

[0021] Further, the number of convolution modules in the encoder structure is three, and the number of upsampling and mixing block combinations in the decoder structure is four;

[0022] The skip structure includes: the output of the first convolution module is connected to the input of the fourth group of upsampling and mixing block combination; the output of the second convolution module is connected to the input of the third group of upsampling and mixing block combination; the output of the GRU layer is connected to the input of the second group of upsampling and mixing block combination.

[0023] Furthermore, the DNA sequence input into the fully convolutional neural network base-resolution primitive mining model is converted into a one-hot encoded matrix.

[0024] Furthermore, in the training process of the fully convolutional neural network base-resolution primitive mining model, the base-resolution MSE loss is used as the loss function, and the expression of the MSE loss is:

[0025]

[0026] Where MSE is the value of the MSE loss function, N is the number of samples in the training data, L is the length of each input sequence, and z i,j and are the predicted signal and the true signal respectively, α is the regularization parameter, and ||*||2 represents L2 regularization.

[0027] Furthermore, the fully convolutional neural network base resolution primitive mining model scales the input base resolution combined signal value through the log10 function, and the corresponding scaling expression is:

[0028] signal=log 10 (1+signals)

[0029] Where signal is the base resolution combined signal value after scaling, and signals is the input base resolution combined signal value.

[0030] Furthermore, the Adam optimization algorithm is used to perform parameter optimization during the training process of the fully convolutional neural network base-resolution primitive mining model.

[0031] Furthermore, the method also includes distinguishing the DNA sequence to be predicted into a binding sequence and a non-binding sequence according to the maximum value of the predicted base resolution signal.

[0032] Furthermore, the method also includes locating a binding region with a length of 100bp according to the position of the maximum value of the predicted base resolution signal, selecting the sub-region with the highest score according to the training weight of the first convolutional layer in the encoder structure, arranging the selected sub-regions, calculating the corresponding position frequency matrix, and finally matching the position frequency matrix with the primitives in the standard database to obtain the prediction result of the primitives.

[0033] Compared with the prior art, the present invention has the following advantages:

[0034] (1) The present invention proposes a fully convolutional neural network base-resolution primitive mining model, which introduces binding signal data in addition to the input sequence data, so as to solve the problem of high-precision prediction of primitive-specific binding sites. Currently, common primitive mining based on deep learning mainly focuses on predicting whether a DNA sequence is a primitive, but such methods cannot predict primitives at base resolution. FCNsignal developed by the present invention uses a fully convolutional neural network combined with a gated recursive unit to model sequence data and binding signal data, thereby realizing the prediction and analysis of primitive-specific binding sites at the base resolution level.

[0035] (2) The experimental results on 53 transcription factor binding datasets and 6 chromatin accessibility datasets show that the model proposed in this embodiment is superior to the currently common deep learning methods. In addition, the above experiments also found that the model has the potential to identify disease-related single nucleotide polymorphisms. BRIEF DESCRIPTION OF THE DRAWINGS

[0036] Figure 1A schematic diagram of a flow chart of a base resolution prediction method based on a fully convolutional neural network provided in an embodiment of the present invention;

[0037] Figure 2 A schematic diagram of the structure of a fully convolutional neural network base-resolution primitive mining model provided in an embodiment of the present invention;

[0038] Figure 3 A schematic diagram of the Pearson correlation coefficient of four TFs including CTCF, YY1, FOSL1 and USF1 provided in an embodiment of the present invention;

[0039] Figure 4 Schematic diagram of a model prediction result provided in an embodiment of the present invention, wherein (A) is that FCNsignal is mainly composed of an encoder structure, a decoder structure and a skip structure, which takes a DNA sequence as input and predicts a base resolution signal; (B) is that FCNsignal can accurately capture the signals of the shifted binding regions located at 250bp, 500bp and 750bp respectively; (C) is that FCNsignal can accurately capture the signals of the binding regions of randomly inserted negative examples;

[0040] Figure 5 The performance comparison of an FCNsignal provided in an embodiment of the present invention and a comparative method on 53 ChIP-seq data sets; wherein, (A) is the MSE and Pearsonr values ​​of BPNet, FCNA and FCNsignal, and the Wilson test p-value (paired) of these two indicators between FCNsignal and the other two methods; (B) is the Auc and PRAuc values ​​of BPNet, DanQ, DeepCNN, FCNA, FCNsignal and LSGKM, and the Wilson test p-value (paired) of the two indicators between FCNsignal and the other five methods, and the red triangle represents the average value;

[0041] Figure 6 A distribution diagram of the maximum signal values ​​in positive and negative samples of CEBPB, CUX1, TCF, and CTCF provided in an embodiment of the present invention;

[0042] Figure 7This is a comparison of the motif prediction performance of an FCNsignal and other comparative methods provided in an embodiment of the present invention; wherein (A) is the log2 (p-value) value of BPNet, DanQ, DeepCNN, FCNA, FCNA*, FCNsignal, MEME and STREME, and the Wilson test p-value (paired) of FCNsignal and other 7 methods. Please note that the red triangle represents the average value, and a value of 0 means that the target is not found; wherein, (B) is a detailed view of the identified motif of K562. The motifs in the inner ring are the target motifs, and the motifs in the outer ring are the discovered motifs. The size of the circle corresponds to the log2 (p-value) value, and different colors are used to represent different transcription factor categories. Motifs marked with pink belong to the same transcription factor family that shares a binding motif, while transcription factor motifs marked with other colors are more likely to be indirect transcription factor motifs;

[0043] Figure 8 A schematic diagram of the performance of an FCNsignal for six TF ChIP-seq data sets under four indicators provided in an embodiment of the present invention;

[0044] Fig. 9 The performance of FCNsignal in locating the DNA binding region of transcription factors provided in an embodiment of the present invention, wherein (A) is the result of FCNsignal locating the potential binding regions of CTCF and YY1 on the entire chromosome 17, and the distribution of H3K27ac signals in supporting and unsupporting sequences, wherein "support" means existing in the real peak, and vice versa means "unsupported"; (B) is the number of motif instances found by FIMO in the localization region and the negative region, and the distribution of the matching log2 (p-value); (C) is the motif enrichment analysis of the localization region and the negative region. The dotted line represents the result of the negative region; (D) is the base resolution signal, the real signal, and the H3K27ac signal of the DNA sequence of any length predicted by FCNsignal; the range of the DNA sequence above is chr17: 8376183-8384180 (7997bp length), and the range of the DNA sequence below is chr17: 41548466-41580329 (31863bp length);

[0045] Fig.10The performance comparison of an FCNsignal provided in an embodiment of the present invention and a comparative method on six ATAC-seq data sets; wherein, (A) is the MSE and Pearsonr values ​​of BPNet and FCNsignal; (B) is the Auc and PRAuc values ​​of LSGKM, DeepEmbed, DeopenBPNet and FCNsignal; (C) is the Pearsonr between the predicted maximum signal and the true maximum signal of MCF7 (Pearsonr: 0.534) and IMR90 (Pearsonr: 0.556); for MCF7, the Pearnsonr of the DNA sequence of chromatin openness is 0.075 and 0.634 respectively; for IMR90, the Pearnsonr of the DNA sequence of chromatin openness is -0.054 and 0.693 respectively;

[0046] Fig.11 Schematic diagram of FCNsignal successfully identifying disease-related causal SNPs from the LD group; (A) is the SNP scores for myeloma, pan-autoimmune disease and CLL by FCNsignal (top) and DeltaSVM (bottom); the risk variants for myeloma, pan-autoimmune disease and CLL are rs4487645, rs6927172 and rs539846, respectively; (B) is the distribution of SNP scores for cancer breast predicted by FCNsignal (left) and DeltaSVM (right); the SNP scores of the LD group predicted by FCNsignal are more concentrated in the low-score area than those predicted by DeltaSVM. DETAILED DESCRIPTION

[0047] In order to make the purpose, technical solutions and advantages of the embodiments of the present invention clearer, the technical solutions 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 part of the embodiments of the present invention, not all of the embodiments. Generally, the components of the embodiments of the present invention described and shown in the drawings here can be arranged and designed in various different configurations.

[0048] Therefore, the following detailed description of the embodiments of the present invention provided in the accompanying drawings is not intended to limit the scope of the invention claimed for protection, but merely represents selected embodiments of the present invention. Based on the embodiments of the present invention, all other embodiments obtained by ordinary technicians in this field without creative work are within the scope of protection of the present invention.

[0049] It should be noted that similar reference numerals and letters denote similar items in the following drawings, and therefore, once an item is defined in one drawing, further definition and explanation thereof is not required in subsequent drawings.

[0050] Example 1

[0051] like Figure 1 As shown, this embodiment provides a base resolution prediction method based on a fully convolutional neural network, comprising the following steps:

[0052] Data preparation stage S1: select the data set, obtain the DNA sequence and the corresponding base-resolution binding signal value;

[0053] FCNsignal model training stage S2: construct a fully convolutional neural network base resolution primitive mining model, which includes: an encoder structure, a decoder structure and a skip structure. The acquired data set is used to train the fully convolutional neural network base resolution primitive mining model, in which the base resolution combined with the signal value is used as a label for supervised learning;

[0054] Base resolution prediction stage S3: The DNA sequence to be predicted is input into the trained fully convolutional neural network base resolution primitive mining model to predict the base resolution signal.

[0055] The following is a detailed description of each stage.

[0056] 1. Data preparation stage

[0057] In this example, 53 ChIP-seq datasets for sequence-specific TFs of HeLa-S3, GM12878 and K562 cell lines, as well as 6 ATAC-seq datasets for A549, GM12878, HepG2, IMR90, K562 and MCF7 cell lines were downloaded from the ENCODE project. When selecting data, first select a dataset with overlapping biological experiments to ensure the quality of the data. Secondly, the peak value and signal are obtained after processing the standard data provided by ENCODE. The signal feature (signal) used in the experiment is calculated using the signal (signals) of the input data, and the specific calculation method is as follows:

[0058] signal=log 10 (1+signals)

[0059] Where signal is the base resolution combined signal value after scaling, and signals is the input base resolution combined signal value.

[0060] The positive sample is the peak data amplified to a length of 1000 bases and the corresponding signal, while the sequence of the same length 3000 bases upstream of the peak and the corresponding signal are selected as negative samples (please note that in order to ensure consistency with the real environment, the signal of the negative sample is not set to 0). When taking positive and negative samples, the experiment first considers maintaining similar background information between positive and negative samples, and secondly ensures that the negative sample does not bind to the target sequence and the binding signal value is also relatively low.

[0061] To test the ability of the FCNsignal model in locating transcription factor DNA binding regions, the entire chromosome 17 of hg38 was downloaded from UCSC, and all non-redundant peaks of two TFs (CTCF and YY1) were collected from ReMap2020. In addition, to verify the ability of the FCNsignal model in identifying causal single nucleotide polymorphisms (SNPs), several SNPs groups with strong linkage disequilibrium effects were collected. Candidate single nucleotide polymorphisms for myeloma were obtained from the corresponding literature, including 10 risk variants that change IRF4 binding at 7p15.3, 7 genetic susceptibility variants that affect 8 TFs binding at 6p23, and 27 risk variants that disrupt RELA binding at 15q15.1 for chronic lymphocytic leukemia (CLL). In addition, this example collected a set of genetic variations from previous studies, which contained 29 single nucleotide polymorphisms associated with breast cancer, and several single nucleotide polymorphisms had strong linkage disequilibrium effects (LD).

[0062] 2. FCNsignal model

[0063] FCNsignal is a symmetric deep neural network similar to U-net. The model takes DNA sequence as input and predicts base-resolution signals. FCNsignal can achieve multiple tasks simultaneously by using the predicted signal value, and the recognition accuracy is very high. The model consists of the following parts: (1) The encoder structure on the left is used to downsample the input and extract specific features of the sequence; (2) The decoder structure on the right is used to upsample the features and model the signal of each base; (3) The skip structure in the middle is responsible for combining the feature information learned in the encoder with the discriminant information in the decoder.

[0064] Specifically, Figure 2As shown in the figure, the encoder structure consists of three convolutional modules, each of which consists of a convolutional layer, a ReLU layer, a maximum pooling layer, a dropout layer, a bidirectional GRU layer, a dropout layer, and a global average pooling layer. Specifically, the convolutional block is used to gradually reduce the feature dimension and encode the specific features of the sequence, the bidirectional GRU layer is used to learn the long-term dependencies in the sequence-specific features, and the global average layer is used to capture the global context of the DNA sequence. Specifically, the calculation process of this part is as follows:

[0065] x (c) =ReLU(S (c) *W (c) +b (c) )

[0066] x (c) =MaxPool(x (c) )

[0067] x (c) = Dropout(x (c) ,ratio=0.2)

[0068] Where S (c) ,x (c) ,W (c) and b (c) They are the input, output, weight matrix, and bias vector of the cth convolutional block in the encoder structure; * represents the convolution operation.

[0069] r t =σ(W r x t +U r h t-1 )

[0070]

[0071] z t =σ(W z x t +U z h t-1 )

[0072]

[0073] where r t and z t They are the update gate and reset gate σ represents the sigmoid function, represents the element-wise multiplication operation; W and U are weight matrices.

[0074] h=h f +h b

[0075] y=GlobalAvePool(h)

[0076] where h f and h b It is the forward output and the backward output.

[0077] like Figure 2 As shown in Figure 1, the decoder structure consists of four upsampling layers and four mixing blocks. Each block consists of a batch normalization (BN) layer, a ReLU layer, and a convolutional layer. Specifically, the three jumps are used to combine the feature information learned in the encoder with the discriminative information in the decoder, and the upsampling layer is used to restore the size of the downsampled features. The mixing layer is used to rescale the value of the upsampled features. The calculation process of this part can be described by the following formula.

[0078] z (c) =Bilinear(y (c) )

[0079] z (c) =z (c) +x (c)

[0080] z (c) =ReLU(BN(z (c) ))*W (c) +b (c)

[0081] Where W (c) and b (c) are the weight matrix and bias vector of the cth mixing block in the decoder structure respectively; the Bilinear function is a bilinear interpolation operation; x (c) is the feature of the encoder from the same level. After the ReLU activation function, the final output z (c) This is the final prediction result of the FCNsignal signal. For the specific implementation of FCNsignal, please refer to the code in the Github link in Appendix C.

[0082] FCNsignal model training adopts a random split strategy to divide the experimental data set into training data, verification data, and test data. Specifically, for each data set, 80% of the peaks with corresponding signals are randomly selected as training data, and then 10% of them are randomly selected as verification data, and the remaining 20% ​​are used as test data. The FCNsignal model converts TFBSs prediction into a base-resolution regression task. The traditional binary classification error is not applicable, so the base-resolution MSE (mean square error) loss is used. The specific calculation method is as follows:

[0083]

[0084] Where MSE is the value of the MSE loss function, N is the number of samples in the training data, L is the length of each input sequence, and z i,j and are the predicted signal and the true signal respectively, α is the regularization parameter, which is used to weigh the fitting objective and the universality objective of the trained model, and ||*||2 represents L2 regularization.

[0085] The loss function is optimized by the Adam optimization algorithm with a batch size of 500. The gradient calculation uses the back-propagation algorithm, and the learning rate uses the exponential decay method with a decay rate of 0.9. The weights of the FCNsignal model are initialized by the Glorot initialization method. The learning rate parameter, Adam's beta parameter, and the regularization parameter are randomly selected from the pre-prepared sets [0.01, 0.001, 0.0001], [0.9, 0.99, 0.999], and [0, 0.001], respectively, and the dropout is randomly selected from 0.2 to 0.5. These hyperparameters are randomly sampled 15 times, and the validation data is used to select the hyperparameter set corresponding to the best performance.

[0086] 3. Other instructions

[0087] In this example, the openness of protein-DNA binding is defined as the number of reads falling into the binding region, which is called openness. To calculate this openness, the filtered bam files were first collected from ENCODE, and then the 'multicov' command in bedtools was used to calculate the number of reads falling into each peak. The openness values ​​of all peaks were further scaled with the log10 function as openness values.

[0088] Traditional prediction of transcription factor and DNA binding sites focuses on identifying whether a DNA sequence contains a motif sequence. This simple binary classification is not conducive to the specific analysis of motifs. In the study of this chapter, it is assumed that to some extent, the maximum value of the ChIP-seq signal can reflect the degree of openness of transcription factor binding to DNA. In order to verify this hypothesis, this example calculates the Pearson correlation coefficient (Pearson r) between the openness value and the maximum signal. Figure 3 As shown in Figure 1, the results of the four TF datasets confirm that the maximum value of the signal is significantly correlated with the openness of the transcription factor binding to DNA. Therefore, the maximum value of the signal can be used to distinguish between bound and non-bound sequences, because the openness of bound sequences is much higher than that of non-bound sequences, and can be used to locate the binding region, because the maximum value is most likely to be located in the binding region.

[0089] For the prediction of motifs, the position of the maximum value is first used to locate the binding region of length 100 bp. Then, the trained weights of the first convolutional layer are used to select the sub-regions with the highest score. Secondly, these selected sub-regions are arranged to calculate the corresponding position frequency matrices (PFMs). Finally, TOMTOM is used to match the PFMs with experimentally verified motifs in the standard database.

[0090] Given a specific cell line, how to calculate the score of SNPs is roughly divided into four steps. First, this chapter trains FCNsignal with ATAC-seq data of related cells. Second, the 1000bp region around the SNP position is determined, and the signal y of the corresponding reference sequence and the changed sequence is predicted respectively. i ref and i ref Third, the absolute value of the difference between the two signals at the SNP position is calculated as the SNP score, which is equal to Where i represents the position of the SNP. Finally, the SNP score is used to identify causal SNPs.

[0091] 4. Results

[0092] The FCNsignal model consists of an encoder structure, a decoder structure, and a skip structure. The model aims to predict the base-resolution signal of transcription factor DNA binding. In essence, FCNsignal is a base-resolution regression model that inputs DNA sequences, where each sequence is converted into a one-hot encoded matrix, A = [1, 0, 0, 0], C = [0, 1, 0, 0], G = [0, 0, 1, 0], T = [0, 0, 0, 1], and uses base-resolution binding signal values ​​as labels for supervised learning. In addition, FCNsignal can also be used to identify whether a DNA sequence is a transcription factor DNA binding region. In addition, the model locates the transcription factor DNA binding region and predicts the binding motif by using the maximum value of the binding signal. In the experiment, a total of 53 ChIP-seq and 6 ATAC-seq datasets were used to build the model and evaluate the prediction performance of FCNsignal.

[0093] like Figure 4 As shown in A in Figure 1, FCNsignal mainly consists of an encoder structure, a decoder structure and a skip structure. The model takes the DNA sequence as input and predicts the base resolution signal. FCNsignal can achieve multiple tasks at the same time by using the maximum value of the predicted signal, and the recognition accuracy is also improved. Figure 4As shown in B in Figure 1, FCNsignal can accurately capture the signals of the translocated binding regions located at 250 bp, 500 bp, and 750 bp. Figure 4 As shown in C, FCNsignal can accurately capture the signal of the binding region after randomly inserting negative examples.

[0094] FCNsignal can accurately capture transcription factor DNA binding motifs. In order to test the ability to capture transcription factor DNA binding patterns, two verification experiments were designed in this embodiment. In the first experiment, first find the maximum signal value of the peak; then, expand the length of the peak to 1000bp with the position of the maximum value as the center; finally, move the position of the maximum value by -250bp and +250bp respectively, and repeat the same process twice. For the second experiment, first find the maximum value of the peak signal; then extract a region with a length of 100bp around the maximum value position; finally, randomly insert the extracted region into the negative sample and keep it at 1000bp in length. Subsequently, use the trained FCNsignal to predict the base resolution signal of these data. Figure 4 B and Figure 4 As shown in C, it was found that the FCNsignal method proposed in this example can accurately identify the sequence-specific patterns of transcription factor DNA binding regions and model their corresponding signals.

[0095] 5. Experimental verification

[0096] 5.1. Overall performance of FCNsignal on ChIP-seq dataset

[0097] In order to test the overall performance of FCNsignal on ChIP-seq datasets, this example collected 53 ChIP-seq datasets of sequence-specific TFs from the ENCODE project, including 21, 20, and 12 TF datasets in GM12878, K562, and HeLa-S3 cell lines, respectively. The overall performance of the proposed method was studied on the test data.

[0098] The MSE and Pearson r between the predicted signal and the true signal are used to evaluate the performance of FCNsignal in modeling base resolution signals. The FCNA and BPNet methods are also introduced in the experiment for comparison. Figure 5 As shown in A in , the average MSE value of FCNsignal is lower than that of other methods, while the average Pearsonr value is higher than that of other methods. FCNsignal outperforms the other two methods in MSE and Pearsonr indicators. These results show that FCNsignal can model base resolution signals well.

[0099] The experiment selected LSGKM, DanQ, DeepCNN, FCNA and BPNet as comparison methods, and used Auc and PRAuc indicators to evaluate the performance of FCNsignal in distinguishing binding or non-binding sequences. Figure 5 As shown in B, FCNsignal outperforms the five comparison methods with a lower Wilson test p-value, and the average Auc and PRAuc values ​​are also higher than them. In addition, the experiment also found that FCNsignal performed much worse than LSGKM on three TFs, including CEBPB (Auc: 0.883 vs. 0.913), CUX1 (Auc: 0.837 vs. 0.882), and TCF7 (Auc: 0.885 vs. 0.944). In order to analyze this phenomenon, for each of them, the maximum value distribution of the signal in the positive and negative samples was calculated. As shown in Figure 6 As shown in the figure, for each of the three TFs, most of the maximum values ​​in the positive and negative samples are intersected, so it is difficult for FCNsignal to accurately distinguish positive or negative samples by the maximum value of the signal. In contrast, for CTCF, there is almost no intersection between the maximum values ​​in the positive and negative samples, so FCNsignal can accurately distinguish positive or negative samples, and its performance is better than LSGKM (Auc: 0.993 vs. 0.955). Overall, these observations indicate that FCNsignal can accurately distinguish the bound sequences, but it may also be affected by the signal quality.

[0100] In order to explore the impact of sample size on the performance of FCNsignal, the experiment first sorted all TF datasets from large to small by number, and then compared the three methods of LSGKM, FCNsignal and BPNet using the Auc indicator. The experiment found that: (1) LSGKM performs similarly to the FCNsignal model on small-scale datasets, or even exceeds FCNsignal, but its performance on large-scale datasets is inferior to that of the FCNsignal model; (2) The performance of BPNet on small-scale datasets is far behind that of LSGKM and FCNsignal models, but the results get better and better as the number of datasets increases. In general, these observations show that FCNsignal is less affected by the number of samples, but its performance improves as the number of samples increases.

[0101] 5.2. FCNsignal accurately identifies transcription factor DNA binding motifs

[0102] FCNsignal can accurately predict transcription factor DNA binding motifs. Seven comparison methods were selected in the experiment, including MEME, STREME, DanQ, DeepCNN, FCNA*, FCNA and BPNet, to study the performance of FCNsignal in predicting transcription factor DNA binding motifs.

[0103] The log2(p-value) generated by TOMTOM was used as the evaluation metric. Figure 7 As shown in A: (1) FCNsingal outperforms MEME, STREME, DanQ, DeepCNN and BPNet at low Wilson test p-value, and has a higher average log2(p-value) than them; (2) FCNA has competitive performance with FCNsignal, but has a lower average log2(p-value) than FCNsignal. (3) FCNA* performs significantly better than FCNsignal, and its average log2(p-value) is the highest among all methods; (4) FCNsignal discovers more motifs (50 / 53) than other comparison methods, and only FCNA* can discover all motifs. Both FCNA and FCNsignal are based on symmetric FCN structures, so their performance is comparable. FCNA* uses base-resolution labels manually annotated by PCMs (position counting matrices), while FCNsignal uses base-resolution signals generated by high-throughput sequencing technology, so FCNA* achieves very impressive performance. However, FCNA* is extremely dependent on PCMs, so it cannot distinguish between positive and negative samples, which means that it can get very good results even with negative samples, while FCNsignal can distinguish between positive and negative samples by the maximum value of the signal. In addition, unlike the signals generated by sequencing technology, real PCMs are not easy to obtain.

[0104] Collectively, these results demonstrate that FCNsignal can accurately and efficiently predict transcription factor DNA binding motifs.

[0105] In addition, FCNsignal can also be used to find indirect transcription factor DNA binding motifs. Taking all TFs in the K562 cell line as an example, in addition to each target TF motif, the top 5 matching motifs are picked out according to log2(p-value). Figure 7As shown in Figure B, the following observations were made experimentally: (1) FCNsignal can easily find similar motifs belonging to the same TF family. For example, the matching TF motifs of NFYA include NFYB and NFYC, which belong to the heteromeric CCAAT binding family. The matching TF motifs of GATA2 include GATA1, GATA3, and GATA4, which belong to the C4-GATA related family; (2) FCNsignal can find indirect TF motifs that interact with the target TF. For example, JUN-related motifs (JUNB and JUND) and Fos-related motifs (FOSL1 and FOSL2) are mutually inclusive, which has been shown that all JUN-FOS isoforms often strongly bind to TPA response elements. The matching TF motifs of GATA2 involve JUN and FOSB, and have significant p-values, which has been shown that GATA2 often occupies the same chromatin sites with the heterodimeric components of AP-1, c-JUN and c-FOS.

[0106] FCNsignal can locate potential transcription factor DNA binding sites across chromosomes. First, the cross-cell line prediction ability of FCNsignal was tested by using a model trained on datasets of GM12878 and K562 cell lines to predict datasets of HeLa-S3 cell lines, which included CTCF, ELK1, and JUND. The performance of FCNsignal in predicting datasets of different cell lines was slightly worse than its performance in predicting datasets of the same cell line, but the overall cross-cell line prediction performance was very high. The above observations mean that FCNsignal is able to predict specific transcription factor DNA binding sites across the entire genome.

[0107] In order to test the ability of FCNsignal in predicting potential transcription factor DNA binding sites on the entire chromosome, the experiment applied a chromosome segmentation strategy to divide the experimental data set into training data, validation data, and test data. First, six TF ChIP-seq data sets were selected from HeLa-S3 and K562 cell lines respectively; then, for each data set, the data of chromosomes 17 and 18 were used as test sets and validation sets, and the rest of the data were used as training sets; finally, these data were used to train and test FCNsignal.

[0108] The results are as follows Figure 8As shown in the figure, FCNsignal can achieve good prediction performance regardless of whether chromosome splitting or random splitting strategy is used. In addition, this section uses the entire chromosome 17 to verify the performance of the FCNsignal model in locating transcription factor DNA binding regions. The experimental process is as follows: (1) The entire chromosome 17 is divided into sequences of 1000 bp in length; (2) The training models of CTCF and YY1 are used to predict the base-resolution binding signals of these sequences; (3) Only when the maximum value of the predicted signal exceeds the manually set threshold (set to 1.5 in this experiment), several sequences with high confidence will be selected; (4) A region of 60 bp in length around the maximum value position is extracted. As a result, 2132 and 1060 potential regions containing CTCF and YY1 TFBSs were found, respectively.

[0109] In order to verify these localized short regions, the number of intersections between the localized regions and the real peaks collected by ReMap2020 was calculated based on the principle that the intersection rate exceeded 50%. Fig. 9 As shown in A, 81% (1725 / 2132) of CTCF and 63% (671 / 1060) of YY1 were supported by true peaks. In addition, the experiment analyzed the distribution of H3K27ac signals in supported and unsupported regions, and found the following: First, the H3K27ac values ​​of supported regions were significantly higher than those of unsupported regions, which means that supported regions mostly come from open chromatin regions, while unsupported regions may come from unopened regions; Second, unsupported regions also contain a small number of outliers with high H3K27ac values, which means that these outliers are similar to CTCF / YY1-related binding sites. Third, FIMO

[35] was used to find motif instances with high significance (p-value <1e-4), and CentriMo

[36] was used for motif enrichment analysis, analyzing the localization regions and negative regions respectively. Different from selecting localization regions, if the maximum value of the predicted signal is less than the threshold, these sequences are selected as negative regions, of which the same number as the localization regions are randomly selected. The results are shown in Fig. 9 As shown in Figure B, the number of motif instances found by FIMO in the positive region (1682 for CTCF and 514 for YY1) is much greater than that found in the negative region (117 for CTCF and 34 for YY1), and the log2(p-value) value of the positive region is much higher than that of the negative region. Fig. 9 As shown in Figure C, the motifs of CTCF and YY1 are also more enriched in the positive regions than in the negative regions. These results confirm that the positive regions are more likely to contain TFBSs. Overall, FCNsignal has the ability to classify and locate potential TF-DNA binding sites across chromosomes, but it is not easy to distinguish similar binding sites (such as those bound by TFs of the same family).

[0110] Since FCNsignal is an FCN-based framework, it can accept DNA sequences of arbitrary length. To account for this, the experiment randomly cut out some sequences of arbitrary length containing multiple peaks from chromosome 17. Then, the trained models of CTCF and YY1 were used to predict the base-resolution binding signals of these sequences. Fig. 9 As shown in D in Figure , FCNsignal can take sequences of arbitrary length as input and accurately predict binding signals with base resolution.

[0111] 5.3. Overall performance of FCNsignal on ATAC-seq dataset

[0112] In order to test the performance of FCNsignal on ATAC-seq datasets, this chapter collected 6 ATAC-seq datasets from the ENCODE project, including A549, GM12878, HepG2, IMR90, K562, and MCF7 cell lines. The ATAC-seq dataset was processed in the same way as the ChIP-seq dataset. As we all know, ATAC-seq data are often applied to the task of predicting chromatin accessibility, which is equivalent to identifying open or unopened regions of the entire genome. In this study, FCNsignal was used to predict chromatin accessibility. Simply put, FCNsignal inputs a DNA sequence and predicts its corresponding base-resolution signal, and then uses the maximum value of the signal to distinguish between open and unopened regions. Four competitive comparison methods, including LSGKM, DeepEmbed, Deopen, and BPNet, were used to compare with the FCNsignal model.

[0113] MSE and Pearsonr indicators are used to compare the performance of FCNsignal and BPNet in base resolution signal modeling. Fig.10 As shown in A, FCNsignal performs better than BPNet on all datasets. The Auc and PRAuc indicators are used to compare the performance of FCNsignal and the four comparison methods in distinguishing whether a sequence is open or not. Fig.10 As shown in B, except that BPNet outperforms FCNsignal on the A549 dataset, FCNsignal outperforms the four comparison methods on all datasets. In addition, when the number of samples is large, the performance of BPNet is comparable to that of FCNsignal, which confirms the above conclusion that BPNet performs poorly on small-scale datasets but obtains good results on large-scale datasets. In general, the above results show that FCNsignal performs better than other comparison methods on the ATAC-seq dataset.

[0114] In the above experiment, only the classification performance of FCNsignal was considered. However, the degree of accessibility of DNA sequences may be different from each other, even if they have the same label. The difference in the degree of DNA accessibility means that the classification model cannot distinguish between potential open regions of DNA with different accessibility. To explore this issue, the experiment still used the maximum value of the predicted signal to simulate the accessibility degree of ATAC-seq peaks. Therefore, the trained FCNsignal was used to predict the openness of DNA sequences from the original test data, which only included positive samples from different cell lines. Fig.10 As shown in Figure C, it is observed that FCNsignal performs well for DNA sequences with higher openness values ​​(average Pearson r: 0.62), but performs poorly for sequences with lower openness values ​​(average Pearson r: 0.07). The possible reason for this phenomenon is that the peak signal is easily affected by noise, especially for those sequences with weak binding signals. Therefore, it is still a huge challenge to identify potential open regions of DNA with different accessibility.

[0115] In addition, the enrichment of different transcription factor binding sites in each cell was analyzed by using the trained FCNsignal. In short, (1) all matching motifs were found through the motif prediction process; (2) the important matching motifs were selected by filtering those with p-value less than 1e-04; (3) the number of occurrences of each motif was calculated as the enrichment value. The heat map of these selected TF motifs in six cell lines was obtained, which clearly showed the enrichment of some TFBSs in different cell lines. For example, (1) CTCF and Sp1-like family (such as Sp1, Sp2) motifs are widely present in the six cell lines; (2) Ets-related family (such as ETV1, ETV6) motifs are enriched in GM12878 cell lines; (3) AP2-related family (such as TFAP2A, TFAP2B) motifs are enriched in MCF7 cell lines; (4) RXR-related family (such as NR2F1, NR2C2) motifs are enriched in HepG2 cell lines; (5) TEF1-related family (such as TEAD1, TEAD3) motifs are enriched in IMR90 cell lines. To explain this, the experiment collected all non-redundant peaks of five TFs including CTCF, ETV6, TFAP2A, NR2C2 and TEAD1 from ReMap2020 and calculated their intersection ratio with the ATAC-seq peaks of the six cell lines. The experiment showed that CTCF had a high intersection ratio with the six cell lines, while the maximum intersection ratios of other TFs (ETV6, TFAP2A, NR2C2, and TEAD1) were located in GM12878, MCF7, HepG2, and IMR90 cell lines, respectively, which is consistent with the above observations. In summary, the motifs learned by FCNsignal can reflect the enrichment degree of TFBSs in different cell lines.

[0116] 5.4. FCNsignal identifies causal SNPs from LD groups

[0117] FCNsignal can identify causal SNPs from LD groups To test the ability of FCNsignal to identify causal SNPs, four causal SNP datasets with strong LD effects were used in the experiment, and DeltaSVM was selected as a comparison method. For myeloma, pan-autoimmune disease, and CLL, FCNsignal trained on the GM12878 cell line was used to calculate the SNP scores. For breast cancer, FCNsignal trained on the MCF7 cell line was used to calculate the SNP scores. DeltaSVM calculated the SNP scores by directly utilizing the output of the LSGKM model on the GM12878 and MCF7 cell lines. Subsequently, these SNPs were prioritized by SNP scores, and the SNP with the highest score in each group was considered to be causal. Fig.11As shown in Figure A, both FCNsignal and DeltaSVM successfully identified the causal SNPs in the myeloma risk variant (rs4487645) and CLL risk variant (rs539846). However, DeltaSVM failed to correctly identify 7 pan-autoimmune genetic susceptibility candidate SNPs (rs6927172), while FCNsignal was able to accurately identify one of them. For breast cancer variants, Fig.11 B in the figure shows that the distribution of SNP scores in the LD group predicted by FCNsignal is more concentrated in the low-score region compared with the variants predicted by DeltaSVM. In addition, rs4784227 is considered to disrupt the binding of FOXA1, while rs3803662 is in obvious LD with it in individuals of European ancestry. Therefore, the ratio of SNP scores predicted by FCNsignal (0.054 / 0.001=54) is much higher than that predicted by DeltaSVM (2.28 / 1.54=1.5). The above results show that the FCNsignal model can identify causal SNPs from the LD group. In addition, BPNet can also successfully identify causal SNPs from the LD group, indicating that the framework for predicting base-resolution signals can be directly used to identify causal SNPs.

[0118] 6. Conclusion

[0119] This example proposes a model for high-precision identification of transcription factor binding sites based on a fully convolutional neural network, named the FCNsignal model. In addition to the input sequence data, the model also introduces binding signal data to achieve base-resolution-level prediction and analysis of motif-specific binding sites. Experimental results on 53 TF ChIP-seq datasets and 6 chromatin accessibility ATAC-seq datasets show that the FCNsignal model is superior to existing common methods. At the same time, this chapter also uses the FCNsignal model to predict DNA sequence motifs and locate all potential transcription factor binding regions on the entire chromosome. The results show that the model can find most of the known binding regions. In addition, experiments have also demonstrated the potential ability of the FCNsignal model in discovering disease-related single nucleotide polymorphisms (SNPs).

[0120] The preferred specific embodiments of the present invention are described in detail above. It should be understood that a person skilled in the art can make many modifications and changes based on the concept of the present invention without creative work. Therefore, any technical solution that can be obtained by a person skilled in the art through logical analysis, reasoning or limited experiments based on the concept of the present invention on the basis of the prior art should be within the scope of protection determined by the claims.

Claims

1. A base resolution prediction method based on a fully convolutional neural network, characterized in that: The following steps are involved: Obtain the DNA sequence to be predicted, input it into the pre-established and trained fully convolutional neural network base resolution primitive mining model, and predict the base resolution signal; The fully convolutional neural network base resolution primitive mining model includes: The encoder structure is used to downsample the input data and extract specific features of the sequence; The decoder structure is used to upsample the features and model the signal of each base; The skip structure is used to combine the features learned in the encoder structure with the discriminative information in the decoder structure; During the training process of the fully convolutional neural network base resolution primitive mining model, base resolution combined with signal value is used as a label for supervised learning.

2. A base resolution prediction method based on a fully convolutional neural network according to claim 1, characterized in that: The encoder structure includes a convolution module, a GRU layer and a global average aggregation layer connected in sequence, the number of the convolution modules is multiple, and each convolution module includes a convolution layer, a ReLU layer, a maximum aggregation layer and a dropout layer connected in sequence, which are used to gradually reduce the dimension of the feature and encode the specific features of the sequence; The GRU layer includes a bidirectional GRU layer and a dropout layer connected to each other, and the bidirectional GRU layer is used to learn long-term dependencies in sequence-specific features; The global average aggregation layer is used to capture the global background of DNA sequences.

3. A base resolution prediction method based on a fully convolutional neural network according to claim 2, characterized in that: The decoder structure includes a plurality of groups of upsampling and mixing block combinations connected in sequence, each group of upsampling and mixing block combinations includes an upsampling layer and a mixing block connected in sequence, and the mixing block includes a batch normalization layer, a ReLU layer and a convolution layer connected in sequence; The upsampling layer is used to restore the size of the downsampled features; The mixing block is used to rescale the values ​​of the upsampled features.

4. A base resolution prediction method based on a fully convolutional neural network according to claim 3, characterized in that: The number of convolutional modules in the encoder structure is three, and the number of upsampling and mixing block combinations in the decoder structure is four; The skip structure includes: the output of the first convolution module is connected to the input of the fourth group of upsampling and mixing block combination; the output of the second convolution module is connected to the input of the third group of upsampling and mixing block combination; the output of the GRU layer is connected to the input of the second group of upsampling and mixing block combination.

5. The method for predicting base resolution based on a fully convolutional neural network according to claim 1, characterized in that: The DNA sequence input to the fully convolutional neural network base-resolution motif mining model is converted into a one-hot encoded matrix.

6. The method for predicting base resolution based on a fully convolutional neural network according to claim 1, characterized in that: The base resolution MSE loss is used as the loss function in the training process of the fully convolutional neural network base resolution primitive mining model. The expression of the MSE loss is: Where MSE is the value of the MSE loss function, N is the number of samples in the training data, L is the length of each input sequence, and z i,j and are the predicted signal and the true signal respectively, α is the regularization parameter, and ||*||2 represents L2 regularization.

7. The method for predicting base resolution based on a fully convolutional neural network according to claim 1, characterized in that: The fully convolutional neural network base resolution primitive mining model scales the input base resolution combined signal value through the log10 function, and the corresponding scaling expression is: signal=log 10 (1+signals) Where signal is the base resolution combined signal value after scaling, and signals is the input base resolution combined signal value.

8. The method for predicting base resolution based on a fully convolutional neural network according to claim 1, characterized in that: During the training process of the fully convolutional neural network base-resolution primitive mining model, the Adam optimization algorithm is used to perform parameter optimization.

9. The method for predicting base resolution based on a fully convolutional neural network according to claim 1, characterized in that: The method further comprises distinguishing the DNA sequence to be predicted into a binding sequence and a non-binding sequence according to the maximum value of the predicted base resolution signal.

10. A method for predicting base resolution based on a fully convolutional neural network according to claim 1, characterized in that: The method also includes locating a binding region with a length of 100bp according to the position of the maximum value of the predicted base resolution signal, selecting the sub-region with the highest score according to the training weight of the first convolutional layer in the encoder structure, arranging the selected sub-regions, calculating the corresponding position frequency matrix, and finally matching the position frequency matrix with the primitives in the standard database to obtain the prediction result of the primitives.

Images