# Example for: model.mutate()
# This will read a PDB file, change its sequence a little, build new
# coordinates for any of the additional atoms using only the internal
# geometry, and write the mutant PDB file. It can be seen as primitive,
# but rapid comparative modeling for substitution mutants. For insertion
# and deletion mutants, follow the standard comparative modeling procedure.
env = environ()
# Read the topology library with non-hydrogen atoms only:
env.libs.topology.read(file='$(LIB)/top_heav.lib')
# To produce a mutant with all hydrogens, uncomment this line:
#env.libs.topology.read(file='$(LIB)/top_allh.lib')
# Read the CHARMM parameter library:
env.libs.parameters.read(file='$(LIB)/par.lib')
# Read the original PDB file and copy its sequence to the alignment array:
code = '1fas'
aln = alignment(env)
mdl = model(env, file=code)
aln.append_model(mdl, atom_files=code, align_codes=code)
# Select the residues to be mutated: in this case all ASP residues:
mdl.pick_atoms(aln, res_types='ASP')
# The second example is commented out; it selects residues '1' and '10'.
#mdl.pick_atoms(aln, selection_segment=('1', '1'),
# selection_status='INITIALIZE', selection_search='SEGMENT',
# selection_from='ALL')
#mdl.pick_atoms(aln, selection_segment=('10', '10'), selection_status='ADD',
# selection_search='SEGMENT', selection_from='ALL')
# Mutate the selected residues into HIS residues (neutral HIS):
mdl.mutate(residue_type='HIS')
# Add the mutated sequence to the alignment arrays (it is now the second
# sequence in the alignment):
aln.append_model(mdl, align_codes='1fas-1')
# Generate molecular topology for the mutant:
mdl.generate_topology(aln, sequence='1fas-1')
# Transfer all the coordinates you can from the template native structure
# to the mutant (this works even if the order of atoms in the native PDB
# file is not standard):
mdl.transfer_xyz(aln)
# Build the remaining unknown coordinates for the mutant:
mdl.build(initialize_xyz=False, build_method='INTERNAL_COORDINATES')
# Write the mutant to a file:
mdl.write(file='1fas-1.atm')