You are here

scoring with omega energy term

3 posts / 0 new
Last post
scoring with omega energy term
#1

Hello,

During my work with score function there was a situation when only one omega backbone angle has a change and I noticed some feature. If you set only one omega angle without changing other backbone torsions and then scoring there will be no change in total energy. Moreover, when I drowned myself a bit deeper I have found out that changing two omegas of different residues do the change in energy but still wrongly. Here complete source and output (first for omega term, then for talaris2014) showing the problem, plot also shows output. This feature appears both when you use set_omega function or chaging via DOF_ID. I'm using 2017.08 build.

#include <devel/init.hh>
#include <core/pose/annotated_sequence.hh>
#include <core/pose/util.hh>
#include <core/pose/Pose.hh>
#include <core/scoring/Energies.hh>
#include <core/scoring/ScoreFunction.hh>
#include <core/scoring/ScoreFunctionFactory.hh>
#include <iostream>

int main(int argc, char * argv[])
{
    devel::init(argc, argv);

    //core::scoring::ScoreFunctionOP scorefx = core::scoring::ScoreFunctionFactory::create_score_function("talaris2014");

    core::scoring::ScoreFunctionOP scorefx(new core::scoring::ScoreFunction);
    scorefx->set_weight(core::scoring::omega, 1.0);

    core::pose::Pose pose1, pose2, pose3, pose4;
    std::string peptide_seq = "AAAAA";
    core::pose::make_pose_from_sequence(pose1, peptide_seq, *core::chemical::ChemicalManager::get_instance()->residue_type_set(core::chemical::FA_STANDARD));

    for (size_t i = 1; i != pose1.total_residue(); i++)
    {
        pose1.set_phi(i, -135.0);
        pose1.set_psi(i, 135.0);
        pose1.set_omega(i, 179.8);
    }

    pose4 = pose3 = pose2 = pose1;

    size_t steps = 6;
    double range = 180;
    double start = 90;

    for (size_t i = 0; i != steps + 1; i++)
    {
        pose1.set_omega(3, start + i * range / steps);
        pose1.set_phi(3, -135);

        pose2.set_omega(3, start + i * range / steps);
        pose2.set_phi(2, -135);

        pose3.set_omega(3, start + i * range / steps);
        pose3.set_omega(4, start + i * range / steps);

        pose4.set_omega(3, start + i * range / steps);

        (*scorefx)(pose1);
        (*scorefx)(pose2);
        (*scorefx)(pose3);
        (*scorefx)(pose4);
        std::cout << std::scientific << start + i * range / steps << '\t'
        << pose1.energies().total_energy() << '\t'
        << pose2.energies().total_energy() << '\t'
        << pose3.energies().total_energy() << '\t'
        << pose4.energies().total_energy() << std::endl;
    }
}

// core::scoring::omega only

9.000000e+01    8.100120e+01    8.100120e+01    1.620008e+02    8.100120e+01
1.200000e+02    3.600120e+01    8.100120e+01    1.170008e+02    8.100120e+01
1.500000e+02    9.001200e+00    8.100120e+01    9.000080e+01    8.100120e+01
1.800000e+02    1.200000e-03    8.100120e+01    8.100080e+01    8.100120e+01
2.100000e+02    9.001200e+00    8.100120e+01    9.000080e+01    8.100120e+01
2.400000e+02    3.600120e+01    8.100120e+01    1.170008e+02    8.100120e+01
2.700000e+02    8.100120e+01    8.100120e+01    1.620008e+02    8.100120e+01

// talaris2014

9.000000e+01    5.534642e+01    5.534642e+01    1.056326e+02    5.534642e+01
1.200000e+02    2.728259e+01    5.540759e+01    7.771807e+01    5.540759e+01
1.500000e+02    1.054679e+01    5.554679e+01    6.108479e+01    5.554679e+01
1.800000e+02    5.057364e+00    5.568236e+01    5.568270e+01    5.568236e+01
2.100000e+02    1.073223e+01    5.573223e+01    6.137627e+01    5.573223e+01
2.400000e+02    2.754141e+01    5.566641e+01    7.804450e+01    5.566641e+01
2.700000e+02    5.557881e+01    5.557881e+01    1.060839e+02    5.557881e+01

So is that mean that I always need to catch such situations and reset other torsions when dealing with omegas? If I understand scoring and omega term correctly pose 1, 2, and 4 energies must have identical values, or I am not right?

Thanks in advance

AttachmentSize
score.png34.91 KB
Category: 
Post Situation: 
Thu, 2017-03-09 12:03
SergeyP

I can confirm that I see the same behavior but it's not immediately obvious to me that this is correct or where the bug is if it is incorrect.  I still haven't decided if there is a pose number range issue (they are numbered from 1, not 0, so total_residue represents a number in the pose, and zero does not) or a degrees versus radians issue.

Thu, 2017-03-09 19:11
smlewis

Thanks for the answer but I don't understand you. What do you mean under a pose number range issue? Is residue numbering influence the energy scoring? It isn't look like a degrees versus radians issue since, as I mentioned, there are same values when working with dof.
 

	size_t steps = 6;
	double range = numeric::NumericTraits<core::Real>::pi();
	double start = numeric::NumericTraits<core::Real>::pi() / 2;

	core::id::TorsionID phi_id2(2, core::id::BB, core::id::phi_torsion);
	core::id::TorsionID phi_id3(3, core::id::BB, core::id::phi_torsion);

	core::id::TorsionID omega_id3(3, core::id::BB, core::id::omega_torsion);
	core::id::TorsionID omega_id4(4, core::id::BB, core::id::omega_torsion);

	for (size_t i = 0; i != steps + 1; i++)
	{
		double angle = start + i * range / steps;

		pose1.set_dof(pose1.conformation().dof_id_from_torsion_id(omega_id3), angle);
		pose1.set_dof(pose1.conformation().dof_id_from_torsion_id(phi_id3),  pose1.dof(pose1.conformation().dof_id_from_torsion_id(phi_id3)));

		pose2.set_dof(pose2.conformation().dof_id_from_torsion_id(omega_id3), angle);
		pose2.set_dof(pose2.conformation().dof_id_from_torsion_id(phi_id2), pose2.dof(pose2.conformation().dof_id_from_torsion_id(phi_id2)));

		pose3.set_dof(pose3.conformation().dof_id_from_torsion_id(omega_id3), angle);
		pose3.set_dof(pose3.conformation().dof_id_from_torsion_id(omega_id4), angle);

		pose4.set_dof(pose4.conformation().dof_id_from_torsion_id(omega_id3), angle);

		(*scorefx)(pose1);
		(*scorefx)(pose2);
		(*scorefx)(pose3);
		(*scorefx)(pose4);
		std::cout << std::scientific << angle << 
		'\t' << pose1.energies().total_energy() << 
		'\t' << pose2.energies().total_energy() << 
		'\t' << pose3.energies().total_energy() << 
		'\t' << pose4.energies().total_energy() << std::endl;

	}

//

1.570796e+00    8.100120e+01    8.100120e+01    1.620008e+02    8.100120e+01
2.094395e+00    3.600120e+01    8.100120e+01    1.170008e+02    8.100120e+01
2.617994e+00    9.001200e+00    8.100120e+01    9.000080e+01    8.100120e+01
3.141593e+00    1.200000e-03    8.100120e+01    8.100080e+01    8.100120e+01
3.665191e+00    9.001200e+00    8.100120e+01    9.000080e+01    8.100120e+01
4.188790e+00    3.600120e+01    8.100120e+01    1.170008e+02    8.100120e+01
4.712389e+00    8.100120e+01    8.100120e+01    1.620008e+02    8.100120e+01

 

Fri, 2017-03-10 11:16
SergeyP