Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

EMthreading parts #5

Open
wants to merge 5 commits into
base: main
Choose a base branch
from
Open

EMthreading parts #5

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
5 commits
Select commit Hold shift + click to select a range
d5e7d08
added extra tests to check connection restraints
dibyendu92 Aug 26, 2021
b786b68
identation problem were checked and more comments to explain the tabl…
dibyendu92 Aug 26, 2021
d10d062
two new tests were added to test chain id mismatch and high mismatch …
dibyendu92 Aug 26, 2021
c09d178
log 0 issue was taken care of by adding -0.1 to the score
dibyendu92 Aug 26, 2021
f033c29
we are not testing this cases for now, have to throughly check the sc…
dibyendu92 Aug 26, 2021
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
13 changes: 9 additions & 4 deletions src/SecondaryStructureParsimonyRestraint.cpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -48,10 +48,12 @@ double SecondaryStructureParsimonyRestraint::unprotected_evaluate(DerivativeAccu

// Get chain
std::string chain = threading::StructureElement(get_model(),ses_[nse]).get_chain();

if (boost::algorithm::contains(chain, prot_chain)) {

// First, get the residue indexes for this SE.
Ints resinds = threading::StructureElement(get_model(), ses_[nse]).get_resindex_list();

// Second, get the probabilities for this SE
Floats ses_ss_prob = atom::SecondaryStructureResidue(get_model(), ses_[nse]).get_all_probabilities();

Expand All @@ -70,15 +72,18 @@ double SecondaryStructureParsimonyRestraint::unprotected_evaluate(DerivativeAccu
// Loop over all the residues
// Compare to sequence ss_probs
for (int j = 0; j < nres; j++){

Floats seq_prob = atom::SecondaryStructureResidue(get_model(), resis[j].get_particle_index()).get_all_probabilities();
std::vector<float> se_prob = ses_ss_probs[j];

float dot = 0;
for (int k = 0; k < 3; k++){
//std::cout << "seq_prob " << seq_prob[k] << " se_prob " << se_prob[k] << std::endl;

dot += seq_prob[k] * se_prob[k];
}
score += -1 * log(dot+0.0001); // to avoid log(0) problem
if (dot < 0.001){
score += -0.1;
} else {score += -1.0 * log(dot);} // to avoid log(0) problem
//score += -1 * log(dot+0.0001); // to avoid log(0) problem
//std::cout << "RES# " << j << " (" << se_prob[0] << se_prob[1] << se_prob[2] << "), (" << seq_prob[0] << seq_prob[1] << seq_prob[2] << ") " << dot << std::endl;

}
Expand Down
48 changes: 15 additions & 33 deletions src/StructureElementConnectivityRestraint.cpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -20,8 +20,15 @@ StructureElementConnectivityRestraint::StructureElementConnectivityRestraint(Mod
std::string name)
: Restraint(m, "SEConnectivityRestraint%1%"), score_func_(score_func), a_(a), b_(b), n_sds_(n_sds) {}

/* Apply the pair score to each particle pair listed in the container.
/*
* Apply the pair score to each particle pair listed in the container.
* The mean and standard deviation of the distance-per-residue for loops based on number of residues and bounding secondary structure element was determined from the CATH S100 database as described
* in Section 2.4.1 in https://doi.org/10.1016/j.pbiomolbio.2019.09.003.
* Helix-loop-helix and sheet-loop-sheet values were computed using 1,668,774 and 614,152 Smotif elements, respectively.
* Sheet-loop-helix and helix-loop-sheet observations were combined, with a total of 1,416,244 Smotif elements (702,557 sheet-loop-helix and 713,687 helix-loop-sheet) used to compute the values.
* The values are taken from Table 1 in https://doi.org/10.1016/j.pbiomolbio.2019.09.003
*/

static float ll_means_[90] = {3.81, 3.036, 2.836, 2.511, 2.275, 2.178, 2.026, 1.876, 1.835, 1.669, 1.658, 1.666, 1.625, 1.53,
1.445, 1.374, 1.292, 1.212, 1.164, 1.133, 1.049, 1.043, 1.074, 0.977, 0.965, 0.938, 0.868, 0.824, 0.805,
0.788, 3.81, 3.19, 1.846, 1.607, 1.274, 1.14, 1.139, 1.198, 1.177, 1.115, 1.029, 1.048, 0.935, 0.91, 0.908,
Expand All @@ -37,30 +44,6 @@ static float ll_sds_[90] = {0.027, 0.284, 0.397, 0.441, 0.483, 0.499, 0.504, 0.5
0.373, 0.395, 0.47, 0.418, 0.36, 0.349, 0.359, 0.312, 0.302, 0.281, 0.279, 0.264, 0.259, 0.346, 0.257,
0.067, 0.278, 0.361, 0.418, 0.45, 0.448, 0.455, 0.436, 0.452, 0.438, 0.416, 0.407, 0.402, 0.411, 0.405,
0.381, 0.378, 0.373, 0.36, 0.372, 0.338, 0.322, 0.308, 0.285, 0.289, 0.296, 0.298, 0.294, 0.286, 0.208};
/*
static float nhh_means_[30] = {3.81, 3.036, 2.836, 2.511, 2.275, 2.178, 2.026, 1.876, 1.835, 1.669, 1.658, 1.666, 1.625, 1.53,
1.445, 1.374, 1.292, 1.212, 1.164, 1.133, 1.049, 1.043, 1.074, 0.977, 0.965, 0.938, 0.868, 0.824, 0.805,
0.788};

static float nss_means_[30] = {3.81, 3.19, 1.846, 1.607, 1.274, 1.14, 1.139, 1.198, 1.177, 1.115, 1.029, 1.048, 0.935, 0.91, 0.908,
0.85, 0.83, 0.852, 0.849, 0.761, 0.722, 0.742, 0.684, 0.677, 0.611, 0.587, 0.596, 0.565, 0.576, 0.532};

static float nsh_means_[30] = {3.809, 3.137, 2.818, 2.482, 2.154, 1.928, 1.749, 1.67, 1.531, 1.428, 1.377, 1.282, 1.261, 1.203,
1.135, 1.045, 1.004, 1.02, 0.977, 0.928, 0.865, 0.834, 0.811, 0.756, 0.761, 0.749, 0.777, 0.74, 0.655,
0.648};


static float nhh_sds_[30] = {0.027, 0.284, 0.397, 0.441, 0.483, 0.499, 0.504, 0.537, 0.534, 0.538, 0.545, 0.507, 0.494, 0.468,
0.447, 0.428, 0.439, 0.415, 0.432, 0.392, 0.382, 0.38, 0.401, 0.381, 0.38, 0.317, 0.328, 0.304, 0.318,
0.273};

static float nss_sds_[30] = {0.027, 0.313, 0.293, 0.469, 0.419, 0.474, 0.49, 0.505, 0.447, 0.501, 0.475, 0.479, 0.417, 0.451, 0.416,
0.373, 0.395, 0.47, 0.418, 0.36, 0.349, 0.359, 0.312, 0.302, 0.281, 0.279, 0.264, 0.259, 0.346, 0.257};

static float nsh_sds_[30] = {0.067, 0.278, 0.361, 0.418, 0.45, 0.448, 0.455, 0.436, 0.452, 0.438, 0.416, 0.407, 0.402, 0.411, 0.405,
0.381, 0.378, 0.373, 0.36, 0.372, 0.338, 0.322, 0.308, 0.285, 0.289, 0.296, 0.298, 0.294, 0.286, 0.208};

*/

float StructureElementConnectivityRestraint::get_sd_distance_per_residue() const {
int nres = get_number_of_residues();
Expand Down Expand Up @@ -112,8 +95,7 @@ float StructureElementConnectivityRestraint::get_model_distance() const {
//std::cout << "chain 1: " << chain_1 << " & chain 2: " << chain_2 << std::endl;
// Compute and return the distance
if (chain_1 == chain_2){
//int ra_i = threading::StructureElement(get_model(),a_).get_start_res();
//int rb_i = threading::StructureElement(get_model(),b_).get_start_res();

algebra::Vector3D a_coords = threading::StructureElement(get_model(), a_).get_coordinates().back();
//std::cout<<" a_coords "<<a_coords<<std::endl;
algebra::Vector3D b_coords = threading::StructureElement(get_model(), b_).get_coordinates().back();
Expand All @@ -140,16 +122,16 @@ float StructureElementConnectivityRestraint::get_mean_distance_per_residue() con
Floats ses_ss_prob_a = atom::SecondaryStructureResidue(get_model(), a_).get_all_probabilities();
Floats ses_ss_prob_b = atom::SecondaryStructureResidue(get_model(), b_).get_all_probabilities();
if (ses_ss_prob_a[0]==1 && ses_ss_prob_b[0]==1) {
sset = 0;
} else if (ses_ss_prob_a[1]==1 && ses_ss_prob_b[1]==1){
sset = 1;
sset = 0;
} else if (ses_ss_prob_a[1]==1 && ses_ss_prob_b[1]==1){
sset = 1;
} else if ((ses_ss_prob_a[1]==0 && ses_ss_prob_b[1]==1) || (ses_ss_prob_a[1]==1 && ses_ss_prob_b[1]==0)){
sset = 2;
sset = 2;
} else {sset = 0;}

if (nres > 29) {
if (nres > 29) {
mean = ll_means_[29+(sset*30)];
}else { mean = ll_means_[nres+(sset*30)];};
}else {mean = ll_means_[nres+(sset*30)];};

return mean;
};
Expand Down
Loading