Lesson 2.10

Gibbs Sampling in Action

4 min read✓ Read

Watch the lectureHow Rolling Dice Helps Us Find Regulatory Motifs (Part 3)

We will illustrate how GibbsSampler works on the same strings Dna that we considered before. Imagine that, at the initial step, the algorithm has chosen the following 4-mers (shown in red) and has randomly selected the third string for removal. To be more precise, GibbsSampler does not really remove the third string; it ignores it at this particular step and may analyze it again in subsequent steps.

This results in the following motif, count, and profile matrices.

Note that the profile matrix deviates only slightly from the uniform profile, making us wonder whether we have any chance to be steered toward the implanted motif. We now use this profile matrix to compute the probabilities of all 4-mers in the deleted string ccgGCGTtag:

Note that all but two of these probabilities are zero. This situation is similar to the one we encountered with GreedyMotifSearch, and as before, we need to augment zero probabilities with small pseudocounts to avoid disastrous results.

Application of Laplace’s Rule of Succession to the count matrix above yields the following updated count and profile matrices:

After adding pseudocounts, the 4-mer probabilities in the deleted string "ccgGCGTtag" are recomputed as follows:

Since these probabilities sum to C = 80/8⁴, our hypothetical seven-sided die is represented by the random number generator

Random(4/8⁴, 8/8⁴, 8/8⁴, 24/8⁴, 12/8⁴, 16/8⁴, 8/8⁴) =
Random((4/8⁴)/(80/8⁴), (8/8⁴)/(80/8⁴), (8/8⁴)/(80/8⁴), (24/8⁴)/(80/8⁴), (12/8⁴)/(80/8⁴), (16/8⁴)/(80/8⁴), (8/8⁴)/(80/8⁴)) =
Random (4/80, 8/80, 8/80, 24/80, 12/80, 16/80, 8/80)

Let's assume that after rolling this seven-sided die, we arrive at the Profile-randomly generated 4-mer GCGT (the fourth 4-mer in the deleted sequence). The deleted string ccgGCGTtag is now added back to the collection of motifs, and GCGT substitutes the previously chosen ccgG in the third string in Dna, as shown below. We then roll a fair five-sided die and randomly select the first string from Dna for removal.

After constructing the motif and profile matrices, we obtain the following:

Note that the profile matrix looks more biased toward the implanted motif than the previous profile matrix did.

We update the count and profile matrices with pseudocounts:

Then, we compute the probabilities of all 4-mers in the deleted string "ttACCTtaac":

When we roll a seven-sided die, we arrive at the Profile-randomly generated k-mer ACCT, which we add to the collection Motifs.

After rolling the five-sided die once again, we randomly select the fourth string for removal.

We further add pseudocounts and construct the resulting count and profile matrices:

We now compute the probabilities of all 4-mers in the deleted string "cactaACGAg":

We need to roll a seven-sided die to produce a Profile-randomly generated 4-mer. Assuming the most probable scenario in which we select ACGA, we update the selected 4-mers as follows:

You can see that the algorithm is beginning to converge. Rest assured that a subsequent iteration will produce all implanted motifs after we select the second string in Dna (when the incorrect 4-mer GTct will likely change into the implanted (4, 1)-motif ACGT).

Try it Watch RandomizedMotifSearch and GibbsSampler at work

Now make the dice roll in your own code: Implement RandomizedMotifSearch and GibbsSampler in the interactive text.

STOP and Think Run GibbsSampler on the Subtle Motif Problem. What do you find?

Although GibbsSampler performs well in many cases, it may converge to a suboptimal solution, particularly for difficult search problems with elusive motifs. A local optimum is a solution that is optimal within a small neighboring set of solutions, which is in contrast to a global optimum, or the optimal solution among all possible solutions. Since GibbsSampler explores just a small subset of solutions, it may "get stuck" in a local optimum. For this reason, similarly to RandomizedMotifSearch, it should be run many times with the hope that one of these runs will produce the best-scoring motifs. Yet convergence to a local optimum is just one of many issues we must consider in motif finding; for more details, see "DETOUR: Complications of Motif Finding" in the print companion or at Cogniterra.

When we run GibbsSampler 2,000 times on the Subtle Motif Problem with implanted 15-mer AAAAAAAAGGGGGGG (each time with new randomly selected k-mers for N = 200 iterations), it returns a collection Motifs with consensus string AAAAAAgAGGGGGGt and Score(Motifs) equal to 38. This score is even lower than the score of 40 expected from the implanted motifs!