SSAGES::AntiBetaRMSDCV Class Reference

Collective variable to measure antiparallel beta secondary structure. More...

#include <AntiBetaRMSDCV.h>

Inheritance diagram for SSAGES::AntiBetaRMSDCV:

Public Member Functions
	AntiBetaRMSDCV (std::vector< int > resids, std::string refpdb, double unitconv, int mode)
	Constructor. More...
void	Initialize (const Snapshot &snapshot) override
	Initialize necessary variables. More...
void	Evaluate (const Snapshot &snapshot) override
	Evaluate the CV. More...
Public Member Functions inherited from SSAGES::CollectiveVariable
	CollectiveVariable ()
	Constructor.
virtual	~CollectiveVariable ()
	Destructor.
virtual void	Initialize (const class Snapshot &)
	Initialize CV. More...
virtual void	Evaluate (const class Snapshot &)=0
	Evaluate CV. More...
double	GetValue () const
	Get current value of the CV. More...
virtual double	GetMinimumImage (double) const
	Returns the minimum image of a CV based on the input location. More...
virtual double	GetPeriodicValue (double location) const
	Apply periodic boundaries to a given value. More...
const std::vector< Vector3 > &	GetGradient () const
	Get current gradient of the CV. More...
const Matrix3 &	GetBoxGradient () const
	Get gradient contribution to box. More...
const std::array< double, 2 > &	GetBoundaries ()
	Get CV boundaries. More...
virtual double	GetDifference (double location) const

Static Public Member Functions
static AntiBetaRMSDCV *	Build (const Json::Value &json, const std::string &path)
	Set up collective variable. More...
Static Public Member Functions inherited from SSAGES::CollectiveVariable
static CollectiveVariable *	BuildCV (const Json::Value &json, const std::string &path)
	Set up collective variable. More...

Private Attributes
std::vector< int >	resids_
	Residue IDs for secondary structure calculation.
std::vector< std::vector< std::string > >	atomids_
	Atom IDs and chain IDs for secondary structure calculation: backbone of resids_.
std::vector< int >	ids_only_
	Atom IDs only, for secondary structure calculation: backbone atoms of resids_.
std::string	refpdb_
	Name of pdb reference for system.
double	unitconv_
	Length unit conversion: convert 1 nm to your internal MD units (ex. if using angstroms use 10)
int	mode_
	Specify whether to calculate beta sheet character in intra or inter mode: 0 for either, 1 for inter, 2 for intra.
std::vector< double >	refdists_
	Pairwise distance between backbone atoms for reference structure.

Additional Inherited Members
Protected Attributes inherited from SSAGES::CollectiveVariable
std::vector< Vector3 >	grad_
	Gradient vector dCv/dxi.
Matrix3	boxgrad_
	Gradient w.r.t box vectors dCv/dHij.
double	val_
	Current value of CV.
std::array< double, 2 >	bounds_
	Bounds on CV.

Detailed Description

Collective variable to measure antiparallel beta secondary structure.

Following treatment in Pietrucci and Laio, "A Collective Variable for the Efficient Exploration of Protein Beta-Sheet Structures: Application to SH3 and GB1", JCTC, 2009, 5(9): 2197-2201.

Check 2 blocks of 3 consecutive protein residues for RMSD from reference "ideal" antiparallel beta sheet structure.

Definition at line 43 of file AntiBetaRMSDCV.h.

Constructor & Destructor Documentation

AntiBetaRMSDCV()

SSAGES::AntiBetaRMSDCV::AntiBetaRMSDCV ( std::vector< int >  resids, std::string  refpdb, double  unitconv, int  mode  )

inline

Constructor.

Parameters

resids	IDs of residues for calculating secondary structure
refpdb	String of pdb filename with atom and residue indices.
unitconv	Conversion for internal MD length unit: 1 nm is equal to unitconv internal units
mode	Specification for inter/intra mode for beta sheet character

Construct an AntibetaRMSD CV – calculates anti-beta sheet character by summing pairwise RMSD to an ideal anti-beta sheet structure for all possible 6 residue segments.

Definition at line 80 of file AntiBetaRMSDCV.h.

                                                                                           :
        resids_(resids), refpdb_(refpdb), unitconv_(unitconv), mode_(mode)
        {
            if(resids_.size() != 2 ){
                throw std::invalid_argument("AntiBetaRMSDCV: Input must designate range of residues with 2 residue numbers.");
            }
 
            resids_.clear();
 
            if(resids[0] >= resids[1]){
                throw std::invalid_argument("AntiBetaRMSDCV: Input must list lower residue index first: please reverse residue range.");
            } else if(resids[1] - resids[0] < 5){
                throw std::invalid_argument("AntiBetaRMSDCV: Residue range must span at least 6 residues for secondary structure calculation.");
            }
 
            std::cout << "AntiBetaRMSDCV: Calculating antiparallel beta sheet character from residue " << resids[0] << " to " << resids[1]  << "." << std::endl;
 
            for(int i = resids[0]; i <= resids[1]; i++){
                resids_.push_back(i);
            }
        }

References resids_.

Referenced by Build().

Here is the caller graph for this function:

Member Function Documentation

Build()

static AntiBetaRMSDCV* SSAGES::AntiBetaRMSDCV::Build ( const Json::Value &  json, const std::string &  path  )

inlinestatic

Set up collective variable.

Parameters

json	JSON input value.
path	Path for JSON path specification.

Returns

Pointer to the CV built by this function. nullptr in case of unknown error.

Builds a CV from a JSON node. Returns a pointer to the built cv. If an unknown error is encountered, this function will return a nullptr, but generally it will throw a BuildException on failure.

Warning

Object lifetime is the caller's responsibility.

Definition at line 256 of file AntiBetaRMSDCV.h.

        {
            Json::ObjectRequirement validator;
            Json::Value schema;
            Json::CharReaderBuilder rbuilder;
            Json::CharReader* reader = rbuilder.newCharReader();
 
            reader->parse(JsonSchema::AntiBetaRMSDCV.c_str(),
                          JsonSchema::AntiBetaRMSDCV.c_str() + JsonSchema::AntiBetaRMSDCV.size(),
                          &schema, nullptr);
            validator.Parse(schema, path);
 
            //Validate inputs
            validator.Validate(json, path);
            if(validator.HasErrors())
                    throw BuildException(validator.GetErrors());
 
            std::vector<int> resids;
            for(auto& s : json["residue_ids"])
                resids.push_back(s.asInt());
            auto reference = json["reference"].asString();
 
            double unitconv = json.get("length_unit", 1).asDouble();
 
            int mode = json.get("mode", 0).asInt();
 
            return new AntiBetaRMSDCV(resids, reference, unitconv, mode);
        }

References AntiBetaRMSDCV(), Json::Requirement::GetErrors(), Json::Requirement::HasErrors(), Json::ObjectRequirement::Parse(), and Json::ObjectRequirement::Validate().

Referenced by SSAGES::CollectiveVariable::BuildCV().

Here is the call graph for this function:

Here is the caller graph for this function:

Evaluate()

void SSAGES::AntiBetaRMSDCV::Evaluate ( const Snapshot &  snapshot )

inlineoverride

Evaluate the CV.

Parameters

snapshot

Current simulation snapshot.

Definition at line 174 of file AntiBetaRMSDCV.h.

        {
            // need atom positions for all atoms
            const auto& pos = snapshot.GetPositions();
            auto& comm = snapshot.GetCommunicator();
 
            unsigned int resgroups = resids_.size() - 2;
            double rmsd, dist_norm, dxgrouprmsd;
            int localidx;
            Vector3 dist_xyz, dist_ref;
            std::vector<Vector3> refxyz;
            std::vector< std::vector< Vector3 > > deriv(30, std::vector<Vector3>(30, Vector3{0,0,0}));
            std::fill(grad_.begin(), grad_.end(), Vector3{0,0,0});
            grad_.resize(snapshot.GetNumAtoms(), Vector3{0,0,0});
            val_ = 0.0;
 
            for(size_t i = 0; i < resgroups - 3; i++){
                for(size_t j = i + 3; j < resgroups; j++){
                    // mode: 0 for all, 1 for inter, 2 for intra
                    if((mode_ == 0) || (mode_ == 1 && atomids_[1][5 * j] != atomids_[1][5 * i]) || (mode_ == 2 && atomids_[1][5 * j] == atomids_[1][5 * i])){
                        rmsd = 0.0;
                        std::fill(refxyz.begin(), refxyz.end(), Vector3{0,0,0});
                        refxyz.resize(30, Vector3{0,0,0});
 
                        for(size_t k = 0; k < 15; k++){
                            localidx = snapshot.GetLocalIndex(ids_only_[5 * i + k]);
                            if(localidx != -1) refxyz[k] = pos[localidx];
                        }
 
                        for(size_t k = 0; k < 15; k++){
                            localidx = snapshot.GetLocalIndex(ids_only_[5 * j + k]);
                            if(localidx != -1) refxyz[k + 15] = pos[localidx];
                        }
 
                        MPI_Allreduce(MPI_IN_PLACE, refxyz.data(), 90, MPI_DOUBLE, MPI_SUM, comm);
 
                        // sum over all pairs to calculate CV
                        for(size_t k = 0; k < 29; k++){
                            for(size_t h = k + 1; h < 30; h++){
                                dist_xyz = refxyz[k] - refxyz[h];
                                dist_norm = dist_xyz.norm() - refdists_[29 * k + h];
                                rmsd += dist_norm * dist_norm;
                                deriv[k][h] = dist_xyz * dist_norm / dist_xyz.norm();
                            }
                        }
 
                        rmsd = pow(rmsd / 435, 0.5) / 0.1;
                        val_ += (1 - pow(rmsd, 8.0)) / (1 - pow(rmsd, 12.0));
 
                        dxgrouprmsd = pow(rmsd, 11.0) + pow(rmsd, 7.0);
                        dxgrouprmsd /= pow(rmsd, 8.0) + pow(rmsd, 4.0) + 1;
                        dxgrouprmsd /= pow(rmsd, 8.0) + pow(rmsd, 4.0) + 1;
                        dxgrouprmsd *= -40. / 435;
 
                        for(size_t k = 0; k < 15; k++){
                            for(size_t h = k + 1; h < 15; h++){
                                localidx = snapshot.GetLocalIndex(ids_only_[5 * i + k]);
                                if(localidx != -1) grad_[localidx] += dxgrouprmsd * deriv[k][h];
                                localidx = snapshot.GetLocalIndex(ids_only_[5 * i + h]);
                                if (localidx != -1) grad_[localidx] += dxgrouprmsd * deriv[k][h];
                            }
                            for(size_t h = 0; h < 15; h++){
                                localidx = snapshot.GetLocalIndex(ids_only_[5 * i + k]);
                                if(localidx != -1) grad_[localidx] += dxgrouprmsd * deriv[k][h + 15];
                                localidx = snapshot.GetLocalIndex(ids_only_[5 * j + h]);
                                if(localidx != -1) grad_[localidx] += dxgrouprmsd * deriv[k][h + 15];
                            }
                        }
                        for(size_t k = 0; k < 14; k++){
                            for(size_t h = k + 1; h < 15; h++){
                                localidx = snapshot.GetLocalIndex(ids_only_[5 * j + k]);
                                if(localidx != -1) grad_[localidx] += dxgrouprmsd * deriv[k + 15][h + 15];
                                localidx = snapshot.GetLocalIndex(ids_only_[5 * j + h]);
                                if(localidx != -1) grad_[localidx] += dxgrouprmsd * deriv[k + 15][h + 15];
                            }
                        }
                    }
                }
            }
        }

References atomids_, SSAGES::Snapshot::GetCommunicator(), SSAGES::Snapshot::GetLocalIndex(), SSAGES::Snapshot::GetNumAtoms(), SSAGES::Snapshot::GetPositions(), SSAGES::CollectiveVariable::grad_, ids_only_, mode_, refdists_, resids_, and SSAGES::CollectiveVariable::val_.

Here is the call graph for this function:

Initialize()

void SSAGES::AntiBetaRMSDCV::Initialize ( const Snapshot &  snapshot )

inlineoverride

Initialize necessary variables.

Parameters

snapshot

Current simulation snapshot.

Definition at line 106 of file AntiBetaRMSDCV.h.

        {
            atomids_ = ReadBackbone::GetPdbBackbone(refpdb_, resids_);
            ids_only_.resize(atomids_[0].size());
 
            // check to find all necessary atoms
            std::transform(atomids_[0].begin(), atomids_[0].end(), ids_only_.begin(), [](std::string val) {return std::stoi(val);});
            using std::to_string;
            std::vector<int> found;
            snapshot.GetLocalIndices(ids_only_, &found);
            size_t nfound = found.size();
            MPI_Allreduce(MPI_IN_PLACE, &nfound, 1, MPI_INT, MPI_SUM, snapshot.GetCommunicator());
            if(nfound != resids_.size() * 5)
                throw BuildException({
                    "AntiBetaRMSDCV: Expected to find " +
                    to_string(resids_.size() * 5) +
                    " atoms, but only found " +
                    to_string(nfound) + "."
                });
 
            // reference 'ideal' antiparallel beta sheet structure in nanometers
            // Reference structure taken from values used in Plumed 2.2.0
            std::vector<Vector3> refalpha;
            refalpha.push_back(unitconv_ * Vector3{ .2263, -.3795,  .1722}); // N    i
            refalpha.push_back(unitconv_ * Vector3{ .2493, -.2426,  .2263}); // CA
            refalpha.push_back(unitconv_ * Vector3{ .3847, -.1838,  .1761}); // CB
            refalpha.push_back(unitconv_ * Vector3{ .1301, -.1517,  .1921}); // C
            refalpha.push_back(unitconv_ * Vector3{ .0852, -.1504,  .0739}); // O
            refalpha.push_back(unitconv_ * Vector3{ .0818, -.0738,  .2917}); // N    i+1
            refalpha.push_back(unitconv_ * Vector3{-.0299,  .0243,  .2748}); // CA
            refalpha.push_back(unitconv_ * Vector3{-.1421, -.0076,  .3757}); // CB
            refalpha.push_back(unitconv_ * Vector3{ .0273,  .1680,  .2854}); // C
            refalpha.push_back(unitconv_ * Vector3{ .0902,  .1993,  .3888}); // O
            refalpha.push_back(unitconv_ * Vector3{ .0119,  .2532,  .1813}); // N    i+2
            refalpha.push_back(unitconv_ * Vector3{ .0683,  .3916,  .1680}); // CA
            refalpha.push_back(unitconv_ * Vector3{ .1580,  .3940,  .0395}); // CB
            refalpha.push_back(unitconv_ * Vector3{-.0394,  .5011,  .1630}); // C
            refalpha.push_back(unitconv_ * Vector3{-.1459,  .4814,  .0982}); // O
            refalpha.push_back(unitconv_ * Vector3{-.2962,  .3559, -.1359}); // N    j - 2
            refalpha.push_back(unitconv_ * Vector3{-.2439,  .2526, -.2287}); // CA
            refalpha.push_back(unitconv_ * Vector3{-.1189,  .3006, -.3087}); // CB
            refalpha.push_back(unitconv_ * Vector3{-.2081,  .1231, -.1520}); // C
            refalpha.push_back(unitconv_ * Vector3{-.1524,  .1324, -.0409}); // O
            refalpha.push_back(unitconv_ * Vector3{-.2326,  .0037, -.2095}); // N    j - 1
            refalpha.push_back(unitconv_ * Vector3{-.1858, -.1269, -.1554}); // CA
            refalpha.push_back(unitconv_ * Vector3{-.3053, -.2199, -.1291}); // CB
            refalpha.push_back(unitconv_ * Vector3{-.0869, -.1949, -.2512}); // C
            refalpha.push_back(unitconv_ * Vector3{-.1255, -.2070, -.3710}); // O
            refalpha.push_back(unitconv_ * Vector3{ .0326, -.2363, -.2072}); // N    j
            refalpha.push_back(unitconv_ * Vector3{ .1405, -.2992, -.2872}); // CA
            refalpha.push_back(unitconv_ * Vector3{ .2699, -.2129, -.2917}); // CB
            refalpha.push_back(unitconv_ * Vector3{ .1745, -.4399, -.2330}); // C
            refalpha.push_back(unitconv_ * Vector3{ .1899, -.4545, -.1102}); // O
 
            // calculate all necessary pairwise distances for reference structure
            std::fill(refdists_.begin(), refdists_.end(), 0.0);
            refdists_.resize(900, 0.0);
            for(size_t k = 0; k < 29; k++){
                for(size_t h = k + 1; h < 30; h++){
                    refdists_[29 * k + h] = (refalpha[k] - refalpha[h]).norm();
                }
            }
        }

References atomids_, SSAGES::Snapshot::GetCommunicator(), SSAGES::Snapshot::GetLocalIndices(), SSAGES::ReadBackbone::GetPdbBackbone(), ids_only_, refdists_, refpdb_, resids_, and unitconv_.

Here is the call graph for this function:

---

The documentation for this class was generated from the following file:

• /home/michael/development/SSAGES/src/CVs/AntiBetaRMSDCV.h