Skip to content
dstoeckel edited this page Mar 16, 2015 · 2 revisions

How can I mutate amino acids of a given protein or peptide?

BALLs FragmentDB allows to mutate a residue while keeping the backbone atom positions.

The following example first creates a small peptide (reading a protein from a pdb file is possible as well), mutates two amino acids of the given peptide, and finally writes the mutated peptide into a PDB file.

Please note :

  • To apply the FragmentDB, the program deletes all side chain atoms.
  • The call to setName() expects ThreeLetterCode, otherwise it fails.
  • This example does not change the atom positions, the side chain atoms are added on reasonable positions, but the backbone remains. Any positioning should be done afterwards (see TranslateAMolecule).
#include <BALL/FORMAT/PDBFile.h>
#include <BALL/STRUCTURE/fragmentDB.h>
#include <BALL/STRUCTURE/residueChecker.h>
#include <BALL/STRUCTURE/peptideBuilder.h>
#include <vector>

using namespace BALL;
using namespace std;

int main(int argc, char** argv)
{
	String filename_ori = "ori.pdb";
	String filename_mut = "mut.pdb";

	FragmentDB fdb = FragmentDB("");
	
	// Create a peptide
	Peptides::PeptideBuilder* pb = new BALL::Peptides::PeptideBuilder;
	pb->setFragmentDB(&fdb);
        pb->addAminoAcid("SER",Angle(M_PI, true),Angle(M_PI, true));
	pb->addAminoAcid("CYS",Angle(M_PI, true),Angle(M_PI, true));
	pb->addAminoAcid("PRO",Angle(M_PI, true),Angle(M_PI, true));
	pb->addAminoAcid("ARG",Angle(M_PI, true),Angle(M_PI, true));

	//Peptides::PeptideBuilder pb("SCPR");

        pb->setChainName("TestChain");
        pb->setProteinName("TestProtein");

        Protein* prot = pb->construct();

        System system;
	system.insert(*prot);

	PDBFile outfile_ori(filename_ori, std::ios::out);
	outfile_ori << system;

	String old_seq = Peptides::GetSequence(*prot);
	cout << "ori sequence: " << old_seq << endl;
	cout << "#residues/atoms: " << prot->countResidues() << " / " << prot->countAtoms() << endl;

	// now create a mutated version: mutate to SPCR
	String mut_seq = "SPCR";
	cout << "mutate sequence to " << mut_seq << endl;

	Size seq_index = 0;
	if ( old_seq.size() != mut_seq.size() )
	{
		// one may have to delete or add residues...
		cout << "ERROR: " << old_seq.size() << " != " <<  mut_seq.size() << endl;
		return 1;
	}

	// collect while iterating over all residues
	vector<pair<Residue*, char> > to_mutate;

	for (ResidueIterator rit = prot->beginResidue(); +rit; ++rit)
	{
		if (Peptides::OneLetterCode(rit->getName()) != mut_seq[seq_index])
		{
			cout << "  mutation at position " << seq_index << " from " << Peptides::OneLetterCode(rit->getName())
					 << " to " << mut_seq[seq_index] << endl;

			to_mutate.push_back(pair<Residue*, char>(&*rit, mut_seq[seq_index]));
		}
		seq_index++;
	}

	AtomIterator at_it_ori;

	//    mutate the residues
	for (Size j=0; j<to_mutate.size(); ++j)
	{
		vector<Atom*> to_delete;

		// delete all side chain atoms
		for ( at_it_ori = to_mutate[j].first->beginAtom();
				 +at_it_ori;
				++at_it_ori)
		{
			// is it a backbone atom? Keep them! 
			if ( !( (at_it_ori->getName()== "C") ||
						(at_it_ori->getName()== "O") ||
						(at_it_ori->getName()== "N") ||
						(at_it_ori->getName()== "HA")
						)
				 )
			{
				to_delete.push_back(&*at_it_ori);
			}
		}
		cout << "delete " << to_delete.size() << " side chain atoms ";

		// really delete them
		for (Size k=0; k<to_delete.size(); ++k)
		{
			to_delete[k]->destroy();
		}

		cout << "and mutate " << to_mutate[j].first->getName() << " to ";

		// rename the residue for the fragmentDB
		to_mutate[j].first->setName(Peptides::ThreeLetterCode(to_mutate[j].second));

		cout << to_mutate[j].first->getName() << " ..." << endl;

		// apply the fragmentDB
		to_mutate[j].first->apply(fdb.normalize_names);
		to_mutate[j].first->apply(fdb.add_hydrogens);
		to_mutate[j].first->apply(fdb.build_bonds);

		to_delete.clear();
	} // end of mutate residues


	cout << " done." <<  endl;
	cout << "mutated sequence: " << Peptides::GetSequence(*prot) << endl;
	cout << "mutated #residues/atoms: " << prot->countResidues() << " / " << prot->countAtoms() << endl;

	system.insert(*prot);

	// write the result	
	PDBFile outfile(filename_mut, std::ios::out);
	outfile << system;

	cout << "Wrote files " << filename_ori << " and " << filename_mut << endl;
	cout << " done." <<  endl;

}
Clone this wiki locally