Skip to content
Snippets Groups Projects
Commit 09267d58 authored by Gerardo Tauriello's avatar Gerardo Tauriello
Browse files

Bugfixed BFT testset gen. and using new copy constructor for rotamers.

parent f6ad5cb8
No related branches found
No related tags found
No related merge requests found
"""
GOAL: extract loop examples from structure DB
OUT: fragments.json: list of fragment infos ([chain_idx, offset, length])
-> fragments are linked to PM3's default structure DB
-> each fragment has at least 50% coils and doesn't touch termini
"""
import random, json
......@@ -12,7 +14,7 @@ from promod3 import loop
###############################################################################
min_coord_size = 50 # min. size for StructureDB entry to be considered
terminal_res_to_skip = 2 # first and last x residues are never chosen for loop
num_chains = 3973 # get x random chains from structure DB
num_chains = 4000 # get x random chains from structure DB
num_loops_per_len = 5000 # number of loops to pick for each loop length
loop_range = range(3,13) # range of loop lengths (uniformly distributed)
random.seed(42) # fixed seed for reproducibilty
......@@ -88,7 +90,7 @@ for loop_length in loop_range:
for start_idx in range(terminal_res_to_skip, range_end):
num_coil = sum(is_coil[start_idx:start_idx+loop_length])
to_skip = any(res_to_skip[start_idx:start_idx+loop_length])
if num_coil >= num_coil_thresh:
if num_coil >= num_coil_thresh and not to_skip:
possible_start_indices.append(start_idx)
# show it
print "LOOP LEN %d: %d choices" % (loop_length, len(possible_start_indices))
......@@ -96,9 +98,23 @@ for loop_length in loop_range:
start_indices = random.sample(possible_start_indices, num_loops_per_len)
for start_idx in start_indices:
coord_idx = back_mapping[start_idx]
fragment = loop.FragmentInfo(coord_idx, start_idx - first_idx[coord_idx],
loop_length)
fragments.append(fragment)
offset = start_idx - first_idx[coord_idx]
fragments.append(loop.FragmentInfo(coord_idx, offset, loop_length))
# check all fragments
for fragment in fragments:
# check termini
coord_info = structure_db.GetCoordInfo(fragment.chain_index)
assert(fragment.offset >= terminal_res_to_skip)
assert(fragment.offset + fragment.length + terminal_res_to_skip\
<= coord_info.size)
# check coils
frag_info = loop.FragmentInfo(fragment.chain_index, 0, coord_info.size)
dssp_states = structure_db.GetDSSPStates(frag_info)
coilness = [1 if el in ['T','C','S'] else 0 for el in dssp_states]
assert(len(coilness) == coord_info.size)
num_coil = sum(coilness[fragment.offset:fragment.offset+fragment.length])
assert(num_coil >= round(fragment.length*0.5))
# store it
with open(out_fragments, "w") as json_file:
......
......@@ -261,9 +261,7 @@ void SidechainReconstructor::CollectRotamerGroups_(
env_->GetRotamerGroup(rot_group, res_idx);
if (rot_group) {
// add a copy (we modify rotamer in ApplySelfEnergyThresh!)
RotamerGroupPtr new_group(
new RotamerGroup(rot_group->GetRotamers(),
rot_group->GetResidueIndex()));
RotamerGroupPtr new_group(new RotamerGroup(*rot_group));
rotamer_groups.push_back(new_group);
res->rotamer_res_indices.push_back(i);
}
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment