Just because the algorithm from the previous section failed miserably does not mean that all brute force approaches are doomed. Can we design a faster brute force algorithm based on a different idea?

Instead of checking all cyclic peptides with a given mass, our new approach to solving the Cyclopeptide Sequencing Problem will "grow" candidate linear peptides whose theoretical spectra are “consistent” with the experimental spectrum.

STOP and Think What should it mean for a linear peptide to be consistent with an experimental spectrum? Would you classify VKF as consistent with the spectrum for Tyrocidine B1, reproduced below? What about VKY?

Tyrocidine b1 theoretical spectrum

Given an experimental spectrum Spectrum, we will form a collection Peptides of candidate linear peptides initially consisting of the empty peptide, which is just an empty string (denoted "") having mass 0. At the next step, we will expand Peptides to contain all linear peptides of length 1. We continue this process, creating 18 new peptides of length k + 1 for each amino acid string Peptide of length k in Peptides by appending every possible amino acid mass to the end of Peptide.

To prevent the number of candidate peptides from increasing exponentially, every time we expand Peptides, we trim it by keeping only those linear peptides that remain consistent with the experimental spectrum. We then check if any of these new linear peptides has mass equal to Mass(Spectrum). If so, we circularize this peptide and check whether it provides a solution to the Cyclopeptide Sequencing Problem.

More generally, brute force algorithms that enumerate all candidate solutions but discard large subsets of hopeless candidates by using various consistency conditions are known as branch-and-bound algorithms. Each such algorithm consists of a branching step to increase the number of candidate solutions, followed by a bounding step to remove hopeless candidates. In our branch-and-bound algorithm for the Cyclopeptide Sequencing Problem, the branching step will extend each candidate peptide of length k into 18 peptides of length k + 1, and the bounding step will remove inconsistent peptides from consideration.

Note that the spectrum of a linear peptide does not contain as many masses as the spectrum of a cyclic peptide with the same amino acid sequence. For instance, the theoretical spectrum of the cyclic peptide NQEL contains 14 masses (corresponding to "", N, Q, E, L, LN, NQ, QE, EL, ELN, LNQ, NQE, QEL, and NQEL). However, the theoretical spectrum of the linear peptide NQEL, shown below, does not contain masses corresponding to LN, LNQ, or ELN, since these subpeptides “wrap around” the end of the linear peptide.

NQEL linear spectrum

Exercise Break How many subpeptides does a linear peptide of given length n have? (Include the empty peptide and the entire peptide.) Solve at Cogniterra.

Input: An integer n.
Output: The number of subpeptides of a linear peptide of length n.

Given an experimental spectrum Spectrum of a cyclic peptide, a linear peptide is consistent with Spectrum if every mass in its theoretical spectrum is contained in Spectrum. If a mass appears more than once in the theoretical spectrum of the linear peptide, then it must appear at least that many times in Spectrum in order for the linear peptide to be consistent with Spectrum. For example, a linear peptide can still be consistent with the theoretical spectrum of NQEL if the peptide’s spectrum contains 242 twice. But it cannot be consistent with the theoretical spectrum of NQEL if its spectrum contains 113 twice.

The key to our new algorithm is that every linear subpeptide of a cyclic peptide Peptide is consistent with Cyclospectrum(Peptide). Thus, to solve the Cyclopeptide Sequencing Problem for Spectrum, we can safely ban all peptides that are inconsistent with Spectrum from the growing set Peptides, which powers the bounding step that we described above.

For example, the linear peptide VKF (with spectrum {0, 99, 128, 147, 227, 275, 374}) will be banned because it is inconsistent with Tyrocidine B1’s spectrum, reproduced below. But the linear peptide VKY will not be banned because every mass in its theoretical spectrum ({0, 99, 128, 163, 227, 291, 390}) is present in Tyrocidine B1's spectrum.

Tyrocidine b1 theoretical spectrum

What about the branching step? Given the current collection of linear peptides Peptides, define Expand(Peptides) as a new collection containing all possible extensions of peptides in Peptides by a single amino acid mass. We can now provide the pseudocode for the branch-and-bound algorithm, called CyclopeptideSequencing.

CyclopeptideSequencing(Spectrum)
  CandidatePeptides ← a set containing only the empty peptide
  FinalPeptides ← empty list of strings
  while CandidatePeptides is nonempty
    CandidatePeptides ← Expand(CandidatePeptides)
    for each peptide Peptide in CandidatePeptides
      if Mass(Peptide) = ParentMass(Spectrum)
        if Cyclospectrum(Peptide) = Spectrum and Peptide is not in FinalPeptides
          append Peptide to FinalPeptides
        remove Peptide from CandidatePeptides
      else if Peptide is not consistent with Spectrum
        remove Peptide from CandidatePeptides
  return FinalPeptides

Code Challenge Implement CyclopeptideSequencing (pseudocode reproduced below). Solve at Cogniterra or Rosalind.Solve it with the autograderFree to start in the interactive text

Note: After the failure of the first brute-force algorithm we considered, you may be hesitant to implement CyclopeptideSequencing for fear that its runtime will be prohibitive. The potential problem with CyclopeptideSequencing is that it may generate incorrect k-mers at intermediate stages (i.e., k-mers that are not subpeptides of a correct solution). In practice, however, this is not a concern. See "Charging Station: How fast is CyclopeptideSequencing?" in the print companion or at Cogniterra

CyclopeptideSequencing(Spectrum)
  CandidatePeptides ← a set containing only the empty peptide 
  FinalPeptides ← empty list of strings        
  while CandidatePeptides is nonempty
    CandidatePeptides ← Expand(CandidatePeptides)
    for each peptide Peptide in CandidatePeptides
      if Mass(Peptide) = ParentMass(Spectrum)
        if Cyclospectrum(Peptide) = Spectrum and Peptide is not in FinalPeptides
          append Peptide to FinalPeptides
        remove Peptide from CandidatePeptides
      else if Peptide is not consistent with Spectrum
        remove Peptide from CandidatePeptides
  return FinalPeptides

It is hard to imagine a worst-case scenario in which CyclopeptideSequencing takes a long time to run, but no one has been able to guarantee that this algorithm will not generate a huge number of incorrect k-mers at intermediate stages. BFCyclopeptideSequencing is exponential, and although in practice, CyclopeptideSequencing is much faster, this algorithm has not been proven to be polynomial. Thus, from the perspective of an introductory algorithms course focusing on theoretical computer science, the practical CyclopeptideSequencing algorithm is just as inefficient as BFCyclopeptideSequencing, since neither algorithm’s running time can be bounded by a polynomial.

A Branch-and-Bound Algorithm for Cyclopeptide Sequencing