ex_bf_by_residue

ex_bf_by_residue reads a PDB file and prints per-residue statistics of B-factors. The output provides: amino acid type (1-letter code), residue ID, and minimum, average and maximum b-factors for that residue

USAGE:
ex_bf_by_residue 2gb1.pdb

Keywords:

Categories:

  • core::data::io::Pdb

Input files:

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
#include <core/data/io/Pdb.hh>
#include <core/data/structural/structure_selectors.hh>
#include <utils/string_utils.hh>
#include <utils/exit.hh>

std::string program_info = R"(

ex_bf_by_residue reads a PDB file and prints per-residue statistics of B-factors. The output provides:
amino acid type (1-letter code), residue ID, and minimum, average and maximum b-factors for that residue
USAGE:
    ex_bf_by_residue 2gb1.pdb

)";


/** @brief Reads a PDB file and per-residue statistics of B-factors
 *
 * CATEGORIES: core::data::io::Pdb;
 * KEYWORDS:   PDB input; B-factors; structure selectors
 * IMG: Bfactor_plot.png
 * IMG_ALT: B-factors of 2GB1 PDB deposit
 */
int main(const int argc, const char *argv[]) {

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

  core::data::io::Pdb reader(argv[1]); // file name (PDB format, may be gzip-ped)
  core::data::structural::Structure_SP strctr = reader.create_structure(0);
  core::data::structural::IsBB is_bb;
  core::data::structural::IsAA is_aa;
  for (auto res_it = strctr->first_residue(); res_it != strctr->last_residue(); ++res_it) {
    if (!is_aa(**res_it)) continue;
    double min = 9999.0, max = -999.0, avg = 0.0, n = 0.0;
    for (auto atom : **res_it) {
//      if (is_bb(*atom)) continue; // --- uncommment that line to compute statistics for side chain only
      core::real bf = atom->b_factor();
      if (min > bf) min = bf;
      if (max < bf) max = bf;
      avg += bf;
      n += 1.0;

    }
    std::cout << (*res_it)->residue_type().code1 << " " <<
    utils::string_format("%4d %5.2f %5.2f %5.2f\n", (*res_it)->id(), min, avg / n, max);
  }
}
../_images/Bfactor_plot.png