mmtbx.building package

Submodules

mmtbx.building.extend_sidechains module

class mmtbx.building.extend_sidechains.conformation_scorer(old_residue, new_residue)

Bases: object

Stand-in for the conformation scoring class in mmtbx.refinement.real_space; instead of calculating fit to the map, this simply uses the change in position of the first atom being moved at each rotation. This allows us to superimpose the conformations for those atoms which are present in both the old and the new residues.

apply_final()
reset(sites_cart, selection)
update(sites_cart, selection)
mmtbx.building.extend_sidechains.correct_sequence(pdb_hierarchy, sequences, truncate_to_cbeta=False, out=<_io.TextIOWrapper name='<stdout>' mode='w' encoding='utf-8'>)

Modify the sequence for the pdb hierarchy to match that of the aligned sequence. This will remove incompatible atoms; the sidechains will still need to be extended separated. For proteins only - mismatches in nucleic acids will only result in a warning.

Parameters:
  • pdb_hierarchy – iotbx.pdb.hierarchy.root object

  • sequences – list of iotbx.bioinformatics.sequence objects

  • trucate_to_cbeta – chop off entire sidechain to C-beta (default: leave common atoms in place)

  • out – output filehandle (default = stdout)

Returns:

number of atom_group objects renamed

class mmtbx.building.extend_sidechains.extend_and_refine(pdb_hierarchy, xray_structure, fmodel, params, cif_objects=(), out=<_io.TextIOWrapper name='<stdout>' mode='w' encoding='utf-8'>, output_model=None, output_map_coeffs=None, prefix=None, write_files=True, reset_segid=True, verbose=True)

Bases: object

Run the combined sidechain extension and real-space fitting, and optionally write final model and map coefficients.

mmtbx.building.extend_sidechains.extend_protein_model(pdb_hierarchy, mon_lib_srv, add_hydrogens=None, selection=None)

Rebuild a sidechain by substituting an ideal amino acid and rotating the sidechain to match the old conformation as closely as possible. Limited functionality:

  1. Amino-acids only,

  2. side chain atoms only.

  3. Not terminii aware.

  4. Not aware of v2.3 vs v3.2 atom names e.g. HB1,HB2 vs HB2,HB3.

  5. Skips altlocs.

mmtbx.building.extend_sidechains.extend_residue(residue, target_atom_group, mon_lib_srv)

Rebuild a sidechain by substituting an ideal amino acid and rotating the sidechain to match the old conformation as closely as possible. Limited functionality:

  1. Amino-acids only, 2) side chain atoms only.

class mmtbx.building.extend_sidechains.prefilter(fmodel, out, backbone_min_sigma=1.0)

Bases: object

Optional filter for excluding residues with poor backbone density from being extended. This is done as a separate callback to enable the main rebuilding routine to be independent of data/maps.

mmtbx.building.extend_sidechains.refit_residues(pdb_hierarchy, cif_objects, fmodel, use_rotamers=True, anneal=False, verbose=True, allow_modified_residues=False, out=<_io.TextIOWrapper name='<stdout>' mode='w' encoding='utf-8'>)

Use real-space refinement tools to fit newly extended residues.

mmtbx.building.make_library module

Tools for assembling an ensemble of related structures for local rebuilding by homology.

class mmtbx.building.make_library.ensembler(reference_hierarchy, related_chains, atom_selection_string=None, backbone_only=<libtbx.AutoType object>, calpha_only=None, nproc=<libtbx.AutoType object>, sieve_fit=False, frac_discard=0.5, log=<libtbx.utils.null_out object>)

Bases: object

Superpose a collection of single-chain models on a reference chain. This is intentionally very limited in function, but flexible with respect to the choice of atoms to superpose.

align_model(i_model)
as_multi_model_hierarchy()
realign(atom_selection_string)
mmtbx.building.make_library.extract_and_superpose(reference_hierarchy, search_directory, sequence, params, nproc=<libtbx.AutoType object>, out=<_io.TextIOWrapper name='<stdout>' mode='w' encoding='utf-8'>)
mmtbx.building.make_library.extract_peptide_fragments_by_sequence(pdb_hierarchy, sequence, renumber_from=None, mainchain_only=False, remove_hydrogens=True, reset_chain_id='A', reset_icode=' ')

Bases: object

Multiprocessing wrapper for pulling matching chains out of a collection of PDB hierarchies (either from local files or from a PDB mirror).

mmtbx.building.make_library.fetch_similar_pdb_ids(sequence, min_identity=0.95, min_resolution=3.0, xray_only=True, data_only=False, expect=1e-07, sort_by_resolution=True, ligands=None, log=None)

Finds closely related structures in the RCSB database, optionally limited to structures with specific ligands.

mmtbx.building.make_library.find_similar_chains(pdb_hierarchy, sequence, source_info=None, first_chain_only=False, remove_alt_confs=True, reset_chain_id=None, min_identity=0.95, min_identity_epsilon=0.02, log=<libtbx.utils.null_out object>)

Find all chains with at least the specified fractional sequence identity to the target, and extract as related_chain objects with new pdb hierarchies.

mmtbx.building.make_library.join_hierarchies(hierarchies, model_ids=None)
mmtbx.building.make_library.load_all_models_in_directory(dir_name, limit_extensions=True, recursive=False)

Load all models in the specified directory, returning a list of file names and iotbx.file_reader objects.

mmtbx.building.make_library.load_pdb_models(pdb_ids, fault_tolerant=False, log=<libtbx.utils.null_out object>)

Given a list of PDB IDs, load hierarchies for all of them (ignoring any unknown atoms).

Combines load_pdb_models and extract_related_models, with multiprocessing support. Functionally equivalent to load_related_models_in_directory but with models from the PDB (local or remote).

Combines load_all_models_in_directory and extract_related_models.

class mmtbx.building.make_library.related_chain(**kwds)

Bases: slots_getstate_setstate_default_initializer

Container for a protein chain (as hierarchy object) and metadata, extracted based on similarity to a target sequence.

chain_id
identity
pdb_hierarchy
source_info
mmtbx.building.make_library.selection_by_symmetry(crystal_symmetries, source_infos, exclude_space_group=None, only_space_group=None, log=<libtbx.utils.null_out object>)

Return a list of indices for those crystal symmetries which either match the only_space_group parameter, or do not match exclude_space_group.

Module contents

mmtbx.building.atom_group_as_hierarchy(atom_group)

Convert an iotbx.pdb.hierarchy.atom_group object to a separate hierarchy (leaving the original hierarchy untouched), which allows us to use the pickling capabilities of iotbx.pdb.hierarchy.root.

class mmtbx.building.box_build_refine_base(pdb_hierarchy, xray_structure, processed_pdb_file, target_map, selection, d_min, out, cif_objects=(), resolution_factor=0.25, selection_buffer_radius=5, box_cushion=2, target_map_rsr=None, geometry_restraints_manager=None, debug=False)

Bases: object

Base class for handling functions associated with rebuilding and refinement of atoms in a box, with the building methods implemented elsewhere.

anneal(simulated_annealing_params=None, start_temperature=None, cool_rate=None, number_of_steps=50)

Run real-space simulated annealing using the target map (not the RSR map, if this is different). In practice, the non-selection atoms in the box should almost always be restrained to their current positions, but the setup is left to the calling code.

cc_model_map(selection=None, radius=1.5)

Calculate the correlation coefficient for the current model (in terms of F(calc) from the xray structure) and the target map, calculated at atomic positions rather than grid points. This will be much less accurate than the CC calculated in the original crystal environment, with full F(model) including bulk solvent correction.

fit_residue_in_box(mon_lib_srv=None, rotamer_manager=None, backbone_sample_angle=20)
geometry_minimization(selection=None, bond=True, nonbonded=True, angle=True, dihedral=True, chirality=True, planarity=True, max_number_of_iterations=500, number_of_macro_cycles=5, rmsd_bonds_termination_cutoff=0, rmsd_angles_termination_cutoff=0)

Perform geometry minimization on a selection of boxed atoms.

get_selected_sites(selection=None, hydrogens=True)
mean_density_at_sites(selection=None)

Compute the mean density of the target map at coordinates for non-hydrogen atoms.

only_residue()

Returns the atom_group object corresponding to the original selection. Will raise an exception if the selection covers more than one atom_group.

real_space_refine(selection=None)

Run simple real-space refinement, returning the optimized weight for the map target (wx).

restrain_atoms(selection, reference_sigma, limit=1.0, reset_first=True, top_out_potential=False)

Apply harmonic reference restraints to the selected atoms, wiping out any previously existing harmonic restraints.

update_coordinates(sites_cart)

Set coordinates for the boxed xray_structure and pdb_hierarchy objects.

update_original_coordinates()

Copies the new coordinates from the box to the selected atoms in the original xray_structure object, and returns the (complete) list of coordinates.

update_sites_from_pdb_atoms()

Update the xray_structure coordinates after manipulating the PDB hierarchy.

mmtbx.building.extract_iselection(pdb_objects)
mmtbx.building.generate_sidechain_clusters(residue, mon_lib_srv)

Extract Chi angle indices (including rotation axis) from the atom_group

mmtbx.building.get_difference_maps(fmodel, resolution_factor=0.25)
mmtbx.building.get_model_map_stats(selection, target_map, model_map, unit_cell, sites_cart, pdb_atoms, local_sampling=False)

Collect basic statistics for a model map and some target map (usually an mFo-DFc map), including CC, mean, and minimum density at the atomic positions.

mmtbx.building.get_nearby_water_selection(pdb_hierarchy, xray_structure, selection, selection_buffer_radius=5)
mmtbx.building.get_non_hydrogen_atom_indices(pdb_object)
mmtbx.building.get_window_around_residue(residue, window_size=2)
mmtbx.building.is_stub_residue(atom_group)
mmtbx.building.iter_residue_groups(hierarchy)
class mmtbx.building.local_density_quality(fofc_map, two_fofc_map, sites_cart=None, atom_names=None, unit_cell=None, atom_selection=None, xray_structure=None, radius=2.5)

Bases: object

atoms_below_fofc_map_level(sigma_cutoff=-2.5)
atoms_below_two_fofc_map_level(sigma_cutoff=-2.5)
atoms_outside_density(**kwds)
density_at_atom(atom_name)

Return a simple object containing density levels for the named atom, or None if no match was found.

fraction_of_nearby_grid_points_above_cutoff(sigma_cutoff=2.5)
max_fofc_value()
number_of_atoms_below_fofc_map_level(*args, **kwds)
number_of_atoms_below_two_fofc_map_level(*args, **kwds)
number_of_atoms_outside_density(*args, **kwds)
show_atoms_outside_density(**kwds)
mmtbx.building.remove_sidechain_atoms(pdb_objects)
mmtbx.building.reprocess_pdb(pdb_hierarchy, crystal_symmetry, cif_objects, out)
mmtbx.building.residue_density_quality(atom_group, unit_cell, two_fofc_map, fofc_map)
mmtbx.building.residues_are_adjacent(residue1, residue2, max_sep=2.5)
mmtbx.building.run_real_space_annealing(xray_structure, pdb_hierarchy, selection, target_map, d_min, processed_pdb_file=None, cif_objects=(), resolution_factor=0.25, params=None, wc=1, target_map_rsr=None, rsr_after_anneal=False, reference_sigma=0.5, out=None, debug=False)

Run simulated annealing with a real-space target. For maximum flexibility, a separate map may be used for the initial real-space refinement step.

mmtbx.building.show_chain_resseq_ranges(residues, out=<_io.TextIOWrapper name='<stdout>' mode='w' encoding='utf-8'>, prefix='')

Given a list of residues (either residue_group or atom_group objects), print a summary of the chains and residue ranges present.