IMP Reference Guide
2.12.0
The Integrative Modeling Platform
|
Fitting atomic structures into a cryo-electron microscopy density map. More...
Fitting atomic structures into a cryo-electron microscopy density map.
MultiFit is a computational method for simultaneously fitting atomic structures of components into their assembly density map at resolutions as low as 25Å. The component positions and orientations are optimized with respect to a scoring function that includes the quality-of-fit of components in the map, the protrusion of components from the map envelope, as well as the shape complementarity between pairs of components. The scoring function is optimized by an exact inference optimizer DOMINO that efficiently finds the global minimum in a discrete sampling space.
For more information please see the MultiFit website.
See the IMP::cnmultifit module for a similar protocol for handling symmetric complexes.
Generally, this module is not used directly; instead, the multifit
command line tool is used. For an example, see Modeling of 3sfd.
A webserver that uses both this module and the IMP::cnmultifit module is also available.
Build an anchor graph from a PDB file.
Author(s): Keren Lasker
Maintainer: benmwebb
License: GPL This library is free software; you can redistribute it and/or modify it under the terms of the GNU General Public License as published by the Free Software Foundation; either version 2 of the License, or (at your option) any later version.
Publications:
Classes | |
class | AnchorsData |
Storage of anchors (points and edges) More... | |
class | AssemblyHeader |
Holds data about the assembly density needed for optimization. More... | |
class | ComplementarityRestraint |
Compute the complementarity between two molecules. More... | |
class | ComponentHeader |
Holds data about a component needed for optimization. More... | |
class | DataPointsAssignment |
class | DensityDataPoints |
Stores density voxels as a vector of Array1D. More... | |
class | DummyRestraint |
A simple Restraint that always returns a score of zero. More... | |
class | Ensemble |
An ensemble of fitting solutions. More... | |
class | FFTFitting |
Fit a molecule inside its density by local or global FFT. More... | |
class | FFTFittingOutput |
Storage of the results from an FFT fit. More... | |
class | FittingSolutionRecord |
A fitting solution record. More... | |
class | FittingStates |
class | GeometricHash |
Geometric Hash table. More... | |
class | MergeTreeBuilder |
Utility class for building merge trees. More... | |
class | ProbabilisticAnchorGraph |
Probabilistic anchor graph. More... | |
class | ProteinsAnchorsSamplingSpace |
stores the anchors sampling space for each protein More... | |
class | ProteomicsData |
Storage of proteomics data. More... | |
class | ProteomicsEMAlignmentAtomic |
Align proteomics graph to EM density map. More... | |
class | RadiusOfGyrationRestraint |
Ensure the radius of gyration of particles fits the predicted one. More... | |
class | RigidLeavesRefiner |
Return all rigid body members that are also hierarchy leaves. More... | |
class | SettingsData |
Holds header data for optimization. More... | |
class | WeightedExcludedVolumeRestraint |
Calculate score based on fit to EM map. More... | |
Functions | |
void | add_surface_index (core::Hierarchy mhd, Float apix, FloatKey shell_key=FloatKey("surf_ind"), FloatKey radius_key=core::XYZR::get_radius_key(), FloatKey weight_key=atom::Mass::get_mass_key()) |
Add shell index to leaves. More... | |
FittingSolutionRecords | convert_em_to_multifit_format (const em::FittingSolutions &em_fits) |
algebra::Transformation3Ds | convert_multifit_format_to_transformations (const FittingSolutionRecords &recs) |
em::FittingSolutions | convert_multifit_to_em_format (const FittingSolutionRecords &multifit_fits) |
FittingSolutionRecords | convert_transformations_to_multifit_format (const algebra::Transformation3Ds &trans) |
atom::Hierarchy | create_coarse_molecule_from_density (em::DensityMap *dmap, float dens_threshold, int num_beads, Model *mdl, float bead_radius) |
Coarsen a density map based on voxels clustering. More... | |
em::DensityMap * | create_hit_map (core::RigidBody rb, Refiner *rb_ref, const FittingSolutionRecords &sols, em::DensityMap *damp) |
IMP::Restraint * | create_weighted_excluded_volume_restraint (core::RigidBody rb1, core::RigidBody rb2, FloatKey shell_key=FloatKey("surf_ind")) |
Create a weighted excluded volume restraint between two rigid bodies. More... | |
multifit::FittingSolutionRecords | fft_based_rigid_fitting (atom::Hierarchy mol2fit, em::DensityMap *dmap, double density_threshold, double angle_sampling_interval_rad) |
FFT fit of a molecule in the density. More... | |
float | get_actual_radius_of_gyration (ParticlesTemp ps) |
IntsList | get_anchor_indices_matching_secondary_structure (const AnchorsData &ad, const atom::SecondaryStructureResidues &ssrs, Float max_rmsd=0.7) |
Get lists of anchors that match a sequence of secondary structures. More... | |
void | get_anchors_for_density (em::DensityMap *dmap, int number_of_means, float density_threshold, std::string pdb_filename, std::string cmm_filename, std::string seg_filename, std::string txt_filename) |
Generate anchors in several formats for a given density map. More... | |
float | get_approximated_radius (int len) |
float | get_approximated_radius_of_gyration (int len) |
FittingSolutionRecords | get_close_to_point (const FittingSolutionRecords &fit_sols, atom::Hierarchy mh, IMP::Particle *ap, Float dist) |
prune solutions by distance to an anchor point More... | |
IntsList | get_connected_components (em::DensityMap *dmap, float threshold, float edge_threshold) |
Return connected components based on density values. More... | |
ProteinsAnchorsSamplingSpace | get_part_of_sampling_space (const ProteinsAnchorsSamplingSpace &prots_ss, const Strings &prot_names) |
Get the sampling space of few of the proteins. More... | |
multifit::SettingsData * | get_partial_assembly_setting_data (multifit::SettingsData *prots_sd, const Strings &prot_names) |
Get the assembly data for a few of the proteins. More... | |
ProteomicsData * | get_partial_proteomics_data (const ProteomicsData *pd, const Strings &prot_names) |
algebra::Vector3Ds | get_points_close_to_molecule (const atom::Hierarchy &mh, const algebra::Vector3Ds points, Float max_dist) |
algebra::Vector3Ds | get_points_far_from_molecule (const atom::Hierarchy &mh, const algebra::Vector3Ds points, Float max_dist) |
algebra::Vector3D | get_segment_maximum (const DataPointsAssignment &dpa, em::DensityMap *dmap, int segment_id) |
algebra::Vector3D | get_segment_maximum (const DataPointsAssignment &dpa, DensGrid *dmap, int segment_id) |
void | get_segmentation (em::DensityMap *dmap, double apix, double density_threshold, int num_means, const std::string pdb_filename, const std::string cmm_filename, const std::string seg_filename, const std::string txt_filename) |
Segment a density map using the anchor graph. More... | |
em::DensityMap * | grid2map (const DensGrid &dg, float spacing) |
Ensemble * | load_ensemble (multifit::SettingsData *sd, Model *mdl, const ProteinsAnchorsSamplingSpace &mapping_data) |
AnchorsData | molecule2anchors (atom::Hierarchy mh, int k) |
em::FittingSolutions | pca_based_rigid_fitting (core::RigidBody rb, Refiner *rb_refiner, em::DensityMap *em_map, Float threshold, FloatKey wei_key=atom::Mass::get_mass_key(), algebra::PrincipalComponentAnalysis dens_pca_input=algebra::PrincipalComponentAnalysis()) |
Compute fitting scores for a given set of rigid transformations. More... | |
em::FittingSolutions | pca_based_rigid_fitting (ParticlesTemp ps, em::DensityMap *em_map, Float threshold, FloatKey wei_key=atom::Mass::get_mass_key(), algebra::PrincipalComponentAnalysis dens_pca_input=algebra::PrincipalComponentAnalysis()) |
Compute fitting scores for a given set of rigid transformations. More... | |
AnchorsData | read_anchors_data (const char *txt_filename) |
std::pair< algebra::Vector3Ds, CEdges > | read_cmm (const std::string &cmm_filename) |
FittingSolutionRecords | read_fitting_solutions (const char *fitting_fn) |
Fitting solutions reader. More... | |
IntsList | read_paths (const char *txt_filename, int max_paths=INT_MAX) |
Read paths. More... | |
ProteinsAnchorsSamplingSpace | read_protein_anchors_mapping (multifit::ProteomicsData *prots, const std::string &anchors_prot_map_fn, int max_paths=INT_MAX) |
ProteomicsData * | read_proteomics_data (const char *proteomics_fn) |
Proteomics reader. More... | |
SettingsData * | read_settings (const char *filename) |
em::DensityMap * | remove_background (em::DensityMap *dmap, float threshold, float edge_threshold) |
Returns a map containing all density without the background. More... | |
void | write_chimera (const std::string &chimera_filename, const DataPointsAssignment &dpa) |
void | write_cmm (const std::string &cmm_filename, const std::string &marker_set_name, const AnchorsData &dpa) |
void | write_connolly_surface (atom::Atoms atoms, TextOutput fn, float density, float probe_radius) |
Write the Connolly surface for a set of atoms to a file. More... | |
void | write_fitting_solutions (const char *fitting_fn, const FittingSolutionRecords &fit_sols, int num_sols=-1) |
Write fitting solutions to a file. More... | |
void | write_markers (const algebra::PrincipalComponentAnalysisD< 3 > &pca, std::ostream &out) |
void | write_paths (const IntsList &paths, const std::string &txt_filename) |
void | write_pdb (const std::string &pdb_filename, const DataPointsAssignment &dpa) |
void | write_protein_anchors_mapping (const std::string &anchors_prot_map_fn, const std::string &anchors_fn, const std::vector< std::pair< String, String > > &prot_paths) |
void | write_protein_anchors_mapping (const std::string &anchors_prot_map_fn, const ProteinsAnchorsSamplingSpace &pa, const Strings &prot_names) |
void | write_segment_as_mrc (em::DensityMap *dmap, const DataPointsAssignment &dpa, int segment_id, Float resolution, Float apix, const std::string &filename) |
void | write_segment_as_pdb (const DataPointsAssignment &dpa, int segment_id, const std::string &filename) |
void | write_segments_as_mrc (em::DensityMap *dmap, const DataPointsAssignment &dpa, Float resolution, Float apix, Float threshold, const std::string &filename) |
Write segments in MRC format. More... | |
void | write_segments_as_pdb (const DataPointsAssignment &dpa, const std::string &filename) |
void | write_settings (const char *filename, const SettingsData *sd) |
void | write_txt (const std::string &txt_filename, const AnchorsData &ad) |
Standard module functions | |
All | |
std::string | get_module_version () |
std::string | get_module_name () |
std::string | get_data_path (std::string file_name) |
Return the full path to one of this module's data files. More... | |
std::string | get_example_path (std::string file_name) |
Return the full path to one of this module's example files. More... | |
typedef IMP::Vector< AlignmentParams > IMP::multifit::AlignmentParamsList |
Pass or store a set of AlignmentParams .
Definition at line 322 of file AlignmentParams.h.
typedef IMP::Vector< AnchorsData > IMP::multifit::AnchorsDataList |
Pass or store a set of AnchorsData .
Definition at line 89 of file anchors_reader.h.
typedef IMP::Vector< ComplementarityParams > IMP::multifit::ComplementarityParamsList |
Pass or store a set of ComplementarityParams .
Definition at line 269 of file AlignmentParams.h.
A vector of reference-counting object pointers.
Definition at line 73 of file SettingsData.h.
A vector of weak (non reference-counting) pointers to specified objects.
Definition at line 73 of file SettingsData.h.
typedef IMP::Vector< ConnectivityParams > IMP::multifit::ConnectivityParamsList |
Pass or store a set of ConnectivityParams .
Definition at line 95 of file AlignmentParams.h.
A vector of reference-counting object pointers.
Definition at line 60 of file DensityDataPoints.h.
typedef IMP::Vector<IMP::WeakPointer< DensityDataPoints > > IMP::multifit::DensityDataPointsListTemp |
A vector of weak (non reference-counting) pointers to specified objects.
Definition at line 60 of file DensityDataPoints.h.
typedef IMP::Vector< DominoParams > IMP::multifit::DominoParamsList |
Pass or store a set of DominoParams .
Definition at line 48 of file AlignmentParams.h.
typedef IMP::Vector<IMP::Pointer< Ensemble > > IMP::multifit::Ensembles |
A vector of reference-counting object pointers.
Definition at line 49 of file ensemble_analysis.h.
typedef IMP::Vector<IMP::WeakPointer< Ensemble > > IMP::multifit::EnsemblesTemp |
A vector of weak (non reference-counting) pointers to specified objects.
Definition at line 49 of file ensemble_analysis.h.
typedef IMP::Vector< EVParams > IMP::multifit::EVParamsList |
Pass or store a set of EVParams .
Definition at line 186 of file AlignmentParams.h.
typedef IMP::Vector< FiltersParams > IMP::multifit::FiltersParamsList |
Pass or store a set of FiltersParams .
Definition at line 209 of file AlignmentParams.h.
typedef IMP::Vector< FittingParams > IMP::multifit::FittingParamsList |
Pass or store a set of FittingParams .
Definition at line 235 of file AlignmentParams.h.
Pass or store a set of FittingSolutionRecord .
Definition at line 84 of file FittingSolutionRecord.h.
typedef IMP::Vector< FragmentsParams > IMP::multifit::FragmentsParamsList |
Pass or store a set of FragmentsParams .
Definition at line 123 of file AlignmentParams.h.
typedef IMP::Vector<IMP::Pointer< ProbabilisticAnchorGraph > > IMP::multifit::ProbabilisticAnchorGraphs |
A vector of reference-counting object pointers.
Definition at line 67 of file anchor_graph.h.
typedef IMP::Vector<IMP::WeakPointer< ProbabilisticAnchorGraph > > IMP::multifit::ProbabilisticAnchorGraphsTemp |
A vector of weak (non reference-counting) pointers to specified objects.
Definition at line 67 of file anchor_graph.h.
typedef IMP::Vector< ProteinRecordData > IMP::multifit::ProteinRecordDataList |
Pass or store a set of ProteinRecordData .
Definition at line 68 of file proteomics_reader.h.
Pass or store a set of ProteinsAnchorsSamplingSpace .
Definition at line 65 of file protein_anchors_mapping_reader.h.
typedef IMP::Vector<IMP::Pointer< ProteomicsEMAlignmentAtomic > > IMP::multifit::ProteomicsEMAlignmentAtomics |
A vector of reference-counting object pointers.
Definition at line 115 of file proteomics_em_alignment_atomic.h.
typedef IMP::Vector<IMP::WeakPointer< ProteomicsEMAlignmentAtomic > > IMP::multifit::ProteomicsEMAlignmentAtomicsTemp |
A vector of weak (non reference-counting) pointers to specified objects.
Definition at line 115 of file proteomics_em_alignment_atomic.h.
typedef IMP::Vector< RogParams > IMP::multifit::RogParamsList |
Pass or store a set of RogParams .
Definition at line 143 of file AlignmentParams.h.
typedef IMP::Vector< XlinkParams > IMP::multifit::XlinkParamsList |
Pass or store a set of XlinkParams .
Definition at line 72 of file AlignmentParams.h.
void IMP::multifit::add_surface_index | ( | core::Hierarchy | mhd, |
Float | apix, | ||
FloatKey | shell_key = FloatKey("surf_ind") , |
||
FloatKey | radius_key = core::XYZR::get_radius_key() , |
||
FloatKey | weight_key = atom::Mass::get_mass_key() |
||
) |
Add shell index to leaves.
[in] | mhd | the Hierarchy to modify |
[in] | apix | the resolution of the surface |
[in] | shell_key | the Key used to store the shell index |
[in] | radius_key | the Key used for particle radius |
[in] | weight_key | the Key used for particle weight (e.g. mass) |
atom::Hierarchy IMP::multifit::create_coarse_molecule_from_density | ( | em::DensityMap * | dmap, |
float | dens_threshold, | ||
int | num_beads, | ||
Model * | mdl, | ||
float | bead_radius | ||
) |
Coarsen a density map based on voxels clustering.
[in] | dmap | the density map to coarsen |
[in] | dens_threshold | use only voxels above this threshold for clustering |
[in] | num_beads | the number of beads |
[in] | mdl | model to add the new molecule to |
[in] | bead_radius | bead radius |
IMP::Restraint* IMP::multifit::create_weighted_excluded_volume_restraint | ( | core::RigidBody | rb1, |
core::RigidBody | rb2, | ||
FloatKey | shell_key = FloatKey("surf_ind") |
||
) |
Create a weighted excluded volume restraint between two rigid bodies.
[in] | rb1 | the first rigid body |
[in] | rb2 | the second rigid body |
[in] | shell_key | the attribute that stored the particles surface level with respect to its molecule |
multifit::FittingSolutionRecords IMP::multifit::fft_based_rigid_fitting | ( | atom::Hierarchy | mol2fit, |
em::DensityMap * | dmap, | ||
double | density_threshold, | ||
double | angle_sampling_interval_rad | ||
) |
FFT fit of a molecule in the density.
[in] | mol2fit | a rigid body molecule to fit |
[in] | dmap | the map to fit into |
[in] | density_threshold | ignore density below this value |
[in] | angle_sampling_interval_rad | sampling internal in radians |
IntsList IMP::multifit::get_anchor_indices_matching_secondary_structure | ( | const AnchorsData & | ad, |
const atom::SecondaryStructureResidues & | ssrs, | ||
Float | max_rmsd = 0.7 |
||
) |
Get lists of anchors that match a sequence of secondary structures.
[in] | ad | The AnchorsData |
[in] | ssrs | The SecondaryStructureResidue particles to match |
[in] | max_rmsd | SecondaryStructureResidue match must be below this value (0.816 is rmsd of known SSE to random) |
void IMP::multifit::get_anchors_for_density | ( | em::DensityMap * | dmap, |
int | number_of_means, | ||
float | density_threshold, | ||
std::string | pdb_filename, | ||
std::string | cmm_filename, | ||
std::string | seg_filename, | ||
std::string | txt_filename | ||
) |
Generate anchors in several formats for a given density map.
FittingSolutionRecords IMP::multifit::get_close_to_point | ( | const FittingSolutionRecords & | fit_sols, |
atom::Hierarchy | mh, | ||
IMP::Particle * | ap, | ||
Float | dist | ||
) |
prune solutions by distance to an anchor point
[in] | fit_sols | initial fitting solutions |
[in] | mh | the molecule the fitting solutions apply for |
[in] | ap | anchor point for which the transformed mh (fit) should be close to |
[in] | dist | all fits such that the distance between ap and the fit center is smaller than dist will be included |
IntsList IMP::multifit::get_connected_components | ( | em::DensityMap * | dmap, |
float | threshold, | ||
float | edge_threshold | ||
) |
Return connected components based on density values.
[in] | dmap | the density map to analyze |
[in] | threshold | consider only voxels above this threshold |
[in] | edge_threshold | an edge is added between two neighboring voxels if their density difference is below this threshold |
std::string IMP::multifit::get_data_path | ( | std::string | file_name | ) |
Return the full path to one of this module's data files.
To read the data file "data_library" that was placed in the data
directory of this module, do something like
This will ensure that the code works both when IMP is installed or if used via the setup_environment.sh
script.
std::string IMP::multifit::get_example_path | ( | std::string | file_name | ) |
Return the full path to one of this module's example files.
To read the example file "example_protein.pdb" that was placed in the examples
directory of this module, do something like
This will ensure that the code works both when IMP is installed or if used via the setup_environment.sh
script.
ProteinsAnchorsSamplingSpace IMP::multifit::get_part_of_sampling_space | ( | const ProteinsAnchorsSamplingSpace & | prots_ss, |
const Strings & | prot_names | ||
) |
Get the sampling space of few of the proteins.
multifit::SettingsData* IMP::multifit::get_partial_assembly_setting_data | ( | multifit::SettingsData * | prots_sd, |
const Strings & | prot_names | ||
) |
Get the assembly data for a few of the proteins.
algebra::Vector3Ds IMP::multifit::get_points_close_to_molecule | ( | const atom::Hierarchy & | mh, |
const algebra::Vector3Ds | points, | ||
Float | max_dist | ||
) |
Given a molecule and a set of points, return the indexes of the points that are close to the molecule (up to max_dist) and the res
void IMP::multifit::get_segmentation | ( | em::DensityMap * | dmap, |
double | apix, | ||
double | density_threshold, | ||
int | num_means, | ||
const std::string | pdb_filename, | ||
const std::string | cmm_filename, | ||
const std::string | seg_filename, | ||
const std::string | txt_filename | ||
) |
Segment a density map using the anchor graph.
All voxels above the threshold are segmented into the given number of clusters, and neighboring clusters are linked.
[in] | dmap | the density map to segment |
[in] | apix | the map spacing in angstroms per pixel |
[in] | density_threshold | consider only voxels over this threshold |
[in] | num_means | the number of segments to generate |
[in] | pdb_filename | write cluster centers as CA atoms in PDB format |
[in] | cmm_filename | if not empty, write clusters in CMM format |
[in] | seg_filename | if not empty, write segments in MRC format |
[in] | txt_filename | if not empty, write anchors in text format |
em::FittingSolutions IMP::multifit::pca_based_rigid_fitting | ( | core::RigidBody | rb, |
Refiner * | rb_refiner, | ||
em::DensityMap * | em_map, | ||
Float | threshold, | ||
FloatKey | wei_key = atom::Mass::get_mass_key() , |
||
algebra::PrincipalComponentAnalysis | dens_pca_input = algebra::PrincipalComponentAnalysis() |
||
) |
Compute fitting scores for a given set of rigid transformations.
Fit a protein to its density by principal component matching
[in] | rb | The rigid body to be fitted |
[in] | rb_refiner | The rigid body refiner |
[in] | em_map | The density map to fit to |
[in] | threshold | Use voxels above this threshold for PCA calculations |
[in] | wei_key | The weight key of the particles in the rigid body |
[in] | dens_pca_input | provide precalculated em_map PCA is available |
em::FittingSolutions IMP::multifit::pca_based_rigid_fitting | ( | ParticlesTemp | ps, |
em::DensityMap * | em_map, | ||
Float | threshold, | ||
FloatKey | wei_key = atom::Mass::get_mass_key() , |
||
algebra::PrincipalComponentAnalysis | dens_pca_input = algebra::PrincipalComponentAnalysis() |
||
) |
Compute fitting scores for a given set of rigid transformations.
Fit a protein to its density by principal component matching
[in] | ps | The particles to fit (treated as rigid) |
[in] | em_map | The density map to fit to |
[in] | threshold | Use voxels above this threshold for PCA calculations |
[in] | wei_key | The weight key of the particles in the rigid body |
[in] | dens_pca_input | precalculated em_map PCA |
FittingSolutionRecords IMP::multifit::read_fitting_solutions | ( | const char * | fitting_fn | ) |
Fitting solutions reader.
IntsList IMP::multifit::read_paths | ( | const char * | txt_filename, |
int | max_paths = INT_MAX |
||
) |
Read paths.
ProteomicsData* IMP::multifit::read_proteomics_data | ( | const char * | proteomics_fn | ) |
Proteomics reader.
em::DensityMap* IMP::multifit::remove_background | ( | em::DensityMap * | dmap, |
float | threshold, | ||
float | edge_threshold | ||
) |
Returns a map containing all density without the background.
[in] | dmap | the density map to segment |
[in] | threshold | consider only voxels above this threshold |
[in] | edge_threshold | consider only voxels above this threshold |
void IMP::multifit::write_connolly_surface | ( | atom::Atoms | atoms, |
TextOutput | fn, | ||
float | density, | ||
float | probe_radius | ||
) |
Write the Connolly surface for a set of atoms to a file.
void IMP::multifit::write_fitting_solutions | ( | const char * | fitting_fn, |
const FittingSolutionRecords & | fit_sols, | ||
int | num_sols = -1 |
||
) |
Write fitting solutions to a file.
[in] | fitting_fn | the fitting filename |
[in] | fit_sols | the fitting solutions to write to file |
[in] | num_sols | optional, only write the first num_sols fits. |
void IMP::multifit::write_segments_as_mrc | ( | em::DensityMap * | dmap, |
const DataPointsAssignment & | dpa, | ||
Float | resolution, | ||
Float | apix, | ||
Float | threshold, | ||
const std::string & | filename | ||
) |
Write segments in MRC format.