ap_docking_crmsd

ap_docking_crmsd calculates crmsd between ligand positions after flexible docking to a receptor and a reference.

The program reads in a native pose and at least one PDB file with a computed pose (i.e. a model), each of them must contain a ligand molecule bound to a protein receptor. The ligand can be a small molecule, peptide or even a protein. The program finds the ligand either by residue ID (a three-letter code, such as CAM) or a chain ID - a single letter code.

USAGE:
ap_docking_crmsd 2m56-ref.pdb CAM 00199.pdb 00963.pdb 04473.pdb
ap_docking_crmsd 2m56-ref.pdb X 00199.pdb 00963.pdb 04473.pdb
ap_docking_crmsd - X 00199.pdb 00963.pdb 04473.pdb

where 2m56-ref.pdb is the native and CAM is the three-letter PDB code of the ligand for which crmsd will be evaluated and 00199.pdb and the two other files are conformation after docking. In the second example, X is the ID of the chain containing a ligand molecule.

The program evaluates crmsd based on ligand cooridinates upon an optimal superimposition of a receptor between the refence and any model of any input file. If the reference structure is not given and dash ‘-‘ character is used instead (as in the last example), the program evaluates pairwise all-vs-all crmsd calculations

Categories:

  • core::protocols::PairwiseCrmsd

Output files:

Program source:

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
#include <iostream>
#include <map>

#include <core/data/io/Pdb.hh>
#include <core/protocols/PairwiseCrmsd.hh>
#include <core/data/structural/structure_selectors.hh>
#include <utils/exit.hh>

std::string program_info = R"(

ap_docking_crmsd calculates crmsd between ligand positions after flexible docking to a receptor and a reference.

The program reads in a native pose and at least one PDB file with a computed pose (i.e. a model), each of them must 
contain a ligand molecule bound to a protein receptor. The ligand can be a small molecule, peptide or even a protein.
The program finds the ligand either by residue ID (a three-letter code, such as CAM) or a chain ID - a single letter code.

USAGE:
    ap_docking_crmsd 2m56-ref.pdb CAM 00199.pdb 00963.pdb 04473.pdb
    ap_docking_crmsd 2m56-ref.pdb X 00199.pdb 00963.pdb 04473.pdb
    ap_docking_crmsd - X 00199.pdb 00963.pdb 04473.pdb

where  2m56-ref.pdb is the native and CAM is the three-letter PDB code of the ligand for which crmsd will be evaluated
and 00199.pdb and the two other files are conformation after docking. In the second example, X is the ID of the chain
containing a ligand molecule.

The program evaluates crmsd based on ligand cooridinates upon an optimal superimposition of a receptor between the refence 
and any model of any input file. If the reference structure is not given and dash '-' character is used instead 
(as in the last example), the program evaluates pairwise all-vs-all crmsd calculations

)";

using namespace core::data::structural;

/** @brief ap_docking_crmsd calculates crmsd between two ligand positions after docking to a receptor.
 *
 * CATEGORIES: core::protocols::PairwiseCrmsd
 * KEYWORDS: PDB input; crmsd; docking; structure selectors
 */
int main(const int argc, const char* argv[]) {

  using namespace core::data::structural;

  if (argc < 4) utils::exit_OK_with_message(program_info); // --- complain about missing program parameter

  const std::string code(argv[2]);    // --- The ligand code is the second parameter of the program
  AtomSelector_SP select_ligand = nullptr; // --- Ligand selector object
  if (code.size() == 3) // --- If the code is 3 characters long, its a residue code
    select_ligand = std::static_pointer_cast<AtomSelector>(std::make_shared<SelectResidueByName>(code));
  else {
    AtomSelector_SP select_chain = std::static_pointer_cast<AtomSelector>(std::make_shared<ChainSelector>(code[0]));
    std::shared_ptr<LogicalANDSelector> select_chain_ca = std::make_shared<LogicalANDSelector>();
    select_chain_ca->add_selector( std::make_shared<IsCA>() );
    select_chain_ca->add_selector(select_chain);
    select_ligand = select_chain_ca;
  }

  std::shared_ptr<LogicalANDSelector> select_receptor = std::make_shared<LogicalANDSelector>(); // --- Receptor selector object
  AtomSelector_SP select_not_ligand = std::static_pointer_cast<AtomSelector>(std::make_shared<InverseAtomSelector>(*select_ligand));
  select_receptor->add_selector(std::make_shared<IsCA>());
  select_receptor->add_selector(select_not_ligand);

  core::protocols::PairwiseCrmsd crmsd_protocol(select_receptor, select_ligand);
 
  for(core::index2 i=3;i<argc;++i) {
      core::data::io::Pdb reader(argv[i], core::data::io::is_not_hydrogen);
      if (reader.count_models()>1) {
        for (core::index2 j = 0; j < reader.count_models(); ++j)
          crmsd_protocol.add_input_structure(reader.create_structure(j), utils::string_format("%s:%4d",argv[i], j));
      } else 
          crmsd_protocol.add_input_structure(reader.create_structure(0), argv[i]);
  }

  crmsd_protocol.crmsd_cutoff(20.0); // crmsd cutoff large anough to get some output
  crmsd_protocol.output_stream( std::shared_ptr<std::ostream>(&std::cout, [](void*) {}) );
  if(std::string(argv[1]) != "-") {
    core::data::io::Pdb reader(argv[1], core::data::io::is_not_hydrogen);
    Structure_SP native = reader.create_structure(0);
    crmsd_protocol.calculate(native);
  } else crmsd_protocol.calculate();

}
../_images/file_icon.png