Lanthipeptide and recon multistate design update

database/chemical/residue_type_sets/fa_standard/patches/Nterm2_Oxobutrylation.txt

@@ -70,6 +70,7 @@ SET_ICOOR  3H4 -60.0    70.0   1.089  C4  C3  C2
 #Unpatched H gets its geometry from LOWER, which we have removed
 
 SET_POLYMER_CONNECT LOWER NONE
+ADD_PROPERTY LOWER_TERMINUS ## implies terminus
 
 # If modified, no longer canonical.
 DELETE_PROPERTY CANONICAL_AA

source/src/protocols/cyclic_peptide/crosslinker/lanthionine_util.hh

@@ -25,10 +25,14 @@ namespace protocols {
 namespace cyclic_peptide {
 namespace crosslinker {
 
-// from M062x_def2TZVP simulation
+// from Methionine
 constexpr core::Real LANTHIONINE_UTIL_LANTHIONINE_BOND_LENGTH = 1.807;
 constexpr core::Real LANTHIONINE_UTIL_LANTHIONINE_BOND_C_ANGLE = 2.011; // 115.230 degrees
 constexpr core::Real LANTHIONINE_UTIL_LANTHIONINE_BOND_S_ANGLE = 1.801; // 103.176 degrees
+// from MP2/cc-pvtz geometry optimization
+constexpr core::Real LANTHIONINE_UTIL_METHYLLANTHIONINE_BOND_CG_ANGLE = 1.864; // 106.790 degrees SCC to CG
+constexpr core::Real LANTHIONINE_UTIL_METHYLLANTHIONINE_BOND_CA_ANGLE = 1.963; // 112.460 degrees SCC to CA
+constexpr core::Real LANTHIONINE_UTIL_METHYLLANTHIONINE_BOND_CB_ANGLE = 1.927; // 110.440 degrees CA-CB-CG
 
 /// @brief Given a pose and two residues, set up the lanthionine variant types.
 /// @details Sidechainres gets SIDECHAIN_CONJUGATION; dalares gets ACETYLATED_NTERMINUS_CONNECTION_VARIANT.

source/src/protocols/cyclic_peptide_predict/SimpleCycpepPredictApplication.cc

@@ -3459,7 +3459,7 @@ SimpleCycpepPredictApplication::add_closebond_logic_lanthipeptide(
 	debug_assert( cyclization_type() == SCPA_lanthipeptide );
 	core::Size ala_res = 0;
 	core::Size cys_res = 0;
-	if ( pose->residue_type(cyclization_point_start).base_name() == "CYS" ) {
+	if ( pose->residue_type(cyclization_point_start).base_name() == "CYS" ||  pose->residue_type(cyclization_point_start).base_name() == "DCYS" ) {
 		ala_res = cyclization_point_end;
 		cys_res = cyclization_point_start;
 	} else {

source/src/protocols/recon_design/FindConsensusSequence.cc

@@ -196,7 +196,7 @@ void FindConsensusSequence::apply( core::pose::Pose & pose ) {
 	parse_resfiles();
 
 	// Get the designable sequences of all poses
-	utility::vector1< std::string > all_pose_sequences;
+	utility::vector1< utility::vector1< std::string > > all_pose_sequences;
 	for ( core::Size ii = 1; ii <= designable_residues_.size(); ++ii ) {
 		all_pose_sequences.push_back(
 			get_designable_sequence ( *poses_[ii], designable_residues_[ii] )
@@ -230,24 +230,41 @@ void FindConsensusSequence::apply_mpi( core::pose::Pose & pose ) {
 	utility::vector1< core::Size > my_designable_residues = get_designable_residues( pose, this_nodes_resfile );
 
 	/// Make a string out of the AAs at my designable positions in the current state
-	std::string my_sequence = get_designable_sequence ( pose, my_designable_residues );
+	utility::vector1< std::string > my_sequence = get_designable_sequence ( pose, my_designable_residues );
 
-	utility::vector1<std::string> all_pose_sequences (n_procs_);
+    std::string pass_seq;
+    for( const std::string& resi_base_name: my_sequence) {
+        pass_seq += resi_base_name + " ";
+    }
+
+    utility::vector1< utility::vector1<std::string> > all_pose_sequences (n_procs_);
 	for ( core::Size ii = 1; ii <= n_procs_; ++ii ) {
 		if ( rank_ == ii ) {
 			for ( core::Size jj = 1; jj <= n_procs_; ++jj ) {
-				if ( rank_!=jj ) utility::send_string_to_node( jj-1, my_sequence ); // node ranks are 0-indexed
+				if ( rank_!=jj ) utility::send_string_to_node( jj-1, pass_seq ); // node ranks are 0-indexed
 				else all_pose_sequences[jj] = my_sequence;
 			}
 		} else {
-			all_pose_sequences[ii] = utility::receive_string_from_node( ii-1 ); // node ranks are 0-indexed
+			//all_pose_sequences[ii] = utility::receive_string_from_node( ii-1 ); // node ranks are 0-indexed
+            std::string passed_seq;
+            passed_seq = utility::receive_string_from_node( ii-1 ); // node ranks are 0-indexed
+            //Need to split passed_seq by spaces
+            std::istringstream iss(passed_seq);
+            std::string resi_name;
+            while (iss >> resi_name){
+                all_pose_sequences[ii].push_back(resi_name);
+            }
 		}
 	}
 
 	if ( debug_ ) {
 		TR << "printing all AAs at designable positions to check indices: "<< std::endl;
 		for ( core::Size i = 1; i <= all_pose_sequences.size(); ++i ) {
-			TR << i << ":" << all_pose_sequences[ i ] << ", ";
+			TR << i << ":";
+            for (const std::string& resi : all_pose_sequences[i]) {
+                TR << resi << "-";
+            }
+            TR << ", ";
 		}
 		TR << std::endl;
 	}

source/src/protocols/recon_design/MSDMover.cc

@@ -42,6 +42,7 @@
 
 #include <protocols/minimization_packing/PackRotamersMover.hh>
 
+#include <sstream>
 #include <utility/string_util.hh>
 #include <utility/mpi_util.hh>
 
@@ -168,7 +169,7 @@ MSDMover::apply( core::pose::Pose & pose ) {
 	setup_mover( pose );
 
 	// Find the sequence of the other poses
-	utility::vector1< std::string > other_pose_sequences;
+    utility::vector1< utility::vector1< std::string > > other_pose_sequences;
 	for  ( core::Size ii = 1; ii <= poses_.size(); ++ii ) {
 		if ( ii != current_pose_ ) {
 			other_pose_sequences.push_back(
@@ -217,7 +218,7 @@ MSDMover::apply( core::pose::Pose & pose ) {
 /// applied so they can be removed later
 utility::vector1< core::scoring::constraints::ConstraintCOP >
 MSDMover::apply_linked_constraints( core::pose::Pose & pose,
-	utility::vector1< std::string > other_pose_sequences,
+	utility::vector1< utility::vector1< std::string > > other_pose_sequences,
 	utility::vector1< core::Size > my_designable_residues ) {
 
 	using namespace core::scoring::constraints;
@@ -280,34 +281,50 @@ MSDMover::apply_mpi( core::pose::Pose & pose ) {
 	/// Find out my designable residues from my pose and my resfile
 	utility::vector1< core::Size > my_designable_residues = get_designable_residues( pose, this_nodes_resfile );
 
-	/// Make a string out of the AAs at my designable positions in the current state
-	std::string my_sequence = get_designable_sequence ( pose, my_designable_residues );
+	/// Make a string vector out of the AAs at my designable positions in the current state
+    utility::vector1< std::string > my_sequence = get_designable_sequence ( pose, my_designable_residues );
+
+    std::string pass_seq;
+    for( const std::string& resi_base_name: my_sequence) {
+        pass_seq += resi_base_name + " ";
+    }
 
 	/// Get the AAs at designable positions of the other states I need to cooperate with
-	utility::vector1<std::string> other_pose_sequences( n_procs );
+	utility::vector1< utility::vector1<std::string> > other_pose_sequences( n_procs );
 	for ( core::Size ii = 1; ii <= n_procs; ++ii ) {
 		if ( rank == ii ) {
 			for ( core::Size jj = 1; jj <= n_procs; ++jj ) {
-				if ( rank!=jj ) utility::send_string_to_node( jj-1, my_sequence ); // node ranks are 0-indexed
+				if ( rank!=jj ) utility::send_string_to_node( jj-1, pass_seq ); // node ranks are 0-indexed
 				else other_pose_sequences[jj] = my_sequence;
 			}
 		} else {
-			other_pose_sequences[ii] = utility::receive_string_from_node( ii-1 ); // node ranks are 0-indexed
+            std::string passed_seq;
+			passed_seq = utility::receive_string_from_node( ii-1 ); // node ranks are 0-indexed
+            //Need to split passed_seq by spaces
+            std::istringstream iss(passed_seq);
+            std::string resi_name;
+            while (iss >> resi_name){
+                other_pose_sequences[ii].push_back(resi_name);
+            }
 		}
 	}
 
 	/// Let the master make sure all the sequences are the same length,
 	/// i.e. all the states have same number of designable residues
 	if ( master ) {
-		for ( std::string const & sequence: other_pose_sequences ) {
+		for ( utility::vector1< std::string > const & sequence: other_pose_sequences ) {
 			if ( sequence.size() != my_sequence.size() ) {
 				utility_exit_with_message( "Error: all states must have the same number of designable residues" );
 			}
 		}
 	}
 
 	if ( debug_ ) {
-		TR << "my sequence is " << my_sequence << std::endl;
+		TR << "my sequence is ";
+        for( core::Size i = 1; i <= my_sequence.size(); ++i ) {
+            TR << my_sequence[i] << " ";
+        }
+        TR << std::endl;
 	}
 
 	/// Run my multistate design

source/src/protocols/recon_design/MSDMover.hh

@@ -73,7 +73,7 @@ public:
 	/// encourage poses to adopt same sequence. Returns a COP of the constraints that were
 	/// applied so they can be removed later
 	utility::vector1< core::scoring::constraints::ConstraintCOP >
-	apply_linked_constraints( Pose & pose, utility::vector1< std::string > other_pose_sequences,
+	apply_linked_constraints( Pose & pose, utility::vector1< utility::vector1< std::string > > other_pose_sequences,
 		utility::vector1< core::Size > my_designable_residues );
 
 	/// @brief Initialize mover by checking that input poses were passed correctly,

source/src/protocols/recon_design/recon_util.cc

@@ -45,31 +45,40 @@ get_designable_residues( core::pose::Pose & pose, std::string resfile ) {
 
 /// Based on a pose and the indices of all designable residues, get a
 /// string of all designable AAs concatenated
-std::string
+utility::vector1< std::string >
 get_designable_sequence ( core::pose::Pose & pose, utility::vector1< core::Size > designable_residues ) {
-	std::string sequence = "";
+	//std::string sequence = "";
+    utility::vector1< std::string > sequence;
 	for ( core::Size seqpos: designable_residues ) {
-		sequence += pose.residue( seqpos ).name1();
+		//sequence += pose.residue( seqpos ).name1();
+        sequence.push_back(pose.residue_type( seqpos ).base_name());
 	}
 	return sequence;
 }
 
 /// Based on a list of sequences from poses, get all the AAs present at
 /// the position given by position_no
 utility::vector1< std::string >
-get_candidate_AAs( utility::vector1< std::string > other_pose_sequences,
+get_candidate_AAs( utility::vector1< utility::vector1< std::string > > other_pose_sequences,
 	core::Size position_no ) {
 
 	utility::vector1< std::string > candidate_AAs;
 
 	// Iterate through all poses in my collection and find their AA at this position
-	for ( std::string const & sequence: other_pose_sequences ) {
-		char current_AA = sequence[ position_no-1 ]; // string is zero indexed
-		std::string current_AA_3letter = core::chemical::name_from_aa( core::chemical::aa_from_oneletter_code( current_AA ) );
-		if ( std::find( candidate_AAs.begin(), candidate_AAs.end(), current_AA_3letter ) == candidate_AAs.end() ) {
-			candidate_AAs.push_back( current_AA_3letter );
-		}
-	}
+	//for ( std::string const & sequence: other_pose_sequences ) {
+	//	char current_AA = sequence[ position_no-1 ]; // string is zero indexed
+	//	std::string current_AA_3letter = core::chemical::name_from_aa( core::chemical::aa_from_oneletter_code( current_AA ) );
+	//	if ( std::find( candidate_AAs.begin(), candidate_AAs.end(), current_AA_3letter ) == candidate_AAs.end() ) {
+	//		candidate_AAs.push_back( current_AA_3letter );
+	//	}
+	//}
+    for ( utility::vector1< std::string > const & sequence: other_pose_sequences ) {
+        std::string current_AA = sequence[ position_no ]; // vector1 is 1 indexed
+        //std::string current_AA_3letter = core::chemical::name_from_aa( core::chemical::aa_from_oneletter_code( current_AA ) );
+        if ( std::find( candidate_AAs.begin(), candidate_AAs.end(), current_AA ) == candidate_AAs.end() ) {
+            candidate_AAs.push_back( current_AA );
+        }
+    }
 
 	return candidate_AAs;
 }

source/src/protocols/recon_design/recon_util.hh

@@ -32,13 +32,13 @@ get_designable_residues( core::pose::Pose & pose, std::string resfile );
 
 /// @brief Based on a pose and the indices of all designable residues, get a
 /// string of all designable AAs concatenated
-std::string
+utility::vector1< std::string >
 get_designable_sequence ( core::pose::Pose & pose, utility::vector1< core::Size > designable_residues );
 
 /// @brief Based on a list of sequences from poses, get all the AAs present at
 /// the position given by position_no
 utility::vector1< std::string > get_candidate_AAs(
-	utility::vector1< std::string > other_pose_sequences,
+	utility::vector1< utility::vector1< std::string > > other_pose_sequences,
 	core::Size position_no );
 
 /// @brief Given a list of poses, find the index of a particular pose

source/test/protocols/recon_design/MSDMoverTest.cxxtest.hh

@@ -128,10 +128,19 @@ public:
 		msd_mover->setup_mover( *pose1 );
 
 		// Find the sequence of the other poses
-		utility::vector1< std::string > other_pose_sequences;
-		other_pose_sequences.push_back( "A" );
-		other_pose_sequences.push_back( "C" );
-		other_pose_sequences.push_back( "D" );
+		utility::vector1< utility::vector1< std::string > > other_pose_sequences;
+		//other_pose_sequences.push_back( "A" );
+		//other_pose_sequences.push_back( "C" );
+		//other_pose_sequences.push_back( "D" );
+        utility::vector1< std::string > ala_seq;
+        ala_seq.push_back("ALA");
+        other_pose_sequences.push_back(ala_seq);
+        utility::vector1< std::string > cys_seq;
+        cys_seq.push_back("CYS");
+        other_pose_sequences.push_back(cys_seq);
+        utility::vector1< std::string > asp_seq;
+        asp_seq.push_back("ASP");
+        other_pose_sequences.push_back(asp_seq);
 
 		utility::vector1< utility::vector1< core::Size > > designable_residues = msd_mover->designable_residues();