GSMC: Combining Parallel Gibbs Sampling with Maximal Cliques for Hunting DNA Motif

Abstract

Regulatory elements are responsible for regulating gene transcription. Therefore, identification of these elements is a tremendous challenge in the field of gene expression. Transcription factors (TFs) play a key role in gene regulation by binding to target promoter sequences. A set of conserved sequence patterns with a highly similar structure that is bound by a TF is called a motif. Motif discovery has been a difficult problem over the past decades. Meanwhile, it is a foundation stone in meeting this challenge. Recent advances in obtaining genomic sequences and high-throughput gene expression analysis techniques have enabled the rapid development of computational methods for motif discovery. As a result, a large number of motif-finding algorithms aiming at various motif models have sprung up in the past few years. However, most of them are not suitable for analysis of the large data sets generated by next-generation sequencing.

To better handle large-scale ChIP-Seq data and achieve better performance in computational time and motif detection accuracy, we propose an excellent motif-finding algorithm known as GSMC (Combining Parallel Gibbs Sampling with Maximal Cliques for hunting DNA Motif). The GSMC algorithm consists of two steps. First, we employ the commonly used Gibbs sampling to generating initial motifs. Second, we utilize maximal cliques to cluster motifs according to Similarity with Position Information Contents (SPIC). Consequently, we raise the detection accuracy in a great degree, in the meantime holding comparative computation efficiency. In addition, we can find much more credible cofactor interacting motifs.

1. Introduction

An important goal of biological research to understand the mechanisms that regulate gene expression (Stegmaier et al., 2013). Identification of regulatory elements, especially the binding sites in DNA for transcription factors (TFs), is a vital task in this goal (Das and Dai, 2007). TFs are proteins that bind to DNA, typically upstream from and close to the transcription start site of a gene. The expression of that gene is regulated by activating or inhibiting the transcription mechanism (Tompa et al., 2005). Moreover, transcription regulation is usually triggered by the binding of TFs to specific DNA segments that are known as TF-binding sites (TFBSs) (Zhang and Chen, 2016). A motif is defined as a set of binding sites recognized by one TF. These binding sites are mutually similar, which happens to explain the general characteristics of the motif. Therefore, motif discovery is a significant step in the process of unraveling the mechanism of gene expression. At the most basic level, the motif-finding problem is the task of identifying recurring patterns of conserved short strings that appear in a set of DNA sequences (promoter region upstream from and close to the transcription start site of a gene). Usually, these patterns are fairly short (5 to 20 bps long).

A genome-wide ChIP study produces thousands or more DNA fragments consisting of several hundred base pairs, which cover the binding sites for a TF (Ikebata and Yoshida, 2015). Providing a set of DNA fragment sequences associated with a TF, by comparing the motifs discovered from these sequences with the known TF-binding motifs in a TFBSs database, for example, JASPAR (Sandelin et al., 2004), TRANSFAC (Wingender et al., 2000), we can not only recognize the binding sites for the ChIPed TF but also recognize the cofactor motifs that regulate the TF activity (Smith et al., 2005; Bailey, 2011; Goi et al., 2013).

The majority of early motif-finding algorithms can be divided into two categories. One is word-based (based on the string) methods that mostly adopt exhaustive enumeration, such as Weeder (Pavesi et al., 2001); the other is probabilistic-based (model-based) algorithms, for instance, MEME [multiple EM (expectation maximization) for motif elicitation] (Bailey and Elkan, 1994), AlignACE (Hughes et al., 2000), and ANN-Spec (Workman and Stormo, 2000). The word-based enumerative methods ensure global optimality, and they are appropriate for short motifs finding. Therefore, they are frequently used for motif finding in eukaryotic genomes where motifs are usually shorter than prokaryotes (Das and Dai, 2007). However, the methods cited earlier are not suitable for handling ChIP-seq data. Some of them have undergone reconstruction, arising from many ChIP-tailored algorithms. STEME (suffix free EM for motif elicitation) (Reid and Wernisch, 2011), a ChIP-tailored version of MEME, utilizes a branch-and-bound technology to remove negligible oligomers with significantly low probabilities effectively. For the sake of reducing the computational load in the counting operation, DREME (motif discovery in transcription factor—ChIP-seq data) (Bailey, 2011) and CisFinder (Sharov and Ko, 2009) adopt similar strategies, taking the risk of missing vital motifs in earlier steps of the recursion. Although Hegma (Ichinose et al., 2012) is competitive among current algorithms, the degradation of its detection accuracy is nonignorable.

The model-based approaches are mostly based on the EM (expectation maximization) algorithm (Bailey and Elkan, 1994) or Gibbs sampling (Lawrence et al., 1993). The RPMCMC (repulsive parallel MCMC algorithm for discovering diverse motifs) algorithm is a parallel version of the commonly used Gibbs sampling, running on a parallel-interacting Gibbs sampler. When the routes of different sampling chains are close to each other, a repulsive force (defined as a function of position probability matrices [PPMs]) will drag them forward in different directions. As a result, different sampling chains are urged to seek various regions so that the RPMCMC is more likely than old methods to discover much more diverse motifs. For old and other recent new motif discovery methods, we make a summary that they commonly exhibit low detection accuracy, are instable, and are time consuming. Thus, we are devoted to overcoming these drawbacks and improving them. As with the RPMCMC algorithm, we use parallel mutually exclusive Gibbs sampling to obtain the initial motif. When another motif similarity measure method and motif clustering method are applied to motif screening, we achieve the same or even better performance as the RPMCMC algorithm. These two methods will be demonstrated later.

The main purpose of this study is to derive a novel motif discovery algorithm, which acquires higher detection accuracy while holding competitive computational efficiency. Moreover, the novel method can detect much more believable diverse motifs. To cater to these needs, we propose a new motif-finding algorithm known as GSMC, based on Parallel Gibbs sampling (Lawrence et al., 1993) and Maximal Cliques (Zhang and Chen, 2016) clustering. We adopt parallel-interacting Gibbs sampling for generating initial motifs and then we utilize Maximal Cliques to group initial motifs into different clusters according to Similarity with Position Information Contents (SPIC) (Zhang et al., 2013) among motifs. Finally, we consider the first motif in each cluster as an output motif. On the basis of the ZOOPS (Zero or one motif occurrences per dataset sequence) sequence model, we design the GSMC algorithm without the OOPS constraint and motif length limit (Yu et al., 2015). We implement the GSMC algorithm with C++. When a set of TF ChIP-seq datasets of the ENCODE project (Dunham et al., 2012) is used to compare the performance of GSMC and two other high-performance algorithms, the results show that GSMC happens to make up for the other two algorithms on 20 ChIP-seq data sets for motif discovering; therefore, GSMC is a practical and alternative method in the field of motif finding.

2. Materials and Methods

2.1. Gibbs sampling model

The parallel Gibbs sampling algorithm applies several parallel-interacting Gibbs samplers to produce PPMs. Here, we provide an overview of the Gibbs sampling model.

Let \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $$S = \left( {{S_1} , \ldots , {S_N}} \right)$$ \end{document} denote a dataset of sequences, where N is equal to the number of sequences in the dataset. Each sequence is over the alphabet set \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $$\left\{ {A , C , G , T} \right\} $$ \end{document} . We look within each sequence for mutually similar segments of specified width W. The algorithm maintains two evolving data structures. The first is the pattern description, in the form of a probabilistic model of residue frequencies for each position i from 1 to W. Residue frequencies consist of the variables \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $${p_{i , a}} , {p_{i , c}} , {p_{i , g}} , {p_{i , t}}$$ \end{document} . This pattern description is accompanied by an analogous probabilistic description of the “background frequencies” \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $${q_a} , {q_c} , {q_g} , {q_t}$$ \end{document} and they denote the occurrence frequency of each residue in background sequences. The second data structure is a set of position \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $${u_k} \ \left( {k = 1 , \ldots , N} \right)$$ \end{document} , for the common pattern within the sequences. Given a set of input sequences, the objective is to identify the “best,” defined as the most probable, common pattern. This pattern is achieved by locating the alignment that maximizes the ratio of the corresponding pattern probability to background probability (Lawrence et al., 1993).

Initialization of the algorithm is to choose random starting positions in all sequences and to write them into starting positions set \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $$SS = \left\{ {s{s_1} , \ldots , s{s_N}} \right\} $$ \end{document} . It then proceeds via the iterative execution of the following two steps:

Step one: predictive update step. One of the N sequences, S_i is selected randomly or sequentially. The pattern description \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $${p_{i , j}}$$ \end{document} and background frequencies q_j are then calculated, as described in Equation (1), from the current position u_k in all sequences excluding S_i.

where B is the sum of the b_j. The q_j is calculated analogously, with the corresponding counts taken over all nonpattern positions.

Step two: sampling step. In the sequence S_i, each segment of width W is considered a possible instance of the pattern. \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $$Scor{e_e} \left( {i , j} \right)$$ \end{document} shows the possibility of the segment x from position j to \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $$j + W - 1$$ \end{document} in sequence S_i to be a pattern instance. \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $$Scor{e_b} \left( {i , j} \right)$$ \end{document} represents the possibility of the segment x from position j to \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $$j + W - 1$$ \end{document} in sequence S_i to be a background sequence. The weight \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $${r_{i , j}}$$ \end{document} is assigned to the segment x, and with each segment so weighted, a random one is selected. Its position then becomes the new \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $$s{s_i}$$ \end{document} . All the \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $${r_{i , j}}$$ \end{document} , \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $$Scor{e_e} \left( {i , j} \right)$$ \end{document} , and \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $$Scor{e_b} \left( {i , j} \right)$$ \end{document} of these segments are calculated according to the pattern description \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $${p_{i , j}}$$ \end{document} and background frequencies q_j. The calculation formulas of \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $${r_{i , j}}$$ \end{document} , \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $$Scor{e_e} \left( {i , j} \right)$$ \end{document} , and \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $$Scor{e_b} \left( {i , j} \right)$$ \end{document} are described by Equations (2), (3), and (4). \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} \begin{align*} Scor{e_e} \left( {i , j} \right) = Scor{e_e} \left( {{S_i} \left[ {j , j + W - 1} \right] } \right) = \prod \nolimits_{t = 0}^{W - 1} {{p_{t , {S_{i , t + j}}}}} \tag{2} \end{align*} \end{document} \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} \begin{align*} Scor{e_b} \left( {i , j} \right) = Scor{e_b} \left( {{S_i} \left[ {j , j + W - 1} \right] } \right) = \prod \nolimits_{t = 0}^{W - 1} {{q_{{S_{i , t + j}}}}} \tag{3} \end{align*} \end{document} \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} \begin{align*} { r_ { i , j } } = { \frac { Scor { e_e } \left( { i , j } \right) } { Scor { e_b } \left( { i , j } \right) } } \tag { 4 } \end{align*} \end{document}

After normalization, \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $${r_{i , j}}$$ \end{document} gives the probability that the pattern instance in sequence S_i belongs at position j. The algorithm finds the most probable alignment by selecting a set of u_i that maximizes the product of these ratios. Likewise, the algorithm finds the one that may maximize F, the sum of the logarithms of these ratios. In the notation developed earlier, F is indicated by Equation (5). \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} \begin{align*} F = \sum \nolimits_ { i = 1 } ^W { \sum \nolimits_ { j \in \left\{ { a , c , g , t } \right\} } { { c_ { i , j } } } } \log { \frac { { p_ { i , j } } } { { q_j } } } \tag { 5 } \end{align*} \end{document}

2.2. GSMC algorithm

The GSMC algorithm is composed of two parts, which include motif generating and postprocessing. GSMC shares some similarities with RPMCMC, especially in the phase of motif generating. First, we utilize parallel-interacting Gibbs sampling to produce various motifs. Second, based on SPIC and information content (IC), we employ maximal clique clustering to group motifs into different categories. Finally, we select the optimal motif in each category as the final output motif.

2.2.1. Motif generation

We use the ZOOPS model (Bailey and Elkan, 1994) in which there is zero or one motif occurring per sequence. As the input sequences, \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $${S^ + } = \left\{ {s_1^ + , \ldots , s_n^ + } \right\} $$ \end{document} and the length of them is \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $${L_i} \ \left( {i = 1 , \ldots n} \right)$$ \end{document} . \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $${S^ - } = \left\{ {s_1^ - , \ldots , s_n^ - } \right\} $$ \end{document} , which is the reverse complement of \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $${S^ + }$$ \end{document} . Our model uses the set of n concatenated sequences, \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $$S = \left\{ {{s_1} , \ldots , {s_n}} \right\} $$ \end{document} , where \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $${s_i} = \left( {s_i^ + , s_i^ - } \right) \ \left( {i = 1 , \ldots n} \right)$$ \end{document} . Whether there is a motif instance in sequence s_i or not, it depends on the value of z_i (one or zero). If \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $${z_i} = 1$$ \end{document} , a K-mer motif exists in sequence s_i, starting at the site \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $${u_i} \in \left\{ {1 , \ldots , {L_i} - K + 1 , {L_i} + 1 , \ldots , 2{L_i} - K + 1} \right\} $$ \end{document} . We utilize PPM to represent the motif. The size of PPM is \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $$4 \times K$$ \end{document} because the length of \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $$\left\{ {a , c , g , t} \right\} $$ \end{document} is four, where K denotes the length of motif. Let \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $${ \theta _k} = { \left( {{ \theta _{k , a}} , { \theta _{k , c}} , { \theta _{k , g}} , { \theta _{k , t}}} \right) ^T}$$ \end{document} indicate the kst position of motif and \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $${ \rm{ \Theta }} = \left( {{ \theta _1} , \ldots , { \theta _k}} \right)$$ \end{document} represent the motif. Let \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $${ \theta _0} = { \left( {{ \theta _{0 , a}} , { \theta _{0 , c}} , { \theta _{0 , g}} , { \theta _{0 , t}}} \right) ^T}$$ \end{document} be probability vector of the background sequences.

Provided the input sequences S, the objective is to detect the PPM \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $${ \rm{ \Theta }}$$ \end{document} without the knowledge of the motif length K and the background probability vector \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $${ \theta _0}$$ \end{document} where the latent variables contain \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $$U = \left\{ {{u_1} , \ldots , {u_n}} \right\} $$ \end{document} and \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $$Z = \left\{ {{z_1} , \ldots , {z_n}} \right\} $$ \end{document} .

As shown in Table 1 and Table 2, motif generation comprised two parts, which include initialization and updating. The stage of updating relies on Gibbs sampling.

Table 1.

Initialization Function for Initializing \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $$K , \theta , { \theta _0} , U , Z$$ \end{document}

Algorithm 1: Initialization function
Input: a set of sequences \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $$S , {k_{max}} , {k_{min}} , markov , n \_motif$$ \end{document}
Output: \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $$K , \theta , { \theta _0} , U , Z$$ \end{document}
1: \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $$K \leftarrow rand \left( \ {} \right) \% \ \left( {{k_{max}} - {k_{min}} + 1} \right) + {k_{min}}$$ \end{document}
2: for \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $$j \leftarrow 0$$ \end{document} to 3 do
3: \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $$tmp \left[ j \right] = gamma \_rand \left( {1 , \ 1} \right)$$ \end{document}
4: \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $$\theta \leftarrow tmp \left[ j \right] / \left( {tmp \left[ {0 \left] { + tmp} \right[ 1 \left] { + tmp} \right[ 2 \left] { + tmp} \right[ 3} \right] } \right)$$ \end{document}
5: \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $${ \theta _0} \leftarrow pow \left( {0.25 , \ markov} \right)$$ \end{document}
6: for \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $$p \leftarrow 0$$ \end{document} to the size of S do
7: for \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $$i \leftarrow 0$$ \end{document} to \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $$n \_motif$$ \end{document} do
8: \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $$U \leftarrow rand \left(\ {} \right) \% \ \left( {S \left[ p \right] .length \left(\ {} \right) - K \left[ i \right] - 1} \right) + 1$$ \end{document}
9: \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $$Z \leftarrow 0$$ \end{document}
10: return \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $$K , \theta , { \theta _0} , U , Z$$ \end{document}

Table 2.

Updating Function for Sampling and Saving Motifs

Algorithm 2: Updating function
Input: \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $$K , \theta , { \theta _0} , U , Z , n \_motif , iter , burnin$$ \end{document}
Output: the set of PPMs \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $${ \rm{ \Theta }}$$ \end{document}
1: for \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $$j \leftarrow 0$$ \end{document} to iter do
2: for \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $$i \leftarrow 0$$ \end{document} to \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $$n \_motif$$ \end{document} do
3: \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $$update \left( {{ \theta _i}} \right)$$ \end{document}
4: \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $$update \left( {{ \theta _0}_i} \right)$$ \end{document}
5: \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $$update \left( {{K_i}} \right)$$ \end{document}
6: \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $$update \left( U \right)$$ \end{document}
7: \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $$update \left( Z \right)$$ \end{document}
8: if \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $$j \ge burnin$$ \end{document} then
9: add \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $$\theta$$ \end{document} to \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $${ \rm{ \Theta }}$$ \end{document}
10: return \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $${ \rm{ \Theta }}$$ \end{document}

PPMs, position probability matrices.

The number of replicas, \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $$n \_motif$$ \end{document} set by the user is equal to the number of Gibbs samplers. When the number of iterations comes to the value of burnin, a number of exclusive Gibbs samplers start to sample and save data. Once on the procedure running for iter times, the procedure of motif generation terminates.

2.2.2. Postprocessing: clustering motifs

After accomplishing the section of motif generation, we can attain massive motif PPMs. Nevertheless, most of them are either false or redundant motifs. It is essential to filter and cluster PPMs for picking out true motifs.

First of all, we get rid of many false motifs according to their IC value less than a cutoff value, which is given by Formula (6): \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} \begin{align*} IC = \sum \nolimits_ { i = 1 } ^K { \sum \nolimits_ { j \in \left\{ { a , c , g , t } \right\} } { { p_ { i , j } } \log { \frac { { p_ { i , j } } } { { q_j } } } } } \tag { 6 } \end{align*} \end{document}

Then, we use the SPIC metric to compute the similarity score between each pair of rest motifs. As is shown later, the SPIC metric is given by Formulas (7), (8), (9), and (10). For a motif M_x containing n_x sequences with length L_x, let \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $${F_x} = { \left( {{f_x} \left( {b , X} \right) } \right) _{4 \times {L_x}}}$$ \end{document} be its position frequency matrix (PFM) and P_x be its position-specific scoring matrix (PSSM) that is defined as \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} \begin{align*} { P_x } = { \left( { { P_x } \left( { b , X } \right) } \right) _ { 4 \times { L_x } } } = { \left( { \log { \frac { { p_x } \left( { b , X } \right) } { { q_x } \left( b \right) } } } \right) _ { 4 \times { L_x } } } \qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad(7), \end{align*} \end{document}

where \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $${f_x} \left( {b , X} \right)$$ \end{document} and \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $${P_x} \left( {b , X} \right)$$ \end{document} are the count and probability of base \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $${ \rm{b}} \in \left\{ {A , C , G , T} \right\} $$ \end{document} appearing at position X of M_x (i.e., column X of P_x), respectively, and \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $${q_x} \left( b \right)$$ \end{document} is the probability of base b appearing in the background sequences. A pseudo-count is added when computing these probabilities. The IC of column X of the PSSM P_x is defined as \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} \begin{align*} I \left( { X , { P_x } } \right) = \sum \nolimits_ { b \in \left\{ { A , C , G , T } \right\} } { { p_x } \left( { b , X } \right) { P_x } \left( { b , X } \right) } = \sum \nolimits_ { b \in \left\{ { A , C , G , T } \right\} } { { p_x } \left( { b , X } \right) \log { \frac { { p_x } \left( { b , X } \right) } { { q_x } \left( b \right) } } } \qquad\qquad\qquad\qquad(8). \end{align*} \end{document}

For two motifs M₁ and M₂ with PSSMs P₁ and P₂, and PFMs F₁ and F₂, respectively, the similarity score between position X of M₁ and position Y of M₂ is defined as \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} \begin{align*} Sim \left( { { M_1 } \left( X \right) , { M_2 } \left( Y \right) } \right) = min \left\{ { 1 , { \frac { max \left\{ { score \left( { { P_1 } \left( X \right) , { F_2 } \left( Y \right) } \right) , score \left( { { P_2 } \left( Y \right) , { F_1 } \left( X \right) } \right) } \right\} } { max \left\{ { score \left( { { P_1 } \left( X \right) , { F_1 } \left( X \right) } \right) , score \left( { { P_2 } \left( Y \right) , { F_2 } \left( Y \right) } \right) } \right\} } } } \right\} \qquad\qquad\qquad\qquad(9), \end{align*} \end{document}

Further, we construct a motif similarity graph on the basis of similarity score between each pair of motifs. In the graph, each node stands for a motif. Two nodes will be linked by an edge if and only if the weight value is greater than a preset threshold. The weight value of two connected motifs is the similarity score between them. More specifically, binding site motifs that belong to the same TF are more likely to form highly connected sub-graphs with high edge weights in the motif similarity graph than are those from different TFs or spurious motifs. Therefore, if one motif has a low similarity score with any additional motifs, there is a strong possibility that this motif is spurious. These motifs will be deleted from the graph.

Finally, we utilize Maximal Cliques to cluster motifs on the basis of the motif similarity graph. Clustering motifs via Maximal Cliques consists of four steps. A brief introduction of the four steps is given next.

Step 1: For each node, find a maximal clique associated with it, as depicted in Figure 1. First, sort the neighbor nodes of v in ascending order by the weights of their edges incident to v to get an array \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $$\left\{ {E , D , C , B , A} \right\} $$ \end{document} . Second, successively delete the nodes in the array as well as their incident edges until v and the remaining nodes have the same degree. Finally, for each deleted node in the reverse order, if the node connects with all the remaining nodes in the initial graph, add the node and incident edges in the initial graph.

FIG. 1.

An example of finding maximal clique associated with node v.

Step 2: Merge cliques into clusters. Initially, sort all the unassigned cliques in descending order by the sum of edge weights. Then, compare the first clique with any other unassigned clique (C_i and C_j), if the overlap ratio of the nodes in C_i appearing in C_j divided by C_i is no less than a preset threshold \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $${ \rm{ \alpha }}$$ \end{document} or if the overlap ratio of the nodes in C_i appearing in C_j divided by C_j is no less than a preset threshold \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $${ \rm{ \beta }}$$ \end{document} , merge C_j into C_i, and repeat the previous steps until all of the cliques in the queue are assigned.

Step 3: Delete redundant nodes. Compare the latter cliques with the former, and remove the nodes in the latter cliques appearing in the fore.

Step 4: Sort clusters. Sort all the clusters in a decreasing order of edge-weight sum. Choose the first motif in each cluster as the final output motif.

2.3. Performance assessment

We report the performance of three motif discovery algorithms on the dataset: 20 TF ChIP-seq datasets of the ENCODE project. ChIP-seq experiments produce a great amount of DNA segments, which contain many TFBSs bound to a certain TF. In addition to finding these TFBSs, we can discover some cofactor motifs that are involved in the regulatory module of the primary TF. The performance evaluation of the algorithm is outlined next.

First of all, to compare motif detection accuracy among three algorithms, each predicted motif (PPM) by them is matched to JASPAR CORE motifs by using the online TOMTOM tool (Gupta et al., 2007). For a given predicted PPM, TOMTOM outputs the matching scores to all annotated TFBSs (the name of TFBSs) in JASPAR with the statistical significance (E-values: the expected number of false positives in the matches). For each algorithm, a diversity of the discovered motifs is evaluated with the number of known motifs in JASPAR CORE that are matched significantly to the produced PPMs with the acceptable level of significance at E-value less than 0.05 (Ikebata and Yoshida, 2015). The less the E-value, the higher the detection accuracy. For advanced options of TOMTOM program, we choose Pearson correlation coefficient as the motif column comparison function and set the significance threshold at 0.05. Then, we compare the total running time between the RPMCMC algorithm and our GSMC algorithm. Finally, we make a comparison of the number of cofactor motifs among three algorithms.

2.4. Programs and parameter selection

The RPMCMC program used in the article was released on January 6, 2015 (http://daweb.ism.ac.jp/yoshidalab/motif). The DREME program is available online (http://meme-suite.org/tools/dreme). Two motif discovery algorithms (RPMCMC and GSMC) were compiled and installed on 64-bit Ubuntu 14.04 in VMware. The DREME program is directly utilized online. To ensure a comparison that is as fair as possible among the three motif discovery algorithms, the values of the adjustable parameters in RPMCMC and GSMC were selected so as to obtain the optimal performance respectively. For the RPMCMC algorithm, we keep a default value of some parameters and change the others. We modify the number of iterations to 70 and the number of replicas to 20. For the DREME algorithm, we keep the value of all parameters unchanged because it is available online from the site http://meme-suite.org/tools/dreme. For the GSMC algorithm, the number of iterations and replicas are the same as the RPMCMC's. As a section of the GSMC program, clustering motifs are involved in several parameters, such as motif similarity cutoff, IC threshold, \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $${ \rm{ \alpha }}$$ \end{document} , and \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $${ \rm{ \beta }}$$ \end{document} . In agreement with Zhang and Chen (2016), the values of parameters \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $${ \rm{ \alpha }}$$ \end{document} and \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $${ \rm{ \beta }}$$ \end{document} are set to (0.5 0.5) and the value of the similarity cutoff parameter is set to 0.6. As the third section describes, the IC threshold parameter is sampled from 6 to 8.5 in steps of 0.5 and we choose 7 as its best value. All parameters of GSMC and RPMCMC are given in Table 3.

Table 3.

Default Parameters of GSMC and RPMCMC That Were Used in All Experiments, DREME Was Executed Using the Default Settings

GSMC
Parameter	Value
Max/min motif length	\documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $${L_{min}} = 6 , \ {L_{max}} = 14$$ \end{document}
No. of replicas	\documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $$M = 20$$ \end{document}
No. of iterations	\documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $${ \rm{{ \rm N}}} = 70$$ \end{document}
Burn-in period	\documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $${ \rm{Burning}} = 20$$ \end{document}
IC threshold of motif	\documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $$I{C_{cutoff}} = 7$$ \end{document}
Similarity threshold of motif	\documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $$Similarit{y_{cutoff}} = 0.6$$ \end{document}

RPMCMC
Parameter	Value
Max/min motif length	\documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $${L_{min}} = 6 , \ {L_{max}} = 14$$ \end{document}
No. of replicas	\documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $$M = 20$$ \end{document}
No. of iterations	\documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $${ \rm{{ \rm N}}} = 70$$ \end{document}
Burn-in period	\documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $${ \rm{Burning}} = 20$$ \end{document}
Motif clustering	\documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $${ \rm{ \lambda }} = 0.3$$ \end{document}

DREME; GSMC, Combining Parallel Gibbs Sampling with Maximal Cliques for hunting DNA Motif; IC, information content; RPMCMC.

3. Results and Discussion

3.1. Parameter optimization and performance comparison

Empirically, we set the value of the similarity cutoff parameter to 0.6 and the values of parameters \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $${ \rm{ \alpha }}$$ \end{document} and \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $${ \rm{ \beta }}$$ \end{document} to (0.5 0.5). For the sake of optimizing the value of the IC threshold parameter, we choose one TF ChIP-seq dataset as the parameter test dataset and make 16 times experiments of the dataset for each IC threshold value from 6 to 8.5 in steps of 0.5 whereas the remaining values of other parameters are unchanged. Motif E-values of different IC thresholds are presented in Figure 2A. Average E-values of discrete IC thresholds are given in Figure 2B. Figure 2C describes the E-values' mean-square deviation of different IC thresholds. The three figures are indicated next.

FIG. 2.

Optimal value selection of the parameter IC threshold. (A) The IC thresholds are selected from 6 to 8.5 in steps of 0.5. For each IC threshold value, a ChIP dataset (wgEncodeSydhTfbsHepg2Nrf1IggrabPk) is tested 16 times. Every time, we choose the minimal E-value of all meeting matches between the predicted motifs and the motif NRF1. Then, we draw them in broken lines. (B) For each IC threshold value, we calculate the average E-value of 16 minimal E-values and draw these average E-values in broken lines. (C) For each IC threshold value, we compute the mean-square deviation of the 16 minimum E-values and plot them in the same way. IC, information content.

According to Figure 2B, we can refer to the relationship between the average detection accuracy and the IC threshold value. The less the average E-value, the higher the average detection accuracy. From Figure 2C, we know the fluctuations of the 16 minimal E-values. The less the E-values' mean-square deviation, the higher the stability. Above all, 7 is chosen as the optimal value of the parameter IC threshold in that there is the average least E-value and the second least mean-square deviation in this case.

3.2. Performance on discovering motifs for ENCODE ChIP-Seq datasets

We use the same ENCODE datasets as the paper (Ikebata and Yoshida, 2015). They are SYDH TFBS narrowPeak files (available from NCBI's Gene Expression Omnibus using the accession number GSE31477), and we can download them from http://hgdownload.cse.ucsc.edu goldenPath/hg19/encodeDCC/wgEncodeSydhTfbs/. We choose 20 datasets with different sizes from the 228 ENCODE ChIP-seq datasets as our experimental data. GSMC is implemented in C++, RPMCMC is available on the authors' websites, and DREME can be used online directly. All the tests are conducted on Intel^® Core™ i3 processor with 4-core CPUs and 4 GB of main memory. On the virtual platform of Ubuntu 14.04, we run GSMC and RPMCMC for 20 datasets, record running time, and save all the produced PPMs. With the help of the TOMTOM, we compare the produced PPMs with all the known motifs in JASPAR and leave the least E-value of matching between the ChIPed motif and the query motifs. Then, we make comparisons among GSMC, RPMCMC, and DREME in terms of detection accuracy, computational speeds, and the number of cofactor motifs discovered. From a performance perspective, our GSMC algorithm makes up for rather than takes the place of existing DNA motif-finding methods for analyzing the large ChIP-seq data.

3.2.1. Detection accuracy

Each dataset corresponds to a ChIPed TF, which has a relevant motif. We run each dataset for 50 times and each time, we make a note of the least E-value of all meeting matches (E-value <0.05) that are aimed at the ChIPed motif. Figure 3 illustrates the comparison of detection accuracy among three algorithms on a ChIP dataset (wgEncodeSydhTfbsHepg2Nrf1IggrabPk) in which the binding sites of NRF1 were studied in cell HepG2. When it comes to detection accuracy, DREME surpasses the others and GSMC is nearly equal to RPMCMC.

FIG. 3.

Comparison of motifs' detection accuracy among RPMCMC, GSMC, and DREME on datasets (wgEncodeSydhTfbsHepg2Nrf1IggrabPk). DREME; GSMC, Combining Parallel Gibbs Sampling with Maximal Cliques for hunting DNA Motif; RPMCMC.

In Figure 3, we have made a comparison of detection accuracy among three algorithms on one dataset. Next, we make comparisons on all chosen datasets. Figure 4 provides the comparison of GSMC with RPMCMC and DREME on 20 chosen datasets. For most datasets, GSMC and RPMCMC are superior in detection accuracy to DREME.

FIG. 4.

Comparison of motifs' detection accuracy among RPMCMC, GSMC, and DREME on 20 chosen datasets.

3.2.2. Computational speeds

Figure 5 displays the computational time for GSMC and RPMCMC on 20 datasets. Each dataset is tested for 16 times, and we take the arithmetic mean as the running time. In terms of computation efficiency, GSMC is comparable to RPMCMC. Motif producing of GSMC basically conforms to RPMCMC. Therefore, the bottleneck calculating the posterior probabilities of the motif start sites u_i in RPMCMC still exists in GSMC. However, the clustering section of GSMC is superior to RPMCMC in terms of flexibility and scalability. With these performances, GSMC is more likely to be a high-efficient method in future studies. Though GSMC cannot surpass RPMCMC in computational efficiency, it is a competitive algorithm in the case of computing time.

FIG. 5.

Comparison of computational efficiency between GSMC and RPMCMC. Running time is the CPU execution time.

3.2.3. Comparison on the number of cofactor motifs discovered

A recent study by Chin Lui Goi and Peter Little called cell-type and TF-specific enrichment of transcriptional cofactor motifs in ENCODE ChIP-seq data claims that TF-specific interactions between TFs and cofactors are essential for transcriptional regulation through recruitment of general transcription machinery to gene promoter regions (Goi et al., 2013). It is also mentioned in the study that some of the cofactor motifs are experimentally verified cofactors and others are potentially novel cofactors. The same can be said of the cofactors discovered in our study. Even with all that, they have evident reference value and directive significance to further studies. Therefore, the ability of detecting cofactor motifs is crucial to a motif-discovering algorithm. Just as shown in Figure 6, the cofactor motifs found by GSMC outnumber those found by RPMCMC on almost all datasets.

FIG. 6.

Comparison on the number of cofactor motifs. Each dataset is tested 16 times, and each data point in the figure is the mean value of 16 tests.

Obviously, GSMC is superior in the number of cofactor motifs to RPMCMC. In other words, the ability of GSMC to detect cofactor motifs is better than RPMCMC and basically equal to DREME.

4. Conclusion

Many older popular motif finders cannot be applied to handling large datasets generated by ChIP-seq (Reid and Wernisch, 2014). Though a vast amount of novel algorithms spring up for handling huge datasets, most of them attend to computation speed and lose sight of the accuracy of motif detection. To overcome the drawbacks of them, we develop a ChIP-tailored motif discovery tool called GSMC, which hunts for DNA motifs by combining Parallel Gibbs Sampling with Maximal Cliques. In terms of computation time and detection accuracy, GSMC can rival the most recent motif discoverer RPMCMC, which was specifically designed to handle large datasets. Moreover, GSMC is capable of detecting much more known and potential cofactor motifs where GSMC far exceeds RPMCMC. Besides, GSMC is superior in accuracy of motif detection to DREME. On the whole, we present a novel motif finder reconciling computation speed with detection ability. Locating these motifs plays a critical role in identifying transcriptional regulation. We expect that GSMC will have a place in the process of exploring gene expression.

Footnotes

Acknowledgments

This work was supported by the grants of the National Science Foundation of China (Grant Nos. 61472467, 61672011, and 61471169) and the Collaboration and Innovation Center for Digital Chinese Medicine of 2011 Project of Colleges and Universities in Hunan Province.

Author Disclosure Statement

No competing financial interests exist.

References

Bailey

T.L.

2011. DREME: Motif discovery in transcription factor ChIP-seq data. Bioinformatics, 27, 1653–1659.

Bailey

T.L.

, and Elkan

1994. Fitting a mixture model by expectation maximization to discover motifs in bipolymers. In Proceedings of the Second International Conference on Intelligent Systems for Molecular Biology, 28–36.

Das

M.K.

, and Dai

H.K.

2007. A survey of DNA motif finding algorithms. BMC Bioinformatics, 8, S21.

Dunham

, Kundaje

, Aldred

S.F.

, et al. 2012. An integrated encyclopedia of DNA elements in the human genome. Nature, 489, 57–74.

Goi

C.L.

, Little

, and Xie

2013. Cell-type and transcription factor specific enrichment of transcriptional cofactor motifs in ENCODE ChIP-seq data. BMC Genomics, 14, 1–11.

Gupta

, Stamatoyannopoulos

J.A.

, Bailey

T.L.

, et al. 2007. Quantifying similarity between motifs. Genome Biol. 8, R24.

Hughes

J.D.

, Estep

P.W.

, Tavazoie

, et al. 2000. Computational identification of cis-regulatory elements associated with groups of functionally related genes in Saccharomyces cerevisiae. J. Mol. Biol. 296, 1205–1214.

Ichinose

, Yada

, and Gotoh

2012. Large-scale motif discovery using DNA Gray code and equiprobable oligomers. Bioinformatics, 28, 25–31.

Ikebata

, and Yoshida

2015. Repulsive parallel MCMC algorithm for discovering diverse motifs from large sequence sets. Bioinformatics, 31, 1561–1568.

10.

Lawrence

C.E.

, Altschul

S.F.

, Boguski

M.S.

, et al. 1993. Detecting subtle sequence signals: A Gibbs sampling strategy for multiple alignment. Science, 262, 208–214.

11.

Pavesi

, Mauri

, and Pesole

2001. An algorithm for finding signals of unknown length in DNA sequences. Bioinformatics, 17, S207–S214.

12.

Reid

J.E.

, and Wernisch

2011. STEME: Efficient EM to find motifs in large data sets. Nucleic Acids Res. 39, e126.

13.

Reid

J.E.

, and Wernisch

2014. STEME: A robust, accurate motif finder for large data sets. PLoS One, 9, e90735.

14.

Sandelin

, Alkema

, Engstrom

, et al. 2004. JASPAR: An open-access database for eukaryotic transcription factor binding profiles. Nucleic Acids Res. 32, D91–D94.

15.

Sharov

A.A.

, and Ko

M.S.

2009. Exhaustive search for over-represented DNA sequence motifs with CisFinder. DNA Res. 16, 261–273.

16.

Smith

A.D.

, Sumazin

, Das

, et al. 2005. Mining ChIP-chip data for transcription factor and cofactor binding sites. Bioinformatics, 21, I403–I412.

17.

Stegmaier

, Kel

, Wingender

, et al. 2013. A discriminative approach for unsupervised clustering of DNA sequence motifs. PLoS Comput.Biol., 9, e1002958.

18.

Tompa

, Li

, Bailey

T.L.

, et al. 2005. Assessing computational tools for the discovery of transcription factor binding sites. Nat. Biotechnol. 23, 137–144.

19.

Wingender

, Chen

, Hehl

, et al. 2000. TRANSFAC: An integrated system for gene expression regulation. Nucleic Acids Res. 28, 316–319.

20.

Workman

, and Stormo

2000. ANN-Spec: A method for discovering transcription factor binding sites with improved specificity. Pac. Symp. Biocomput. 5, 467–478.

21.

, Huo

H.W.

, Chen

X.Y.

, et al. 2015. An efficient algorithm for discovering motifs in large DNA data sets. IEEE Trans. Nanobiosci., 14, 535–544.

22.

Zhang

S.Q.

, and Chen

2016. CLIMP: Clustering motifs via maximal cliques with parallel computing design. PLoS One, 11, e0160435.

23.

Zhang

S.Q.

, Zhou

X.G.

, Du

C.B.

, et al. 2013. SPIC: A novel similarity metric for comparing transcription factor binding site motifs based on information contents. BMC Syst. Biol., 7 Suppl 2, S14.