File: test_sidechain_construction.cc

package info (click to toggle)
promod3 3.4.2%2Bds-1
  • links: PTS, VCS
  • area: main
  • in suites: trixie
  • size: 966,596 kB
  • sloc: cpp: 55,820; python: 18,058; makefile: 85; sh: 51
file content (153 lines) | stat: -rw-r--r-- 5,038 bytes parent folder | download | duplicates (3)
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
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
// Copyright (c) 2013-2020, SIB - Swiss Institute of Bioinformatics and
//                          Biozentrum - University of Basel
// 
// Licensed under the Apache License, Version 2.0 (the "License");
// you may not use this file except in compliance with the License.
// You may obtain a copy of the License at
// 
//   http://www.apache.org/licenses/LICENSE-2.0
// 
// Unless required by applicable law or agreed to in writing, software
// distributed under the License is distributed on an "AS IS" BASIS,
// WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
// See the License for the specific language governing permissions and
// limitations under the License.



#include <promod3/loop/sidechain_atom_rule_lookup.hh>
#include <promod3/loop/sidechain_atom_constructor.hh>
#include <promod3/loop/all_atom_positions.hh>
#include <promod3/core/message.hh>
#define BOOST_TEST_DYN_LINK
#include <boost/test/unit_test.hpp>

#include <ost/mol/bond_handle.hh>
#include <ost/io/mol/pdb_reader.hh>
#include <ost/conop/heuristic.hh>

BOOST_AUTO_TEST_SUITE( loop );

using namespace promod3::loop;

namespace {

ost::mol::EntityHandle LoadTestStructure(const String& pdb_name) {
  ost::io::PDBReader reader(pdb_name, ost::io::IOProfile());
  ost::mol::EntityHandle test_ent = ost::mol::CreateEntity();
  reader.Import(test_ent);
  ost::conop::HeuristicProcessor heu_proc;
  heu_proc.Process(test_ent);

  return test_ent;
}

} // anon ns

BOOST_AUTO_TEST_CASE(test_construct_sidechain_atoms) {
  // get data
  // the idea of the whole test is to transform test_ent into an 
  // AllAtomPositions object, extract the sidechain dihedrals,
  // construct a second AllAtomPositions object, reconstruct
  // the sidechains with the extracted dihedral angles and compare.
  ost::mol::EntityHandle test_ent = LoadTestStructure("data/all_aa.pdb");
  ost::mol::ResidueHandleList res_list = test_ent.GetResidueList();
  AllAtomPositions initial_all_pos(res_list);

  // extract the dihedral angles
  std::vector<std::vector<Real> > dihedral_angles;
  for(uint i = 0; i < initial_all_pos.GetNumResidues(); ++i){
    std::vector<Real> current_angles;

    try{
      for(uint j = 0; j < 4; ++j){
        current_angles.push_back(GetChiAngle(initial_all_pos, i, j));
      }
    } catch(promod3::Error& e){
      //don't do anythin...
    }
    dihedral_angles.push_back(current_angles);
  }

  // generate a new AllAtomPositions object only containing backbone atoms
  BackboneList bb_list = initial_all_pos.ExtractBackbone(0, initial_all_pos.GetNumResidues());
  AllAtomPositions all_pos(bb_list);

  // reconstruct the sidechains 
  std::vector<Real> angles;
  for(uint i = 0; i < all_pos.GetNumResidues(); ++i){

    angles.assign(4, std::numeric_limits<Real>::quiet_NaN());

    for(uint j = 0; j < dihedral_angles[i].size(); ++j){
      angles[j] = dihedral_angles[i][j];
    }

    ConstructSidechainAtoms(all_pos, i,
                            angles[0],
                            angles[1],
                            angles[2],
                            angles[3]);
  }

  // compare
  
  // are all atoms set?
  for(uint i = 0; i < all_pos.GetNumResidues(); ++i){
    BOOST_CHECK(all_pos.IsAllSet(i));
  }

  // are they close?
  for(uint i = 0; i < all_pos.GetNumResidues(); ++i){
    uint start_idx = all_pos.GetFirstIndex(i);
    uint end_idx = all_pos.GetLastIndex(i);
    for(uint j = start_idx; j <= end_idx; ++j){
      Real d = geom::Distance(initial_all_pos.GetPos(j), all_pos.GetPos(j));
      BOOST_CHECK(d < 0.2);
    }
  } 
}

BOOST_AUTO_TEST_CASE(test_proline_fallback) {
  // get data
  ost::mol::EntityHandle test_ent = LoadTestStructure("data/all_aa.pdb");
  ost::mol::ResidueHandleList res_list = test_ent.GetResidueList();
  AllAtomPositions all_pos(res_list);

  // search for the proline
  int proline_idx = -1;
  for(uint i = 0; i < all_pos.GetNumResidues(); ++i){
    if(all_pos.GetAA(i) == ost::conop::PRO){
      proline_idx = i;
      break;
    }
  }

  if(proline_idx == -1){
    throw promod3::Error("There must be a proline in the test structure!");
  }

  // extract the dihedral angles
  Real chi1 = GetChi1Angle(all_pos, proline_idx);
  Real chi2 = GetChi2Angle(all_pos, proline_idx);
  geom::Vec3 cd_pos = all_pos.GetPos(proline_idx, PRO_CD_INDEX);

  // reconstruct proline by only providing chi1
  ConstructSidechainAtoms(all_pos, proline_idx, chi1);
  geom::Vec3 new_cd_pos = all_pos.GetPos(proline_idx, PRO_CD_INDEX);
  Real new_chi1 = GetChi1Angle(all_pos, proline_idx);
  Real new_chi2 = GetChi2Angle(all_pos, proline_idx);

  // the chi1 angle should be superclose (read same)
  BOOST_CHECK(std::abs(chi1 - new_chi1) < 0.00001);

  // the chi2 angle and the cd position have changed
  BOOST_CHECK(std::abs(chi2 - new_chi2) > 0.00001);
  BOOST_CHECK(geom::Distance(cd_pos, new_cd_pos) > 0.1);

  // but not too much
  BOOST_CHECK(std::abs(chi2 - new_chi2) < 0.4); // roughly 23 degrees
  BOOST_CHECK(geom::Distance(cd_pos, new_cd_pos) < 0.6);
}

BOOST_AUTO_TEST_SUITE_END();