TADs structural domain prediction method and system based on attribute graph node clustering

Abstract

The invention belongs to the technical field of biological computing, and relates to a TADs structural domain prediction method and system based on attribute graph node clustering. The method comprises the following steps: comparing Hi-C sequencing data with a reference genome, constructing and preprocessing a Hi-C contact matrix, constructing a chromatin attribute graph model, taking chromatin equal-length segments as graph nodes, taking interaction among the nodes as edges, and splicing node attributes by A / B compartment structures and DNA sequence features; the node attribute values serve as new nodes to be inserted into the graph, node amplification is carried out, and an extended attribute graph is formed; performing node clustering on the extended attribute graph, and dividing closely connected structural nodes into the same cluster to serve as a preliminary TADs structural domain; and according to biological prior knowledge, carrying out continuity and size constraint processing on the preliminary TADs, and finally obtaining a predicted TADs structural domain. According to the method, various characteristic information of chromatin is fused through a graph model, and biological priori knowledge is combined, so that the TADs structural domain is effectively predicted.

Description

Technical Field

[0001] The present invention belongs to the field of biological computing technology and relates to a TADs domain prediction method and system based on attribute graph node clustering. Background Art

[0002] High-throughput chromosome conformation capture is a sequencing technology used to study the spatial conformation of chromatin. It cross-links, cuts, connects, and sequences chromatin to obtain the Hi-C map of each chromatin in the entire genome. Generally speaking, the Hi-C map is organized into a symmetrical square matrix, called the Hi-C contact matrix. The subscripts of the matrix elements represent different chromatin fragments, and the values ​​of the matrix elements represent the interaction frequency of the chromatin fragments at the corresponding positions. Chromatin fragments with higher interaction frequencies are closer in three-dimensional space. On the Hi-C map, people have discovered many key chromatin domains, including A / B compartment structures at the Mb (megabase pair) level, subcompartment structures, and topologically associated domains (TADs) at the sub-Mb (sub-megabase pair) level, from large to small, which together constitute the complex domain hierarchy of chromatin. According to the organizational form of TADs domains, the prediction methods can be divided into non-intersecting TADs domain prediction and hierarchical TADs domain prediction.

[0003] Existing TADs domain prediction methods only use the interaction information of chromatin fragments in Hi-C maps, and lack the use of the A / B compartment structure information of the chromatin fragments and the DNA sequence information contained in the chromatin fragments. Their impact on TADs domain prediction is reflected in the following two aspects: (1) The A / B compartment is a higher-order chromatin domain, and the TADs domain, as its internal fine structure, is closely related to the A / B compartment boundary. Existing methods only use the local interaction information of the Hi-C contact matrix, which will cause the loss of the global characteristics of the TADs domain.

[0004] (2) DNA sequences are associated with many transcription factors, and biological experiments have shown that many transcription factors (such as CTCF) tend to be enriched at the boundaries of TADs domains. However, existing methods do not fully consider the impact of DNA sequence characteristics on TADs domains.

[0005] In view of the above reasons, the accuracy of TADs prediction methods in existing Hi-C data is still not high. Therefore, it is necessary to develop a TADs domain prediction method based on Hi-C data that simultaneously utilizes A / B compartment structure information and DNA sequence information. Summary of the invention

[0006] The technical problem to be solved by the present invention is to provide a TADs domain prediction method and system based on attribute graph node clustering in view of the deficiencies in the above-mentioned prior art, so as to solve the technical problem of low accuracy of current TADs domain prediction.

[0007] The objective of the present invention is achieved by the following technical solutions: In a first aspect, the present invention provides a TADs domain prediction method based on attribute graph node clustering, comprising: Generate a Hi-C contact matrix using the paired-end sequencing sequences obtained by Hi-C sequencing and the reference genome, and pre-process the Hi-C contact matrix; Constructing a chromatin attribute graph model; in the attribute graph model, the data structure represents chromatin, chromatin equal-length segments are used as nodes of the attribute graph, the interactions between the nodes are used as edges between the nodes, and the A / B compartment structural features of the chromatin segments and the DNA sequence features are spliced ​​as the attributes of the nodes; Inserting the attribute value of the node as a new node into the attribute graph model, expanding the nodes of the attribute graph, and forming an extended attribute graph; Cluster the nodes in the extended attribute graph and divide them into different clusters. The tightly connected structural nodes are divided into the same cluster as the preliminary TADs structural domain. According to the continuity of the preliminary TADs domain and the minimum allowed TADs domain size, the TADs domain is further processed to obtain the final predicted TADs domain.

[0008] As a further improvement of the present invention, the double-end sequencing sequence measured by Hi-C sequencing and the reference genome are used to generate a Hi-C contact matrix, which specifically includes: Use Hi-C sequencing technology to obtain paired-end sequencing sequences of chromatin; The double-end sequencing data were compared with the reference genome using a sequence alignment device. Low-quality sequences and PCR duplications were filtered out based on the quality scores of the comparison. The number of sequencing sequences in each segment was counted according to the set chromatin segment length to generate a Hi-C contact matrix.

[0009] As a further improvement of the present invention, the Hi-C contact matrix is ​​pre-processed, specifically comprising: The Hi-C contact matrix is ​​normalized using a Hi-C contact matrix normalization method; The values ​​in the normalized Hi-C contact matrix were logarithmically transformed to obtain the logarithmically transformed Hi-C contact matrix.

[0010] As a further improvement of the present invention, chromatin equal-length segments are used as nodes of the attribute graph, and the interactions between the nodes are used as edges between the nodes. The process of extracting the A / B compartment structure and DNA sequence features specifically includes: The equal-length segments of chromatin are taken as nodes of the attribute graph, and the set of all equal-length segments corresponds to the set of nodes of the attribute graph; The node pairs whose contact frequency in the Hi-C contact matrix is ​​higher than the 80% quantile of all contact frequencies are retained as the edge set of the graph.

[0011] As a further improvement of the present invention, the A / B compartment structural features and DNA sequence features of the chromatin segment are spliced ​​as attributes of the node, specifically including: The method for extracting the structural features of the A / B compartment is as follows: first, the average interaction frequency of all chromatin segments at each genomic distance is calculated to form a contact frequency attenuation vector along the genomic distance, and the attenuation vector is used to normalize the Hi-C contact matrix by genomic distance; then, the Pearson correlation coefficient of the two column vectors of the normalized matrix is ​​calculated to form a correlation coefficient matrix; the correlation coefficient matrix is ​​subjected to principal component analysis, and the principal component of each chromatin segment is extracted as the structural features of the A / B compartment; The method for extracting DNA sequence features is as follows: using the K-mer feature extraction method to divide the DNA sequence into lengths K, counting the frequencies of all subsequences of length K in the chromatin segment; sorting the frequencies to form a frequency vector; adding the frequency vectors, and selecting the feature with the highest frequency as the DNA sequence feature.

[0012] As a further improvement of the present invention, the attribute value of the node is inserted into the attribute graph model as a new node, and the node of the attribute graph is expanded to form an extended attribute graph, which specifically includes: Normalize all node attribute values ​​to between 0 and 1, and discretize and sample the normalized node attribute values ​​with a set step size according to the settings; An attribute node is created for each sampled discrete value, the attribute node is inserted into the attribute graph model, attribute edges are inserted based on the attribute nodes, and an extended attribute graph is formed according to the added attribute nodes and attribute edges; the attribute edge is used to connect an attribute node and a structure node, and the structure node refers to an original node in the attribute graph model.

[0013] As a further improvement of the present invention, clustering the nodes in the extended attribute graph specifically includes: Initialize edge weights; the edge weights are obtained according to the interaction frequency between nodes; Calculate the random walk distance matrix of nodes based on the attribute graph structure and edge weights; According to the density of nodes, the k structural nodes with the highest node density are determined as the starting clustering centers; Iteratively update the cluster to which each node belongs, the cluster center, the weight of the attribute edge, and the random walk distance matrix until convergence.

[0014] As a further improvement of the present invention, according to the continuity of the preliminary TADs domain and the minimum allowed TADs domain size, the TADs domain is further processed, specifically including: If the nodes in the same cluster are continuously distributed on the genome, and the length of the continuous fragments on the genome is not less than the set minimum size of the TADs domain, the initial TADs domain will be retained; If nodes in the same cluster are isolated by other clusters on the genome and present segmented distribution, and the length of the continuous segments on the genome is not less than the set minimum size of the TADs domain, the separated cluster will be split into two non-adjacent TADs domains; The nodes in the same cluster are isolated by other clusters on the genome and present a segmented distribution. If there is a continuous fragment whose length is less than the set minimum size of the TADs domain, and the two sides of the continuous fragment less than the set minimum size of the TADs domain belong to the same cluster, then the fragment will be merged into the clusters on both sides to form a larger TADs domain. The nodes in the same cluster are isolated by other clusters on the genome and present a segmented distribution. If there is a continuous fragment whose length is less than the set minimum size of the TADs domain, and the two sides of the continuous fragment less than the set minimum size of the TADs domain belong to different clusters, the average random walk distance between the fragment and the clusters on both sides is evaluated, and it is merged into the cluster with a larger random walk distance to form a larger TADs domain.

[0015] In a second aspect, the present invention provides a TADs domain prediction system based on attribute graph node clustering, which is used to implement the above-mentioned TADs domain prediction method based on attribute graph node clustering, comprising: Hi-C contact matrix extraction module, which is used to compare paired-end sequencing sequences with the reference genome, filter low-quality sequences and PCR duplications, count DNA fragments, and generate Hi-C contact matrices; Data preprocessing module, used to preprocess Hi-C contact matrix; An attribute graph model generation module is used to map chromatin segments and their interaction information onto an attribute graph model; in the attribute graph model, the data structure represents chromatin, chromatin equal-length segments are used as nodes of the attribute graph, the interactions between nodes are used as edges between nodes, and the A / B compartment structural features of the chromatin segments and the DNA sequence features are spliced ​​as attributes of the nodes; the attribute values ​​of the nodes are inserted into the attribute graph model as new nodes, and the attribute graph is amplified to form an extended attribute graph; The node clustering module is used to cluster densely connected structural nodes into the same category according to the structure of the generated attribute graph model to obtain preliminary TADs domains; based on the continuity of the preliminary TADs domains and the minimum allowed TADs domain size, the TADs domains are further processed to obtain the final predicted TADs domains.

[0016] As a further improvement of the present invention, it also includes a TADs domain post-processing and visualization module, which is used to further process the TADs domain according to the biological priors of TADs continuity and the minimum size of TADs, and visualize the predicted TADs domain.

[0017] The beneficial effects of the present invention are as follows: the TADs domain prediction method based on attribute graph node clustering of the present invention models chromatin as a graph with node attributes, makes full use of data related to chromatin structure, and reduces information loss in the modeling process, that is, the Hi-C contact matrix is ​​modeled as an adjacency matrix of the graph, and the high-order A / B compartment structural characteristics and the DNA sequence characteristics of the chromatin fragments are modeled as node attributes of the graph; the TADs domain prediction problem is reasonably converted into a clustering problem of graph nodes with attributes, comprehensively considers the structure and attribute similarity of the nodes, effectively handles the complex internal topology of the TADs domain, and improves the accuracy of TADs domain prediction. BRIEF DESCRIPTION OF THE DRAWINGS

[0018] In order to more clearly illustrate the embodiments of the present invention or the technical solutions in the prior art, the drawings required for use in the embodiments or the description of the prior art will be briefly introduced below. Obviously, the drawings described below are some embodiments of the present invention. For ordinary technicians in this field, other drawings can be obtained based on these drawings without paying creative work.

[0019] Figure 1 It is a flowchart of a TADs domain prediction method based on attribute graph node clustering provided in an embodiment of the present invention.

[0020] Figure 2 It is a flowchart of node attribute extraction of the chromatin graph model provided by an embodiment of the present invention.

[0021] Figure 3 It is a schematic diagram of the TADs domain post-processing method provided in an embodiment of the present invention.

[0022] Figure 4 It is a system structure diagram of a TADs domain prediction system based on attribute graph node clustering provided in an embodiment of the present invention.

[0023] Figure 5 It is a heat map of TADs domains predicted by the Hi-C data of chromatin No. 1 of IMR90 cell line applied in the present invention. DETAILED DESCRIPTION

[0024] In order to make the purpose and technical solution of the present invention clearer and easier to understand, the present invention is further described in detail below in conjunction with the accompanying drawings and embodiments. The specific embodiments described herein are only used to explain the present invention and are not used to limit the present invention.

[0025] The technical solution of the present invention will be clearly and completely described below in conjunction with the accompanying drawings and specific embodiments, wherein the described embodiments are only a part of the embodiments of the present invention, rather than all the embodiments.

[0026] Embodiment 1 like Figure 1-Figure 3 As shown, this embodiment provides a TADs domain prediction method based on attribute graph node clustering. The following is a specific implementation method.

[0027] Step 1: Generate a Hi-C contact matrix using the paired-end sequencing sequences obtained by Hi-C sequencing and the reference genome, and pre-process the Hi-C contact matrix.

[0028] Specifically, the Hi-C contact matrix generation process is used to convert the paired-end sequencing sequences measured by Hi-C into a Hi-C contact matrix that can be processed and visualized by users. The preprocessing process of the Hi-C contact matrix includes normalization and logarithmic transformation. Normalization is used to remove systematic biases introduced during the Hi-C sequencing process, such as enzyme cutting preference, uneven distribution of GC content and length of paired-end short sequences, etc.

[0029] In this implementation, bowtie2 software was used to establish an index of the reference genome, and bowtie was used to compare the paired-end sequencing sequences to the reference genome; HiCPro software was used to filter out sequences with low quality scores and PCR duplicates; the length of the equal-length segment was specified to be 50 Kb, and the number of sequencing sequences in each position segment was counted according to the length of the equal-length segment to generate a Hi-C contact matrix.

[0030] In this implementation, explicit modeling methods such as HiCNorm or implicit modeling methods such as Knight-Ruiz (KR) can be used to standardize the contact matrix to remove systematic bias. The values ​​of the Hi-C contact matrix are transformed to a logarithmic scale, and the base of the logarithmic transformation used is 10, so that the distribution of the values ​​in the matrix is ​​closer to a normal distribution, which is convenient for subsequent TADs domain prediction.

[0031] In this implementation, the HiCNorm normalization method is based on Poisson regression, which uses the effective length of the paired-end short sequences and the GC content of the reference genome to estimate the deviation of the contact matrix, and uses the residual of the Poisson regression as the normalized Hi-C contact matrix.

[0032] In this implementation, KR normalization is an implicit normalization that does not require the GC content characteristics of the genome to be known in advance. It is based on a matrix balancing strategy and aims to obtain a double random matrix in which the sum of the rows and columns of the matrix is ​​equal to 1.

[0033] Step 2: Construct an attribute graph model of chromatin; in the attribute graph model, the data structure represents chromatin, and equal-length segments of chromatin are used as nodes of the attribute graph. The interactions between nodes are used as edges between nodes, and the A / B compartment structural characteristics of the chromatin segments and the DNA sequence characteristics are spliced ​​as the attributes of the nodes.

[0034] The set of nodes and edges of the graph corresponds to the chromatin in the following way: the equal-length segments of the chromatin are used as the nodes of the graph, and the set of all equal-length segments corresponds to the set of nodes of the graph V Considering that TADs domains are mainly distributed near the diagonal of the Hi-C contact matrix, their interaction frequency is significantly higher than the contact frequency at positions far from the diagonal. Therefore, the node pairs with contact frequencies higher than the 80% quantile of all contact frequencies in the Hi-C contact matrix are retained as the edge set of the graph. E .

[0035] The node attributes of the graph include A / B compartment structure features and DNA sequence features. The extraction process of A / B compartment structure and DNA sequence features is as follows: Figure 2 shown.

[0036] For the A / B compartment structural features, we first calculated the average interaction frequencies of all chromatin segments at each genomic distance to form a vector of contact frequency decay along the genomic distance, and used the decay vector to normalize the contact matrix by genomic distance; then we calculated the Pearson correlation coefficients of the two column vectors of the normalized matrix to form a correlation coefficient matrix; finally, we performed principal component analysis on the correlation coefficient matrix, and extracted the principal components of each chromatin segment as the A / B compartment structural features. In this implementation, the number of principal components was 10.

[0037] For DNA sequence features, a fixed-length nucleotide string (K-mer) feature extraction method is used. In order to achieve a balance between computational complexity and prediction accuracy, this embodiment takes K=5. The 5-mer feature extraction method counts the frequency of all possible subsequences of length 5 in the chromatin segment, so the dimension of the feature is 4. 5 = 1024. In order to further reduce the feature dimension of the DNA sequence, the 5-mer feature vectors of all segments are summed, and the features with the highest frequency of occurrence are selected as the DNA sequence features after dimensionality reduction. In this implementation, the feature dimension of the DNA sequence is consistent with the feature dimension of the A / B compartment structure.

[0038] Furthermore, the reduced-dimensional A / B compartment structural features are concatenated with DNA sequence features as node attributes of the chromatin graph model, and the attribute set of all nodes is recorded as A , thus obtaining a complete chromatin attribute graph model G =( V , E , A ),in , m It is the node attribute dimension.

[0039] Step 3: Insert the attribute value of the node as a new node into the attribute graph model, expand the nodes of the attribute graph, and form an extended attribute graph.

[0040] In order to observe the similarity between two chromatin fragments (i.e., nodes in the graph) under a unified metric, the attribute values ​​of the nodes are inserted as new nodes into the attribute graph structure established in step 2. In order to compress the number of amplified nodes and reduce the computational complexity, all attribute values ​​are normalized to between 0 and 1 using the min-max normalization method and discretely sampled with a step size of 0.1. i Dimension properties a i The set of all possible values ​​is Dom ( a i )= , then the enlarged image is G =( V ∪ V a , E ∪ E a ),in V a = is the set of augmented attribute nodes, E a It is a set of augmented attribute edges, used to connect an attribute node and a structure node. The specific operation is: if a structure node vi The value of a jk ,but( v i , a jk ) node pair is an attribute edge. It can be seen that the more attribute values ​​two nodes share, the more likely they are to be connected through the newly added attribute nodes. Therefore, through the expansion of attribute nodes, the attributes and structures in the original graph are organized in one graph, and the nodes corresponding to chromatin fragments with similar A / B compartment structural characteristics and DNA sequence characteristics are more likely to be connected by edges.

[0041] After attribute node amplification, a new graph is obtained. The nodes of the graph include structure nodes corresponding to chromatin segments one by one and attribute nodes corresponding to all node attribute values. The edges of the graph include structure edges connecting structure nodes and attribute edges connecting structure nodes and attribute nodes. Attribute edges are used to connect an attribute node and a structure node. Among them, structure nodes refer to the original nodes in the attribute graph model.

[0042] Step 4: Cluster the nodes in the extended attribute graph and divide them into different clusters. The tightly connected structural nodes are divided into the same cluster as the preliminary TADs domain.

[0043] Specifically, the TADs domain prediction problem can be transformed into a node clustering problem of the graph after attribute node expansion, that is, the structural nodes with more dense connections are divided into the same cluster as the preliminary TADs domain. The TADs domain prediction algorithm based on node clustering is shown in Table 1: Table 1: TADs domain prediction algorithm based on node clustering

[0044] Probability transfer matrix in the algorithm P According to the figure G a The edge set E ∪ E a and weight w 0, w 1, w 2,…, w m The construction is divided into the following four cases: Case (1): From a structure node v i Transfer to another structure node v j Probability Calculated by the following formula:

[0045] Case (2): From a structure node v i Move to an attribute node v jk Probability Calculated by the following formula:

[0046] Case (3): From an attribute node v jk Move to a structure node v i Probability Calculated by the following formula:

[0047] Case (4): There is no edge connection between attribute nodes, so the transition probability between attribute nodes is 0.

[0048] Calculate the probability transfer matrix one by one according to situations (1)-(4) P The element values ​​in , further, calculate l The distance matrix of a random walk:

[0049] Nodes in the algorithm v i The node density D ( i ) refers to the sum of the influence of the node on all nodes in the graph, which can be calculated by the following formula:

[0050] The method for selecting the cluster center point is: sort all the nodes according to the node density, and select the node with the highest node density. k Nodes are used as the initial cluster centers.

[0051] The node allocation method is: v i Assigned to the cluster center point with the maximum random walk distance, that is, v i Assign to the center point c * Cluster C * , c * The following conditions must be met:

[0052] No. t +1 repeat midpointv i The updating method of the attribute weight is:

[0053] in m is the dimension of the node attribute in this implementation, vote i For attributes a i The voting function indicates whether two nodes are in the attribute a i share the same value, vote i is defined as:

[0054] Step 5: Based on the continuity of the preliminary TADs domain and the minimum allowed TADs domain size, the TADs domain is further processed to obtain the final predicted TADs domain.

[0055] The preliminary TADs domain obtained by node clustering may appear discontinuous on the genome, or the TADs domain may be too small. In order to reduce the false positive rate, the clustering results are further processed according to biological prior knowledge to obtain the final TADs domain. Figure 3 As shown, each symbol occupies an equal-length chromatin segment, the minimum size of the TADs domain is set to 3 segments, and the post-processing methods are given in the following four cases: (1) If all nodes in the same cluster are distributed continuously on the genome, and the length of the continuous fragments on the genome is not less than the set minimum size of the TADs domain, the initial TADs domain is retained, such as Figure 3 As shown in situation 1.

[0056] (2) If nodes in the same cluster are isolated by other clusters on the genome and present a segmented distribution, and the length of the continuous fragments on the genome is not less than the set minimum size of the TADs domain, the separated cluster will be split into two non-adjacent TADs domains, such as Figure 3 In case 2, TAD2 and TAD4, TAD3 and TAD5.

[0057] (3) If the nodes in the same cluster are isolated by other clusters on the genome and present a segmented distribution, and there is a continuous fragment whose length is less than the set minimum size of the TADs domain, and the two sides of the continuous fragment less than the set minimum size of the TADs domain belong to the same cluster, then the fragment is merged into the clusters on both sides to form a larger TADs domain, such as Figure 3 TAD2 in case three.

[0058] (4) If the nodes in the same cluster are isolated by other clusters on the genome and present a segmented distribution, and there is a continuous fragment whose length is less than the set minimum size of the TADs domain, and the two sides of the continuous fragment less than the set minimum size of the TADs domain belong to different clusters, then the average random walk distance between the fragment and the clusters on both sides is evaluated, and it is merged into the cluster with a larger random walk distance to form a larger TADs domain, such as Figure 3 As shown in case 4, if the segment has a larger average random walk distance with the left cluster, it is merged into TAD2, otherwise, it is merged into TAD3.

[0059] Embodiment 2 This embodiment provides a TADs domain prediction device based on attribute graph node clustering. The system structure diagram of the TADs domain prediction device is as follows: Figure 4 As shown, it includes the following modules: The Hi-C contact matrix extraction module is used to compare paired-end sequencing sequences with the reference genome, filter low-quality sequences and PCR duplicates, count DNA fragments, and generate Hi-C contact matrices.

[0060] The data preprocessing module is used to preprocess the Hi-C contact matrix; specifically, the Hi-C contact matrix is ​​explicitly or implicitly normalized to remove systematic biases, and the elements of the Hi-C contact matrix are logarithmically transformed.

[0061] An attribute graph model generation module is used to map chromatin segments and their interaction information onto an attribute graph model; in the attribute graph model, the data structure represents chromatin, chromatin equal-length segments are used as nodes of the attribute graph, the interactions between nodes are used as edges between nodes, and the A / B compartment structural features of the chromatin segments and the DNA sequence features are spliced ​​as attributes of the nodes; the attribute values ​​of the nodes are inserted into the attribute graph model as new nodes, and the attribute graph is amplified to form an extended attribute graph; The node clustering module is used to cluster the densely connected structural nodes into the same category according to the structure of the generated attribute graph model, and obtain the preliminary TADs domain; further process the TADs domain according to the continuity of the preliminary TADs domain and the minimum allowed TADs domain size, and obtain the final predicted TADs domain. Specifically, this module is used to cluster the densely connected structural nodes into the same category according to the structure of the generated graph model, and obtain the preliminary TADs domain; calculate the probability transfer matrix and the random walk distance matrix of the nodes according to the graph structure and weights; assign clusters to nodes according to the distance from the node to the cluster center node; and update the attribute weight according to the contribution of each attribute to the clustering result.

[0062] The system also includes a TADs domain post-processing and visualization module, which is used to further process the preliminary TADs domain according to the biological priors of TADs continuity and the minimum size of TADs to obtain the final predicted TADs domain; the Hi-C heat map is visualized using the heatmap visualization function in python or matlab, and the identified TADs domains are annotated on the heat map with squares or isosceles right triangles.

[0063] The prediction device of TADs domain based on attribute graph node clustering of this embodiment corresponds to the prediction method of TADs domain based on attribute graph node clustering of the first embodiment mentioned above; wherein the Hi-C contact matrix extraction module and the data preprocessing module correspond to step 1, the graph model generation module corresponds to steps 2 and 3, the node clustering module corresponds to step 4, and the TADs domain post-processing part in the TADs domain post-processing and visualization module corresponds to step 5, so it will not be repeated here.

[0064] Taking the TADs domain prediction based on the published Hi-C data of chromatin No. 1 of the human IMR90 cell line as an example, the TADs prediction method and effect of the present invention are introduced.

[0065] (1) The short sequences of paired-end sequencing were aligned with the human reference genome hg19, and low-quality alignment sequences and PCR duplicates were filtered out. A Hi-C contact matrix with a segment length of 50 Kb was extracted, and the HiCNorm normalization method was used to remove systematic biases in the matrix. A logarithmic transformation was used to make the data distribution closer to a normal distribution.

[0066] (2) An attribute graph model of chromatin No. 1 of the IMR90 cell line was established. Each 50Kb long chromatin fragment was mapped as a node of the graph, and the fragment pairs with interaction frequency greater than the 80% quantile of all interaction frequencies were mapped as edges of the graph. The Hi-C contact matrix was normalized by genomic distance, and the Pearson correlation coefficient matrix was calculated and principal component analysis was performed to extract the A / B compartment structure information. The K-mer features of the DNA sequence were calculated, and the first 10 dimensions were extracted according to the frequency sorting. The A / B compartment structure features and the DNA sequence structure features were combined to form the attribute nodes of the graph.

[0067] (3) The graph model is expanded with attribute nodes. The nodes that share the same node attributes in the initial graph model are connected through corresponding attribute nodes, so that the attributes and structures in the original graph are organized in a unified graph.

[0068] (4) Perform node clustering on the amplified graph. For the chromatin of IMR90 cell line No. 1, the total length is 249 Mb, and the average size of the TADs domain is 500 Kb, so the expected number of TADs is =498;Set k =498 pairs of amplified graph nodes were clustered to obtain 498 initial TADs domains.

[0069] (5) The clustering results were further post-processed based on biological prior knowledge, and a total of 457 final predicted TADs domains were obtained after post-processing. Since the dimension of the Hi-C contact matrix is ​​too high, we selected a section (150Mb-155Mb) for visualization and marked the predicted TADs domains with isosceles right triangles. Figure 5 As shown, it can be seen that the interaction frequency within the TADs domain on the Hi-C contact matrix predicted by the method disclosed in this patent is significantly higher than that outside the TADs domain, and it has a certain robustness to the regions with missing alignments in Hi-C sequencing (white strips in the heat map), which illustrates the effectiveness of the present invention.

Claims

1. A TADs domain prediction method based on attribute graph node clustering, characterized in that: include: Generate a Hi-C contact matrix using the paired-end sequencing sequences obtained by Hi-C sequencing and the reference genome, and pre-process the Hi-C contact matrix; Construct a property graph model of chromatin; The data structure in the attribute graph model represents chromatin, and equal-length segments of chromatin are used as nodes of the attribute graph, the interactions between nodes are used as edges between nodes, and the A / B compartment structural features of the chromatin segments and the DNA sequence features are spliced ​​as the attributes of the nodes; Inserting the attribute value of the node as a new node into the attribute graph model, expanding the node of the attribute graph, and forming an extended attribute graph; Cluster the nodes in the extended attribute graph and divide them into different clusters. The tightly connected structural nodes are divided into the same cluster as the preliminary TADs domain. According to the continuity of the preliminary TADs domain and the minimum allowed TADs domain size, the TADs domain is further processed to obtain the final predicted TADs domain.

2. The TADs domain prediction method based on attribute graph node clustering according to claim 1, characterized in that: The double-end sequencing sequence measured by Hi-C sequencing and the reference genome are used to generate a Hi-C contact matrix, specifically including: Use Hi-C sequencing technology to obtain paired-end sequencing sequences of chromatin; The double-end sequencing data were compared with the reference genome using a sequence alignment device. Low-quality sequences and PCR duplications were filtered out based on the quality scores of the comparison. The number of sequencing sequences in each segment was counted according to the set chromatin segment length to generate a Hi-C contact matrix.

3. The TADs domain prediction method based on attribute graph node clustering according to claim 1 is characterized in that: Preprocess the Hi-C contact matrix, including: The Hi-C contact matrix is ​​normalized using a Hi-C contact matrix normalization method; The values ​​in the normalized Hi-C contact matrix were log-transformed to obtain the log-transformed Hi-C contact matrix.

4. The TADs domain prediction method based on attribute graph node clustering according to claim 3 is characterized in that: The chromatin segments of equal length are used as nodes of the attribute graph, and the interactions between nodes are used as edges between nodes. The process of extracting A / B compartment structure and DNA sequence features specifically includes: The equal-length segments of chromatin are taken as nodes of the attribute graph, and the set of all equal-length segments corresponds to the set of nodes of the attribute graph; The node pairs whose contact frequency in the Hi-C contact matrix is ​​higher than the 80% quantile of all contact frequencies are retained as the edge set of the graph.

5. The TADs domain prediction method based on attribute graph node clustering according to claim 4 is characterized in that: The A / B compartment structural characteristics of the chromatin segment and the DNA sequence characteristics are spliced ​​as node attributes, including: The method for extracting the structural features of the A / B compartment is as follows: first, the average interaction frequency of all chromatin segments at each genomic distance is calculated to form a contact frequency attenuation vector along the genomic distance, and the attenuation vector is used to normalize the Hi-C contact matrix by genomic distance; then, the Pearson correlation coefficient of the two column vectors of the normalized matrix is ​​calculated to form a correlation coefficient matrix; the correlation coefficient matrix is ​​subjected to principal component analysis, and the principal component of each chromatin segment is extracted as the structural features of the A / B compartment; The method for extracting DNA sequence features is as follows: using the K-mer feature extraction method to divide the DNA sequence into lengths K, counting the frequencies of all subsequences of length K in the chromatin segment; sorting the frequencies to form a frequency vector; adding the frequency vectors, and selecting the feature with the highest frequency as the DNA sequence feature.

6. The TADs domain prediction method based on attribute graph node clustering according to claim 1, characterized in that: Inserting the attribute value of the node as a new node into the attribute graph model, and performing node expansion on the attribute graph to form an extended attribute graph, specifically includes: Normalize all node attribute values ​​to between 0 and 1, and discretize and sample the normalized node attribute values ​​according to the set step size; An attribute node is created for each sampled discrete value, the attribute node is inserted into the attribute graph model, attribute edges are inserted based on the attribute nodes, and an extended attribute graph is formed according to the added attribute nodes and attribute edges; the attribute edge is used to connect an attribute node and a structure node, and the structure node refers to an original node in the attribute graph model.

7. The TADs domain prediction method based on attribute graph node clustering according to claim 1, characterized in that: Clustering nodes in the extended attribute graph includes: Initialize edge weights; the edge weights are obtained according to the interaction frequency between nodes; Calculate the random walk distance matrix of nodes based on the attribute graph structure and edge weights; According to the density of nodes, the k structural nodes with the highest node density are determined as the starting clustering centers; Iteratively update the cluster to which each node belongs, the cluster center, the weight of the attribute edge, and the random walk distance matrix until convergence.

8. The TADs domain prediction method based on attribute graph node clustering according to claim 1, characterized in that: According to the continuity of the preliminary TADs domain and the minimum allowed TADs domain size, the TADs domain is further processed, including: If the nodes in the same cluster are continuously distributed on the genome, and the length of the continuous fragments on the genome is not less than the set minimum size of the TADs domain, the initial TADs domain will be retained; If nodes in the same cluster are isolated by other clusters on the genome and present segmented distribution, and the length of the continuous segments on the genome is not less than the set minimum size of the TADs domain, the separated cluster will be split into two non-adjacent TADs domains; The nodes in the same cluster are isolated by other clusters on the genome and present a segmented distribution. If there is a continuous fragment whose length is less than the set minimum size of the TADs domain, and the two sides of the continuous fragment less than the set minimum size of the TADs domain belong to the same cluster, then the fragment will be merged into the clusters on both sides to form a larger TADs domain. The nodes in the same cluster are isolated by other clusters on the genome and present a segmented distribution. If there is a continuous fragment whose length is less than the set minimum size of the TADs domain, and the two sides of the continuous fragment less than the set minimum size of the TADs domain belong to different clusters, the average random walk distance between the fragment and the clusters on both sides is evaluated, and it is merged into the cluster with a larger random walk distance to form a larger TADs domain.

9. A TADs domain prediction system based on attribute graph node clustering, used to implement the TADs domain prediction method based on attribute graph node clustering according to any one of claims 1 to 8, characterized in that: include: Hi-C contact matrix extraction module, which is used to compare paired-end sequencing sequences with the reference genome, filter low-quality sequences and PCR duplications, count DNA fragments, and generate Hi-C contact matrices; Data preprocessing module, used to preprocess Hi-C contact matrix; The attribute graph model generation module is used to map the chromatin fragments and their interaction information onto the attribute graph model; The data structure in the attribute graph model represents chromatin, and equal-length segments of chromatin are used as nodes of the attribute graph, the interactions between nodes are used as edges between nodes, and the A / B compartment structural features of the chromatin segment and the DNA sequence features are spliced ​​as the attributes of the node; the attribute value of the node is inserted into the attribute graph model as a new node, and the node of the attribute graph is amplified to form an extended attribute graph; The node clustering module is used to cluster densely connected structural nodes into the same category according to the structure of the generated attribute graph model to obtain preliminary TADs domains; based on the continuity of the preliminary TADs domains and the minimum allowed TADs domain size, the TADs domains are further processed to obtain the final predicted TADs domains.

10. The TADs domain prediction system based on attribute graph node clustering according to claim 9, characterized in that: It also includes a TADs domain post-processing and visualization module, which is used to further process the TADs domain according to the biological priors of TADs continuity and the minimum size of TADs, and to visualize the predicted TADs domain.

Images