protomslib package

Submodules

Command

class protomslib.command.DualTopology(protein='protein.pdb', solutes=['solute1.pdb', 'solute2.pdb'], solvent='water.pdb', templates=['solute1.tem', 'solute2.tem'], nequil=5000000.0, nprod=40000000.0, dumpfreq=100000.0, ranseed=None, lambdaval=None, outfolder='out', restrained=[], softcore='mixed', spec_softcore='')

Bases: protomslib.command.ProteinLigandSimulation

This a command file for a dual-topology protein-ligand simulation Generates input for proteins, solutes and solvent write dual-toplogy parameters, lambda-replica exchange, dump results, pdb and restart files equilibrate and production run

Methods

clear() Clear the simulation description
readCommandFile(filename) Read in a simulation from a command file
run(exe[, cmdfile]) Run this simulation
setChunk(chunk) This adds the simulation chunk to ‘chunk’.
setDump(dump, freq) This adds the simulation dump to ‘dump’ with frequnecy freq
setForceField(filename) Set the filename for parfile
setParameter(parameter, value) Set the parameter to value
setProtein(n, filename) Set the filename for proteinN.
setSolute(n, filename) Set the filename for soluteN.
setSolvent(n, filename) Set the filename for solventN.
setStream(stream, filename) Set the direction of the output stream to a filename
writeCommandFile(filename) Write the contents of this simulation object to disc
make_softcore  
spec_string_from_user_input  
make_softcore(softcore, spec_softcore, templates, solutes)
spec_string_from_user_input(N)
class protomslib.command.Equilibration(protein='protein.pdb', solutes=['solute.pdb'], solvent='water.pdb', templates=['solute.tem'], outfolder='', ranseed=None, nsteps=100000.0, pdbfile='equilibrated.pdb')

Bases: protomslib.command.ProteinLigandSimulation

This a command file for equilibration of a protein-ligand system Generates input for proteins, solutes and solvent and write an equilibration and a pdb chunk

Methods

clear() Clear the simulation description
readCommandFile(filename) Read in a simulation from a command file
run(exe[, cmdfile]) Run this simulation
setChunk(chunk) This adds the simulation chunk to ‘chunk’.
setDump(dump, freq) This adds the simulation dump to ‘dump’ with frequnecy freq
setForceField(filename) Set the filename for parfile
setParameter(parameter, value) Set the parameter to value
setProtein(n, filename) Set the filename for proteinN.
setSolute(n, filename) Set the filename for soluteN.
setSolvent(n, filename) Set the filename for solventN.
setStream(stream, filename) Set the direction of the output stream to a filename
writeCommandFile(filename) Write the contents of this simulation object to disc
class protomslib.command.GCAPDual(protein='protein.pdb', solutes=['solute1.pdb', 'solute2.pdb'], solvent='water.pdb', templates=['solute1.tem', 'solute2.tem'], nequil=5000000.0, nprod=80000000.0, dumpfreq=100000.0, ranseed=None, lambdaval=None, outfolder='out', restrained=[], softcore='mixed', gcmcbox=None, gcmcwater='gcmc_water.pdb', adamval=None, watmodel='tip4p', spec_softcore='')

Bases: protomslib.command.ProteinLigandSimulation

This a command file for a dual-topology protein-ligand simulation Generates input for proteins, solutes and solvent write dual-toplogy parameters, lambda-replica exchange, dump results, pdb and restart files equilibrate and production run

Methods

clear() Clear the simulation description
readCommandFile(filename) Read in a simulation from a command file
run(exe[, cmdfile]) Run this simulation
setChunk(chunk) This adds the simulation chunk to ‘chunk’.
setDump(dump, freq) This adds the simulation dump to ‘dump’ with frequnecy freq
setForceField(filename) Set the filename for parfile
setParameter(parameter, value) Set the parameter to value
setProtein(n, filename) Set the filename for proteinN.
setSolute(n, filename) Set the filename for soluteN.
setSolvent(n, filename) Set the filename for solventN.
setStream(stream, filename) Set the direction of the output stream to a filename
writeCommandFile(filename) Write the contents of this simulation object to disc
make_softcore  
spec_string_from_user_input  
make_softcore(softcore, spec_softcore, templates, solutes)
spec_string_from_user_input(N)
class protomslib.command.GCAPSingle(protein='protein.pdb', solutes=['solute1.pdb'], solvent='water.pdb', templates=['solute.tem'], nequil=5000000.0, nprod=80000000.0, dumpfreq=100000.0, lambdaval=None, ranseed=None, outfolder='out', gcmcbox=None, gcmcwater='gcmc_water.pdb', adamval=None, watmodel='tip4p', restrained=[])

Bases: protomslib.command.ProteinLigandSimulation

This a command file for a single-topology protein-ligand simulation Generates input for proteins, solutes and solvent write lambda-replica exchange, dump results, pdb and restart files equilibrate and production run

Methods

clear() Clear the simulation description
readCommandFile(filename) Read in a simulation from a command file
run(exe[, cmdfile]) Run this simulation
setChunk(chunk) This adds the simulation chunk to ‘chunk’.
setDump(dump, freq) This adds the simulation dump to ‘dump’ with frequnecy freq
setForceField(filename) Set the filename for parfile
setParameter(parameter, value) Set the parameter to value
setProtein(n, filename) Set the filename for proteinN.
setSolute(n, filename) Set the filename for soluteN.
setSolvent(n, filename) Set the filename for solventN.
setStream(stream, filename) Set the direction of the output stream to a filename
writeCommandFile(filename) Write the contents of this simulation object to disc
class protomslib.command.GCMC(protein='protein.pdb', solutes=[], solvent='water.pdb', templates=[], gcmcwater='gcmc_water.pdb', gcmcbox=None, nequil=5000000.0, nprod=80000000.0, dumpfreq=200000.0, adamval=None, ranseed=None, watmodel='tip4p', outfolder='out_gcmc')

Bases: protomslib.command.ProteinLigandSimulation

This a command file for a GCMC simulation Generates input for proteins, solutes and solvent write GCMC parameters dump results, pdb and restart files equilibrate and production run

Methods

clear() Clear the simulation description
readCommandFile(filename) Read in a simulation from a command file
run(exe[, cmdfile]) Run this simulation
setChunk(chunk) This adds the simulation chunk to ‘chunk’.
setDump(dump, freq) This adds the simulation dump to ‘dump’ with frequnecy freq
setForceField(filename) Set the filename for parfile
setParameter(parameter, value) Set the parameter to value
setProtein(n, filename) Set the filename for proteinN.
setSolute(n, filename) Set the filename for soluteN.
setSolvent(n, filename) Set the filename for solventN.
setStream(stream, filename) Set the direction of the output stream to a filename
writeCommandFile(filename) Write the contents of this simulation object to disc
class protomslib.command.Jaws1(protein='protein.pdb', solutes=[], solvent='water.pdb', templates=[], outfolder='', ranseed=None, jawswater='jaws_water.pdb', jawsbox=None, nequil=5000000.0, nprod=40000000.0, watmodel='tip4p', dumpfreq=100000.0)

Bases: protomslib.command.ProteinLigandSimulation

This a command file for a JAWS-1 simulation Generates input for proteins, solutes and solvent write JAWS-1 parameters dump results, pdb and restart files equilibrate and production run

Methods

clear() Clear the simulation description
readCommandFile(filename) Read in a simulation from a command file
run(exe[, cmdfile]) Run this simulation
setChunk(chunk) This adds the simulation chunk to ‘chunk’.
setDump(dump, freq) This adds the simulation dump to ‘dump’ with frequnecy freq
setForceField(filename) Set the filename for parfile
setParameter(parameter, value) Set the parameter to value
setProtein(n, filename) Set the filename for proteinN.
setSolute(n, filename) Set the filename for soluteN.
setSolvent(n, filename) Set the filename for solventN.
setStream(stream, filename) Set the direction of the output stream to a filename
writeCommandFile(filename) Write the contents of this simulation object to disc
class protomslib.command.Jaws2(protein='protein.pdb', solutes=[], solvent='water.pdb', templates=[], outfolder='', ranseed=None, jawswater='jaws_water.pdb', jbias=6.5, nequil=5000000.0, nprod=40000000.0, watmodel='tip4p', dumpfreq=100000.0)

Bases: protomslib.command.ProteinLigandSimulation

This a command file for a JAWS-2 simulation Generates input for proteins, solutes and solvent write JAWS-2 parameters dump results, pdb and restart files equilibrate and production run

Methods

clear() Clear the simulation description
readCommandFile(filename) Read in a simulation from a command file
run(exe[, cmdfile]) Run this simulation
setChunk(chunk) This adds the simulation chunk to ‘chunk’.
setDump(dump, freq) This adds the simulation dump to ‘dump’ with frequnecy freq
setForceField(filename) Set the filename for parfile
setParameter(parameter, value) Set the parameter to value
setProtein(n, filename) Set the filename for proteinN.
setSolute(n, filename) Set the filename for soluteN.
setSolvent(n, filename) Set the filename for solventN.
setStream(stream, filename) Set the direction of the output stream to a filename
writeCommandFile(filename) Write the contents of this simulation object to disc
class protomslib.command.ProteinLigandSimulation(protein='protein.pdb', solutes=['solute.pdb'], solvent='water.pdb', templates=['solute.tem'], outfolder=None, ranseed=None)

Bases: protomslib.command.ProtoMSSimulation

This a generic command file for a protein-ligand simulation Generates input for proteins, solutes and solvent but does not write any chunks

Methods

clear() Clear the simulation description
readCommandFile(filename) Read in a simulation from a command file
run(exe[, cmdfile]) Run this simulation
setChunk(chunk) This adds the simulation chunk to ‘chunk’.
setDump(dump, freq) This adds the simulation dump to ‘dump’ with frequnecy freq
setForceField(filename) Set the filename for parfile
setParameter(parameter, value) Set the parameter to value
setProtein(n, filename) Set the filename for proteinN.
setSolute(n, filename) Set the filename for soluteN.
setSolvent(n, filename) Set the filename for solventN.
setStream(stream, filename) Set the direction of the output stream to a filename
writeCommandFile(filename) Write the contents of this simulation object to disc
class protomslib.command.ProtoMSSimulation(filename=None)

This is a ProtoMS command file object.This object holds all of the information about the simulation and can either run the simulation with this information, or can write this information to a command file.

Attributes:
lines :

list of strings all the lines of the ProtoMS command file

Methods

clear() Clear the simulation description
readCommandFile(filename) Read in a simulation from a command file
run(exe[, cmdfile]) Run this simulation
setChunk(chunk) This adds the simulation chunk to ‘chunk’.
setDump(dump, freq) This adds the simulation dump to ‘dump’ with frequnecy freq
setForceField(filename) Set the filename for parfile
setParameter(parameter, value) Set the parameter to value
setProtein(n, filename) Set the filename for proteinN.
setSolute(n, filename) Set the filename for soluteN.
setSolvent(n, filename) Set the filename for solventN.
setStream(stream, filename) Set the direction of the output stream to a filename
writeCommandFile(filename) Write the contents of this simulation object to disc
clear()

Clear the simulation description

readCommandFile(filename)

Read in a simulation from a command file

Parameters:
filename : string

the filename of the command file to read from disc

run(exe, cmdfile=None)

Run this simulation

Parameters:
exe : string

the filename of the executable

cmdfile : string, optional

the filename of the command file

setChunk(chunk)

This adds the simulation chunk to ‘chunk’.

Parameters:
chunk : string

the chunk to be executed, including parameters

setDump(dump, freq)

This adds the simulation dump to ‘dump’ with frequnecy freq

Parameters:
dump : string

the dump to be executed, including parameters

freq : int

the dump frequency

setForceField(filename)

Set the filename for parfile

Parameters:
filename : string

the filename of the force field file

setParameter(parameter, value)

Set the parameter to value

Parameters:
parameter : string

the name of the parameter

value : string

the value of the parameter

setProtein(n, filename)

Set the filename for proteinN. It then returns n+1

Parameters:
n : int

the protein serial number

filename : string

the filename of the protein pdb file

Returns:
int

n + 1

setSolute(n, filename)

Set the filename for soluteN. It then returns n+1

Parameters:
n : int

the solute serial number

filename : string

the filename of the solute pdb file

Returns:
int

n + 1

setSolvent(n, filename)

Set the filename for solventN. It then returns n+1

Parameters:
n : int

the solvent serial number

filename : string

the filename of the solvent pdb file

Returns:
int

n + 1

setStream(stream, filename)

Set the direction of the output stream to a filename

Parameters:
stream : string

the name of the stream

filename : string

the filename to direct the stream to

writeCommandFile(filename)

Write the contents of this simulation object to disc

Parameters:
filename : string

the filename to write the commands to

class protomslib.command.RestraintRelease(protein='protein.pdb', solutes=['solute1.pdb'], solvent='water.pdb', templates=['solute.tem'], nequil=5000000.0, nprod=40000000.0, dumpfreq=100000.0, lambdaval=None, ranseed=None, outfolder='out', restrained=[])

Bases: protomslib.command.ProteinLigandSimulation

Command file to calculate the free energy of introducing a harmonic restraint during the bound leg of a simulation.

Methods

clear() Clear the simulation description
readCommandFile(filename) Read in a simulation from a command file
run(exe[, cmdfile]) Run this simulation
setChunk(chunk) This adds the simulation chunk to ‘chunk’.
setDump(dump, freq) This adds the simulation dump to ‘dump’ with frequnecy freq
setForceField(filename) Set the filename for parfile
setParameter(parameter, value) Set the parameter to value
setProtein(n, filename) Set the filename for proteinN.
setSolute(n, filename) Set the filename for soluteN.
setSolvent(n, filename) Set the filename for solventN.
setStream(stream, filename) Set the direction of the output stream to a filename
writeCommandFile(filename) Write the contents of this simulation object to disc
class protomslib.command.Sampling(protein='protein.pdb', solutes=['solute.pdb'], solvent='water.pdb', templates=['solute.tem'], outfolder='', ranseed=None, nequil=5000000.0, nprod=40000000.0, dumpfreq=100000.0, outprefix='', tune=False)

Bases: protomslib.command.ProteinLigandSimulation

This a command file for MC sampling of a protein-ligand system Generates input for proteins, solutes and solvent and write chunks for equilibration and simulation, as well as dumps

Methods

clear() Clear the simulation description
readCommandFile(filename) Read in a simulation from a command file
run(exe[, cmdfile]) Run this simulation
setChunk(chunk) This adds the simulation chunk to ‘chunk’.
setDump(dump, freq) This adds the simulation dump to ‘dump’ with frequnecy freq
setForceField(filename) Set the filename for parfile
setParameter(parameter, value) Set the parameter to value
setProtein(n, filename) Set the filename for proteinN.
setSolute(n, filename) Set the filename for soluteN.
setSolvent(n, filename) Set the filename for solventN.
setStream(stream, filename) Set the direction of the output stream to a filename
writeCommandFile(filename) Write the contents of this simulation object to disc
class protomslib.command.SingleTopology(protein='protein.pdb', solutes=['solute1.pdb'], solvent='water.pdb', templates=['solute.tem'], nequil=5000000.0, nprod=40000000.0, dumpfreq=100000.0, lambdaval=None, ranseed=None, outfolder='out', restrained=[])

Bases: protomslib.command.ProteinLigandSimulation

This a command file for a single-topology protein-ligand simulation Generates input for proteins, solutes and solvent write lambda-replica exchange, dump results, pdb and restart files equilibrate and production run

Methods

clear() Clear the simulation description
readCommandFile(filename) Read in a simulation from a command file
run(exe[, cmdfile]) Run this simulation
setChunk(chunk) This adds the simulation chunk to ‘chunk’.
setDump(dump, freq) This adds the simulation dump to ‘dump’ with frequnecy freq
setForceField(filename) Set the filename for parfile
setParameter(parameter, value) Set the parameter to value
setProtein(n, filename) Set the filename for proteinN.
setSolute(n, filename) Set the filename for soluteN.
setSolvent(n, filename) Set the filename for solventN.
setStream(stream, filename) Set the direction of the output stream to a filename
writeCommandFile(filename) Write the contents of this simulation object to disc
protomslib.command.generate_input(protein, ligands, templates, protein_water, ligand_water, ranseed, settings)

Generates ProtoMS command files If protein is a valid filename a protein(-ligand) simulation is setup, and if ligands contains at least one valid filename a ligand simulation is setup.

Parameters:
protein : string

the filename of a protein pdb file

ligands : list of strings

filenames of solute pdb files

templates : list of strings

filenames of template files to be included

protein_water : string

the filename of a solvent pdb file for the protein

ligand_water : string

the filename of a solvent pdb file for the ligands

settings : Namespace object (from argparse)

additional settings

Returns:
free_cmd : ProtoMSSimulation

the command file for the ligand simulation

bnd_cmd : ProtoMSSimulation

the command file for the protein simulation

gas_cmd : ProtoMSSimulation

the command file for the gas simulation

FreeEnergy

class protomslib.free_energy.free_energy_argument_parser.FEArgumentParser(*args, **kwargs)

Bases: argparse.ArgumentParser

Thin wrapper around argparse.ArgumentParser designed to allow specification of individual arguments that cannot be used at the same time as one another.

Methods

add_argument(*args, **kwargs) Thin wrapper around ArgumentParser.add_argument that uses additional keyword argument ‘clashes’.
add_subparsers(**kwargs)
check_clashes(parsed) Check that arguments in Namespace parsed are compatible, according to provided argument clashes.
error(message) Prints a usage message incorporating the message to stderr and exits.
exit([status, message])
format_usage()
parse_args(*args, **kwargs) Thin wrapper around ArgumentParser.parse_args that checks for clashing arguments.
print_usage([file])
register(registry_name, value, object)
set_defaults(**kwargs)
add_argument_group  
add_mutually_exclusive_group  
convert_arg_line_to_args  
format_help  
format_version  
get_default  
parse_known_args  
print_help  
print_version  
add_argument(*args, **kwargs)

Thin wrapper around ArgumentParser.add_argument that uses additional keyword argument ‘clashes’. Clashes should be a list that contains the names of other arguments which cannot be provided at the same time as this one.

check_clashes(parsed)

Check that arguments in Namespace parsed are compatible, according to provided argument clashes.

parse_args(*args, **kwargs)

Thin wrapper around ArgumentParser.parse_args that checks for clashing arguments.

Collection of classes to form the basis of a replacement for the current free energy calculation framework.

class protomslib.free_energy.free_energy_base.BAR(lambdas, results_name='results', subdir_glob='./', **kwargs)

Bases: protomslib.free_energy.free_energy_base.Estimator

Estimate free energy differences using Bennett’s Acceptance Ratio.

Methods

add_data(series) Save data from a SnapshotResults.series object for later calculation.
calculate([temp]) Calculate the free energy difference and return a PMF object.
result_class alias of Result
subset([low_bound, high_bound, step]) Return a new class instance containing truncated data series.
add_data(series)

Save data from a SnapshotResults.series object for later calculation. Return the length of the data series added.

Parameters:
series: SnapshotResults.series

free energy simulation data for a particular lambda value

calculate(temp=300.0)

Calculate the free energy difference and return a PMF object.

Parameters:
temp: float, optional

temperature of calculation

class protomslib.free_energy.free_energy_base.BaseResult(*args, **kwargs)

Bases: object

A base class for Result objects.

class protomslib.free_energy.free_energy_base.Estimator(lambdas, results_name='results', subdir_glob='./', **kwargs)

Bases: object

Base class for free energy estimators.

Methods

add_data(series) Save data from a SnapshotResults.series object for later calculation.
calculate([temp]) Calculate the free energy difference and return a PMF object.
result_class alias of Result
subset([low_bound, high_bound, step]) Return a new class instance containing truncated data series.
add_data(series)

Save data from a SnapshotResults.series object for later calculation. Return the length of the data series added.

Parameters:
series: SnapshotResults.series

free energy simulation data for a particular lambda value

calculate(temp=300.0)

Calculate the free energy difference and return a PMF object.

Parameters:
temp: float, optional

temperature of calculation

result_class

alias of Result

subset(low_bound=0.0, high_bound=1.0, step=1)

Return a new class instance containing truncated data series.

Parameters:
low_bound: float, default=0.

Starting position of truncated data series.

high_bound: float, default=1.

Finishing position of truncated data series.

step: int, default=1

Subsampling frequency of data. 1=all data points, 2=every other data point, etc.

class protomslib.free_energy.free_energy_base.FreeEnergyCalculation(root_paths, temperature, estimators=[<class 'protomslib.free_energy.free_energy_base.TI'>, <class 'protomslib.free_energy.free_energy_base.BAR'>, <class 'protomslib.free_energy.free_energy_base.MBAR'>], subdir='', extract_data=True, **kwargs)

Bases: object

Class for performing a free energy calculation from one or more ProtoMS simulation outputs.

Methods

calculate([subset]) For each estimator return the evaluated potential of mean force.
run(args) Public method for execution of calculation.
calculate(subset=(0.0, 1.0, 1))

For each estimator return the evaluated potential of mean force.

Parameters:
subset: tuple of two floats and an int

specify the subset of data to use in the calculation

run(args)

Public method for execution of calculation. This is usually the desired way to use execute calculations. Handles logic for printing of output tables and saving of output pickles and figures. Calls plt.show() to display figures. Runs self._body(args) and uses the return value as the calculation result.

Parameters:
args
class protomslib.free_energy.free_energy_base.MBAR(lambdas, results_name='results', subdir_glob='./', **kwargs)

Bases: protomslib.free_energy.free_energy_base.TI

Estimate free energy differences using the Multistate Bennett’s Acceptante Ratio.

Methods

add_data(series) Save data from a SnapshotResults.series object for later calculation.
calculate([temp]) Calculate the free energy difference and return a PMF object.
result_class alias of Result
subset([low_bound, high_bound, step]) Return a new class instance containing truncated data series.
add_data(series)

Save data from a SnapshotResults.series object for later calculation. Return the length of the data series added.

Parameters:
series: SnapshotResults.series

free energy simulation data for a particular lambda value

calculate(temp=300.0)

Calculate the free energy difference and return a PMF object.

Parameters:
temp: float, optional

temperature of calculation

class protomslib.free_energy.free_energy_base.PMF(coordinate, *args)

Bases: protomslib.free_energy.free_energy_base.Series

A Potential of Mean Force, describing a free energy profile as a function of some coordinate.

Attributes:
dG

The free energy difference at the PMF end points

Methods

as_floats() Return a numpy array containing the values of the data series as floats without error values.
plot(axes[, show_error, fmt, xlabel, ylabel]) Plot this PMF onto the provided figure axes.
dG

The free energy difference at the PMF end points

class protomslib.free_energy.free_energy_base.Quantity(value, error)

Bases: object

Stores both the value and error associated with a free energy estimate.

Methods

fromData(data) Alternative initialiser that uses of a set of data points rather than providing an explicit value and error.
static fromData(data)

Alternative initialiser that uses of a set of data points rather than providing an explicit value and error. Returns a Quantity object with a value and error determined from the mean of and standard error of the provided data. If Quantity objects are provided as data their individual error attributes are ignored

Parameters:
data: sequence of numbers or Quantity objects

data from which to determine value and error

class protomslib.free_energy.free_energy_base.Result(*args, **kwargs)

Bases: protomslib.free_energy.free_energy_base.BaseResult

The result of a free energy calculation. Contains multiple PMF objects that are combined together to give a meaningful free energy difference.

Attributes:
dG

The free energy difference at the PMF end points

lambdas

Lambda values at which the calculation has been performed

pmf

The potential of mean force for this result

dG

The free energy difference at the PMF end points

lambdas

Lambda values at which the calculation has been performed

pmf

The potential of mean force for this result

class protomslib.free_energy.free_energy_base.Series(coordinate, *args)

Bases: object

Object for storing data series, i.e. some value that changes as a function of some coordinate. Designed to act as a base class for PMF classes.

Methods

as_floats() Return a numpy array containing the values of the data series as floats without error values.
plot(axes[, show_error, fmt, xlabel, ylabel]) Plot this PMF onto the provided figure axes.
as_floats()

Return a numpy array containing the values of the data series as floats without error values.

plot(axes, show_error=True, fmt='-', xlabel='Lambda Value', ylabel='Free Energy (kcal/mol)', **kwargs)

Plot this PMF onto the provided figure axes.

Parameters:
axes : matplotlib Axes object

axes onto which the figure will be drawn

show_error : boolean

if True, plot error bars

xlabel : string

label for x-axis

ylabel : string

label for y-axis

**kwargs :

additional keyword arguments to be passed to axes.plot

class protomslib.free_energy.free_energy_base.TI(lambdas, results_name='results', subdir_glob='./', **kwargs)

Bases: protomslib.free_energy.free_energy_base.Estimator

Estimate free energy differences using Thermodynamic Integration.

Methods

add_data(series) Save data from a SnapshotResults.series object for later calculation.
calculate([temp]) Calculate the free energy difference and return a PMF object.
result_class alias of Result
subset([low_bound, high_bound, step]) Return a new class instance containing truncated data series.
add_data(series)

Save data from a SnapshotResults.series object for later calculation. Return the length of the data series added.

Parameters:
series: SnapshotResults.series

free energy simulation data for a particular lambda value

calculate(temp=300.0)

Calculate the free energy difference and return a PMF object.

Parameters:
temp: float, optional

temperature of calculation

protomslib.free_energy.free_energy_base.get_alchemical_arg_parser()

Returns an argument parser suitable for use with alchemical simulation methods e.g. free energy perturbation.

protomslib.free_energy.free_energy_base.get_base_arg_parser()

Returns a generic argparser for all free energy calculation scripts

class protomslib.free_energy.gcmc_free_energy_base.GCI(B_values, volume, results_name='results', nsteps=None, nmin=None, nmax=None, nfits=10, **kwargs)

Bases: protomslib.free_energy.free_energy_base.Estimator

Estimate insertion free energies using Grand Canonical Integration.

Attributes:
N_max
N_min

Methods

add_data(series) Save data from a SnapshotResults.series object for later calculation.
calculate(temp) Calculate the free energy difference and return a GCMCPMF object.
result_class alias of protomslib.free_energy.free_energy_base.Result
subset([low_bound, high_bound, step]) Return a new class instance containing truncated data series.
N_max
N_min
add_data(series)

Save data from a SnapshotResults.series object for later calculation. Return the length of the data series added.

Parameters:
series: SnapshotResults.series

free energy simulation data for a particular B value

calculate(temp)

Calculate the free energy difference and return a GCMCPMF object.

Parameters:
temp: float

Temperature of calculation

class protomslib.free_energy.gcmc_free_energy_base.GCMCMBAR(lambdas, volume, **kwargs)

Bases: protomslib.free_energy.free_energy_base.MBAR

Estimate free energy differences using the Multistate Bennett Acceptance ratio for alchemical simulations performed in the Grand Canonical ensemble at multiple B values.

Attributes:
Bs

Return the B value used for this calculation.

Methods

B_to_chemical_potential(B, temperature) Convert a B value to the equivalent chemical potential for this simulation setup.
add_data(series) Save data from a SnapshotResults.series object for later calculation.
calculate([temp]) Calculate the free energy difference and return a PMF object.
equilibrium_B(temperature) Returns the value of B equating to equilibrium with bulk water.
result_class alias of protomslib.free_energy.free_energy_base.Result
subset([low_bound, high_bound, step]) Return a new class instance containing truncated data series.
B_to_chemical_potential(B, temperature)

Convert a B value to the equivalent chemical potential for this simulation setup. Depends on GCMC volume of the simulation.

Parameters:
B: float

the input B value

temperature: float

the simulated temperature

Bs

Return the B value used for this calculation.

add_data(series)

Save data from a SnapshotResults.series object for later calculation. Return the length of the data series added.

Parameters:
series: SnapshotResults.series

free energy simulation data for a particular simulation

calculate(temp=300.0)

Calculate the free energy difference and return a PMF object.

Parameters:
temp: float, optional

temperature of calculation, default 300 K

equilibrium_B(temperature)

Returns the value of B equating to equilibrium with bulk water.

class protomslib.free_energy.gcmc_free_energy_base.GCMCPMF(coordinate, values, volume, temperature, model, pmf)

Bases: protomslib.free_energy.free_energy_base.Series

A potential of mean force object for Grand Canonical calculations. Stores occupancy data and PMF for insertion free energies.

Attributes:
equilibrium_B

Returns the value of B equating to equilibrium with bulk water.

Methods

as_floats() Return a numpy array containing the values of the data series as floats without error values.
plot(axes[, show_error, fmt, xlabel, ylabel]) Plot this PMF onto the provided figure axes.
equilibrium_B

Returns the value of B equating to equilibrium with bulk water.

class protomslib.free_energy.gcmc_free_energy_base.GCMCResult(*args)

Bases: protomslib.free_energy.free_energy_base.BaseResult

Attributes:
B_values

Return the simulated B values

equilibrium_B

Return the B value representing equilibrium with bulk water

insertion_pmf

Return a PMF object containing insertion free energies

model

Return a Series object containing an average fitted model

occupancies

Return a Series object containing occupancy data

B_values

Return the simulated B values

equilibrium_B

Return the B value representing equilibrium with bulk water

insertion_pmf

Return a PMF object containing insertion free energies

model

Return a Series object containing an average fitted model

occupancies

Return a Series object containing occupancy data

class protomslib.free_energy.gcmc_free_energy_base.Slp(x=array(0), y=array(0), size=1)

Bases: object

Class for fitting and storing a monotonically increasing single layer perceptron with one input and one output.

Fitting is achieved by iteratively optimising one unit at a time, rather than the more typical back-propigation algorithm, to make use of SciPy’s monotonic optimisation routines.

Attributes:
size : integer

the number logsistic terms, equilivalent to the number of ‘hidden units’ in an artificial neural network

x : numpy array

the explanatory variable, which is the Adams value or chemical potential in GCMC

y : numpy array

the response variable, which is the average number of inserted grand canonical molecules

weights : numpy array

the matrix of coefficients of the logistic terms

offset : numpy array

the intercept of the artificial neural network. Grouped together with the weights when fitting

streams : list

contains the output of each logistic function or ‘unit’

predicted : numpy array

the predicted output of the artificial neural network

error: float

the error (either mean-squared, Huber, or absolute loss) of the artificial neural network with respect to the response variable.

Methods

absolute() Calculates the absolute loss of the model and updates the error
epoch([grad_tol, pin_min, pin_max, …]) Fits each unit and intercept of the artificial neural network once
fit_offset([grad_tol, cost, c]) Fits the intercept of the artificial neural network
fit_unit(unit[, grad_tol, pin_max, …]) Fits an individual unit of the artificial neural network
forward() Evaluates the output from ALL units for a given set of weights and updates the prediction of the whole network
forward_unit(unit) Evaluates the output from ONE unit for a given set of weights and updates the prediction of the whole network
huber([c]) Calculates the pseudo Huber loss function and updates the error of the model
msd() Calculates the mean-squared error of the model and updates the error
predict(x) Returns the prediction for the artificial neural network
randomise([pin_max]) Produces initial values of weights by random sampling
randsearch(samples) Performs a crude initial random search of the space of weights, and saves the best set of weights
train([iterations, grad_tol_low, …]) Iterates over a specified amount of training epochs.
absolute()

Calculates the absolute loss of the model and updates the error

epoch(grad_tol=0.001, pin_min=None, pin_max=None, monotonic=True, cost='msd', c=2.0)

Fits each unit and intercept of the artificial neural network once

Parameters:
grad_tol : float

the tolorance of the gradient that determines when optimisation stops

pin_min : float

the value of the fixed intercept. If None, it will be fitted.

pin_max : float

the value of the maximum value to model can produce. If None, it will be fitted

monotonic : boolean

whether to insure the artificial neural network is monotonically increasing

cost : string

the type of cost/loss function that will be used can be ‘msd’, ‘abs’ or ‘huber’

c : float

the paramater of the huber loss function

fit_offset(grad_tol=0.001, cost='msd', c=2.0)

Fits the intercept of the artificial neural network

Parameters:
grad_tol : float

the tolorance of the gradient that determines when optimisation stops

cost : string
the type of cost/loss function that will be used

can be ‘msd’, ‘abs’ or ‘huber’

c : float

the paramater of the Huber loss function

fit_unit(unit, grad_tol=0.001, pin_max=None, monotonic=True, cost='msd', c=2.0)

Fits an individual unit of the artificial neural network

Parameters:
unit : int

the index of the logistic function whose weights (parameters) will be optimised

grad_tol : float

the tolorance of the gradient that determines when optimisation stops

pin_max : float

the value of the maximum value to model can produce. If None, it will be fitted

monotonic : boolean

whether to insure the artificial neural network is monotonically increasing

cost : string

the type of cost/loss function that will be used can be ‘msd’, ‘abs’ or ‘huber’

c : float

the paramater of the Huber loss function

forward()

Evaluates the output from ALL units for a given set of weights and updates the prediction of the whole network

forward_unit(unit)

Evaluates the output from ONE unit for a given set of weights and updates the prediction of the whole network

Parameters:
unit : the index of the logisitic function
huber(c=2)

Calculates the pseudo Huber loss function and updates the error of the model

Parameters:
c : number

the pseudo Huber loss function parameter

msd()

Calculates the mean-squared error of the model and updates the error

predict(x)

Returns the prediction for the artificial neural network

Parameters:
x : the vector of Adams or chemical potentials
randomise(pin_max=True)

Produces initial values of weights by random sampling

Parameters:
pin_max : whether to constrain the random weights so that the sum

of the units cannot exceed the maximum value of the response variable

randsearch(samples)

Performs a crude initial random search of the space of weights, and saves the best set of weights

Parameters:
samples : how many initial guesses of the weights shall be attempted
train(iterations=100, grad_tol_low=-3, grad_tol_high=1, pin_min=None, pin_max=None, monotonic=True, cost='msd', c=2.0)

Iterates over a specified amount of training epochs. After each epoch, the fitting tolerance will decrease to gradually refine the fits.

Parameters:
iterations : int

the maximum number of training epochs that will be tested.

grad_tol_low : int

ten to the power of this value is the initial tolerance of the gradient

grad_tol_high : int

ten to the power of this value is the final tolerance of the gradient

pin_min : float

the value of the fixed intercept. If left blank, it will be fitted.

pin_max : float

the value of the maximum value to model can produce. If None, it will be fitted

monotonic : boolean

whether to insure the artificial neural network is monotonically increasing

cost : string

the type of cost/loss function that will be used can be ‘msd’, ‘abs’ or ‘huber’

c : float

the paramater of the huber loss function

class protomslib.free_energy.gcmc_free_energy_base.TitrationCalculation(root_paths, temperature, volume, nsteps=None, nmin=None, nmax=None, nfits=5, pin_min=None, **kwargs)

Bases: protomslib.free_energy.free_energy_base.FreeEnergyCalculation

Calculate free energy differences using the Grand Canonical Integration method. Collates water occupancy data from simulations carried out at diferent chemical potentials and gives NVT binding free energies.

Methods

calculate([subset]) For each estimator return the evaluated potential of mean force.
run(args) Public method for execution of calculation.
calculate(subset=(0.0, 1.0, 1))

For each estimator return the evaluated potential of mean force.

Parameters:
subset: tuple of two floats and an int

specify the subset of data to use in the calculation

protomslib.free_energy.gcmc_free_energy_base.fit_ensemble(x, y, size, repeats=20, randstarts=1000, iterations=100, grad_tol_low=-3, grad_tol_high=1, pin_min=False, pin_max=None, monotonic=True, cost='msd', c=2.0, verbose=True)

Fits a collection of ANN models to GCMC titration data. Each ANN is fitted with different set of initial ANN weights.

Parameters:
x : numpy array

the B/Adams values or chemical potentials used in a titration

y : numpy array

the average value of water at each chemical potential

repeats : int

the number of repeats of the fitting. This gives the “ensemble” of models

others

As specified in Slp.

protomslib.free_energy.gcmc_free_energy_base.get_gci_arg_parser()
protomslib.free_energy.gcmc_free_energy_base.insertion_pmf(N, gcmc_model, volume, T=298.15)

Calculates the free energy to insert water from ideal gas to the GCMC volume.

Parameters:
N : numpy array

vector of bound water numbers between which relative free energies will be calculated

gcmc_model : Slp object

an ANN model (single layer perceptron) that will be used to calculate the free energy

volume : float

the volume of the GCMC region in Angstroms^3

T : float

the temperature of the simulation in Kelvin

protomslib.free_energy.gcmc_free_energy_base.integrated_logistic(params, Bi, Bf)

The integral of the logistic function, for use in the function “insertion_pmf”. The problem is that it is prone to numerical instabilities when the exponents are large and positive. While one could fix the problem, I’ve opted for numerical integration instead in “insertion_pmf”. This function has been used to verify the numerical integration.

Parameters:
params : list or numpy array

the parameters of the logistic function

Bi : float

the lower value of the interval of integration

Bf : float

the higher value of the interval of integration

protomslib.free_energy.gcmc_free_energy_base.inverse_slp(model, y)

For a given ANN model, this function returns the value of x that corresponds to the value of y. It’s basically the inverse function of an ANN. Used to determine which B value produces a given average of waters.

Parameters:
model : Slp object

the ANN (single layer perceptron) to return the inverse function of

y : float

the value which you want to find the corresponding value of x

protomslib.free_energy.gcmc_free_energy_base.logistic(params, x)

Logistic function

Parameters:
params : list or numpy array

the three parameters of the logistic function

x : numpy array

the explanatory variable

protomslib.free_energy.gcmc_free_energy_base.plot_insertion_pmf(results, title='')

Convenience function to plot the insertion pmf data from repeats

protomslib.free_energy.gcmc_free_energy_base.plot_titration(results, ax, dot_fmt='b')

Convenience function to plot the titration data from repeats.

protomslib.free_energy.gcmc_free_energy_base.pseudohuber(r, c)

The pseudo Huber loss/cost function. Smooth version of Huber loss function.

Parameters:
r : float

residual, i.e. difference between predicted value and its target

c : float

parameter that determines roughly where residuals will have less weight

protomslib.free_energy.gcmc_free_energy_base.simplex_sample(size)

Samples uniformly over a simplex, such that the output sums to 1.

Parameters:
size : integer

the number of random numbers you want

Collection of classes to handle output table formatting

class protomslib.free_energy.table.Column(header='', value_fmt='%.4f', error_fmt='%.4f')

Bases: object

Table column entry.

Attributes:
width

Width of the column in number of characters

Methods

add_data(data) Add a entry to this column.
add_data(data)

Add a entry to this column.

Parameters:
data: Quantity object or number

entry to add

width

Width of the column in number of characters

class protomslib.free_energy.table.SubColumn(fmt='%.4f')

Bases: object

The left or right side of a Column.

Methods

add_data(data) Add a entry to this subcolumn.
add_data(data)

Add a entry to this subcolumn.

Parameters:
data: number

entry to add

class protomslib.free_energy.table.Table(title, fmts, headers=[], space=3)

Bases: object

Output table for displaying data

Methods

add_blank_row() Add a blank row to the table.
add_row(data) Add a row to the table.
add_blank_row()

Add a blank row to the table.

add_row(data)

Add a row to the table.

Parameters:
data: list of objects

table entries for this run, objects must be convertable to strings according to the corresponding column format

GCMC

protomslib.gcmc.clear_gcmcbox(gcmcbox, waters)

Removes solvent molecules from the GCMC/JAWS-1 box as they will not be able to move during the simulation

Parameters:
gcmcbox : string or PDBFile object

the gcmcbox

waters : string or PDBFile object

the water molecule to check

Returns:
int

the number of removed water molecules

PDBFile

the cleared solvation box

protomslib.gcmc.distribute_particles(box, particles, watermodel='t4p', resname='WA1', partnumb=None)

Randomly distribute molecules in a box

Parameters:
box : a list or dictionary

the origin and length of the box where the molecules will be placed

particles : string or PDB object

the number of waters to include in the box or the pdb object or filename with the molecules to include in the box

watermodel : string, optional

(only used when particles is a number) for instance, “t4p”, “tip4p”, “t3p” “tip3p”, the water model for the generated waters

partnumb : string,optional

(only used when particles is a file) the number of particles. If not specified it is set to be as many as there are in the file

Returns:
pdb object

the pdb object with the molecules distributed randomly in the box

protomslib.gcmc.make_gcmcbox(pdb, filename, padding=2.0)

Make a GCMC/JAWS-1 simulation box around a PDB-structure

Parameters:
pdb : PDBFile object

the PDB structure

filename : string

the name of the output file

padding : float, optional

the amount of extra space around the ligand to add

protomslib.gcmc.print_bequil(boxlen)

Prepare

protomslib.prepare.make_dummy(pdbfile)

Make a dummy PDB structure to match a molecule

It will place the dummy at the centre of the molecule

Parameters:
pdbfile : string or PDBFile

the pdb structure of the molecule

Returns:
PDBFile instance

the dummy structure

protomslib.prepare.pdb2pms(pdb_in, forcefield, conversion_file)

Convert atom names to ProtoMS standard

The protein object is not modified by this routine, but the Residue and Atom objects are.

Parameters:
pdb_in : PDBFile

the pdb file to modify inline

forcefield : string

the style of the input PDB-file

conversion_file :

a file with conversion rules

Returns:
a PDBFile instance with ProtoMS names
protomslib.prepare.read_convfile(file=None, inmode=None, outmode=None)

Reads a conversion file into a dictionary

Parameters:
file : string

the filename of the conversion file

inmode : string

the style of the input PDB-file

outmode : string

the style of the output PDB-file

Returns:
a dictionary of the conversion
Raises:
ValueError

if any of the parameters are none

protomslib.prepare.run_antechamber(lig, charge, resnam=None)

Wrapper for antechamber with default parameters and AM1-BCC charges

Parameters:
lig : PDBFile or string

the ligand to run Antechamber on

charge : int

the net charge of the ligand

resnam : string, optional

the residue name of the ligand

Returns:
string

the filename of the created prepi file

protomslib.prepare.run_parmchk(lig)

Wrapper for parmcheck with default parameters

Parameters:
lig : PDBFile or string

the ligand to run Parmcheck on

Returns:
string

the filename of the created frcmod file

protomslib.prepare.scoop.scoop(protein, ligand, innercut=16, outercut=20, flexin='full', flexout='sidechain', terminal='neutralize', excluded=[], added=[], scooplimit=10)

Generates a scoop from protein structure

The protein object is not modified by this routine.

Parameters:
protein : PDBFile

pdb instance for the protein to be scooped

ligand : PDBFile or string

Either pdb instance for ligand to be scooped around, or string giving the name of a file containing 3 floats to act as coords for scoop centre, or a string with three floating point numbers

innercut : float, optional

Maximum distance from ligand defining inner region of the scoop

outercut : float, optional

Maximum distance from ligand defining outer region of the scoop

flexin : string, optional

Gives the degree of flexibility for residues of the inner region Can be ‘rigid’, ‘sidechain’ or ‘flexible’

flexout : string, optional

As flexin but for residues of the outer scoop

terminal : string, optional

Determines what to do with terminal residues Can be ‘keep’, ‘doublekeep’,’neutralize’

excluded : list of int

List of indices for residues to be excluded from scoop

added : list of int

List of indices for residues to be included in outer scoop

Returns:
PDBFile

an object representing the scooped protein

protomslib.prepare.water.alignhydrogens(watcoords, template, tol=1e-06)
Rotates a template water molecule such that it has the same orientation as a reference water molecule. Useful for converting between different water models in cases when one wishes to retain the same orientations of hydrogens. It assumes that both the template and the reference water have maching oxygen locations.
Parameters:
watcoords : numpy array

the coordinates of the water molecule that will be aligned to. The first row must be the x,y,z coordinates of the oxygen atom.

Hydrogens must be present in the proceeding two rows

template : numpy array

the template water molecule whose hydrogens will be aligned to match that of watcoords. Does not have to be same water model (e.g. tip4p versus spc) as watcoords.

tol : float

the tolerance level for the BFGS solver. At a tolerance level of 0.1, orientating 1000 water molecules takes about 10s. With a lower tolerance, alignment is more accurate, but can take up to 17s per 1000 waters.

Returns:
numpy array

the rotated set of water coordinates.

protomslib.prepare.water.convertwater(pdb_in, watermodel, ignorH=False, watresname=None)

Converts water in a pdb object to ideal water geometries of an input model

The protein object is not modified by this routine, but the Residue and Atom objects are.

Parameters:
pdb_in : PDBFile

the PDB object containing the water molecules that will be converted

watermodel : string

the name of the water model that will be used in the transformtation, e.g. tip4p, t3p.

ignorH : bool

whether to ignore hydrogens in the input water.

Returns:
PDBFile

a pdb file whose solvent elements have the geometry of the desired solvent model

protomslib.prepare.water.rotatesolute(solutecoords, alpha, beta, gamma)
Rotates a molecule about the centre of geometry (the mean of x,y,z of each atom in the molecule).
Parameters:
solutecoords : numpy array

the coordinates of the molecule.

alpha : float or integer

the angle (in radians) by which rotation is performed about the x-axis.

beta : float or integer

the angle (in radians) by which rotation is performed about the y-axis.

gamma : float or integer

the angle (in radians) by which rotation is performed about the z-axis.

Returns:
numpy array

the rotated set of molecule coordinates.

protomslib.prepare.water.rotatewat(watcoords, alpha, beta, gamma)
Rotates a water molecule about the oxygen atom, which is assumed to be the first element in the coordinate array.
Parameters:
watcoords : numpy array

the coordinates of the water molecule. The first row must be the x,y,z coordinates of the oxygen atom. Hydrogens must be present in the proceeding two rows

alpha : float or integer

the angle (in radians) by which rotation is performed about the x-axis.

beta : float or integer

the angle (in radians) by which rotation is performed about the y-axis.

gamma : float or integer

the angle (in radians) by which rotation is performed about the z-axis.

Returns:
numpy array

the rotated set of water coordinates.

protomslib.prepare.water.set_jaws2_box(water_files, length=3)

Sets the header of the pdbs containing one single molecule to define a box around the first atom

water_files : PDBSet
a pdb set containing the pdb files with one molecule
length : int, optional
the dimensions of the box around the molecule
protomslib.prepare.water.split_waters(waters)

Split waters in a PDB file to many PDB-files for JAWS-2

Parameters:
waters : string or PDBFile object

the PDB structure

Returns:
PDBSet

single waters

PDBSet

excluding single waters

Simulation Objects

Useful classes for setup and analysis of ProtoMS simulations. Contain classes to

  1. Read, modify and write structures in PDB format
  2. Store parameter collections
  3. Read and store ProtoMS results files
  4. Read, modify and write ProtoMS template files
class protomslib.simulationobjects.Atom(index=0, name='?', resname='???', resindex=0, coords=[])

Bases: object

Class for holding a PDB atom record

Attributes:
index : integer

the number of the atom in the pdb file

resindex : integer

the number of the residue in the pdb file

name : string

the name of the atom

resname : string

the name of the residue

coords : NumpyArray

Cartesian coordinates

type : string

either atom or hetatm

element : string

the element of the atom

Methods

getElement() Set the element of this atom
getElement()

Set the element of this atom

class protomslib.simulationobjects.AtomSet(record=None)

Class to hold a set of atoms and their associated parameter

Attributes:
atoms : list of strings

the atoms in this set

param : int or ForceFieldParameter

the parameter associated with the atoms

Methods

parse(record) Parse a line from a ProtoMS file
parse(record)

Parse a line from a ProtoMS file

Parameters:
record : string

the line to parse from

class protomslib.simulationobjects.DihParameter(index, ats, terms)

Class to hold a dihedral parameter from a ProtoMS template file

A separate class is needed to handle dihedral parameters due to their composition of multiple terms.

Attributes:
index : integer

the serial number of the parameter

ats : list of string

the name of the atoms associated with the parameter

b0 : list of DihTerm

the terms that make up this dihedral parameter

class protomslib.simulationobjects.DihTerm(index, k1, k2, k3, k4)

Class to hold individual terms making up a Dihedral

Attributes:
index : integer

the serial number of the parameter

k1 : float

The term barrier height

k2 : float

OPLS inversion constant (should be either 0.0 or 1.0 for AMBER ffs)

k3 : float

The term periodicity

k4 : float

The term phase

class protomslib.simulationobjects.EnergyResults(line=None)

Class to hold energy results

Attributes:
curr : float

the energy at the current lambda

back : float

the energy at the previous lambda

forward : float

the energy at the next lambda

type : string

a label

Methods

parse_line(line) Parse energies from the result file
parse_line(line)

Parse energies from the result file

Parameters:
line : string

the line to parse

class protomslib.simulationobjects.ForceFieldParameter(record=None)

Class to hold a general parameter set

Attributes:
key : string

a label

index : integer

a serial number

params : list

either a list of strings or a list of ForceFieldParameter, containing the actual parameters

Methods

parse(record) Parse a line from a ProtoMS file
parse(record)

Parse a line from a ProtoMS file

Parameters:
record : string

the line to parse from

class protomslib.simulationobjects.MolTemplate

Class to hold a ProtoMS template

Attributes:
name : string

the name of the template

type : string

the kind of template, e.g. solute

translate : float

the translation displacement

rot : float

the rotational displacement

atoms : list of TemplateAtom objects

the atoms in this template

connectivity : list of TemplateConnectivity objects

the connectivity in this template

variables : list of string

the variable geometries

atomclass : TemplateAtom class

the class to create objects of atoms

Methods

parse_from(fileobj) Parse a full template from a file
write_to(fileobj) Write the template to disc
write_zmat(filename) Write the z-matrix of this template to disc
parse_from(fileobj)

Parse a full template from a file

Terminates when found another mode or end of file

Parameters:
fileobj : file object

the file to parse from

write_to(fileobj)

Write the template to disc

Parameters:
fileobj : file object

the file to write to

write_zmat(filename)

Write the z-matrix of this template to disc

Parameters:
filename : string

the name of the file to write to

class protomslib.simulationobjects.MyArgumentParser(prefix_chars='-', add_fullhelp=True, **kwargs)

Bases: argparse.ArgumentParser

Extention to ArgumentParser to allow separation of of help in to help for most common options and full help that list all arguments

Methods

add_argument(dest, …[, name, name])
add_subparsers(**kwargs)
error(message) Prints a usage message incorporating the message to stderr and exits.
exit([status, message])
format_usage()
parse_args([args, namespace])
print_usage([file])
register(registry_name, value, object)
set_defaults(**kwargs)
add_argument_group  
add_mutually_exclusive_group  
convert_arg_line_to_args  
format_help  
format_version  
get_default  
parse_known_args  
print_fullhelp  
print_help  
print_version  
format_help(fullhelp=False)
print_fullhelp(file=None)
print_help(file=None)
class protomslib.simulationobjects.PDBFile(filename=None)

Class for holding the atoms and residues of a PDB file

Attributes:
center : NumPy array

the center of coordinates of all atoms

header : string

additional lines to print before ATOM records

name : string

the name of the file

residues : dictionary of Residue objects

all non-solvent residues

solvents : dictionary of Residue objects

all solvent residues

Methods

copy() Make a copy of the residues and solvents dictionaries, not the Residue and Atom objects themselves
getBox([atomlist, reslist]) Calculate the smallest box to encompass the atoms
getCenter() Calculate the center of geometry
read(filename) Read a pdb-file from disc
read_from(f[, resname]) Read a pdb structure from a fileobject
write(filename, **kwargs) Write the PDB file to disc
copy()

Make a copy of the residues and solvents dictionaries, not the Residue and Atom objects themselves

Returns:
PDBFile

the newly created instance

getBox(atomlist='all', reslist='all')

Calculate the smallest box to encompass the atoms

Parameters:
atomlist = list, optional

the atoms to be taken into acount to get the box

reslist = list, optional

the residues to be taken into acount to get the box

Returns:
dictionary of Numpy arrays

the center, length and origin of the box

getCenter()

Calculate the center of geometry

Returns:
NumPy array

the center

read(filename)

Read a pdb-file from disc

It will overwrite existing residues and renumber all residues from 1

Parameters:
filename : string

the name of the pdb file

read_from(f, resname=None)

Read a pdb structure from a fileobject

Terminates reading when found a record starting with END

Parameters:
f : file object

the object to read from

resname : string, optional

the name of the residue to read, only read this

write(filename, **kwargs)

Write the PDB file to disc

Parameters:
filename : string

the name of the file

renumber : boolean, optional

indicate if residues should be renumbered

header : string, optional

additional lines to print before the atom records

solvents : boolean, optional

if to write the solvents

class protomslib.simulationobjects.PDBSet

Hold a collection of PDBFile objects

Attributes:
pdbs : list of PDBFile objects

the pdb files

Methods

from_residues(pdbfile) Split a PDB file into a set by residue
read(filename[, resname, skip, readmax]) Read a set of pdb structures from a file
write(filenames[, solvents]) Write the set to disc
from_residues(pdbfile)

Split a PDB file into a set by residue

Parameters:
pdbfile : PDBFile object

the pdb-file to split

read(filename, resname=None, skip=0, readmax=None)

Read a set of pdb structures from a file

Parameters:
filename : string

the name of the file to read

resname : string, optional

the name of the residue to read, only read this

skip : int, optional

number of snapshot at the beginning to skip

readmax : int, optional

maximum number of snapshots to read

write(filenames, solvents=True)

Write the set to disc

If one filename is given all PDB-files are written to one file

Parameters:
filenames : list of strings

the name of the file(s) to write to

solvents : boolean

if to write the solvents

class protomslib.simulationobjects.Parameter(index, ats, k, b0)

Class to hold a parameter from a ProtoMS template file

Not valid for dihedral parameters but contains sufficient information to check that a parameter exists

Attributes:
index : integer

the serial number of the parameter

ats : list of string

the name of the atoms associated with the parameter

k : float

the force constant

b0 : float

the equilibrium value

class protomslib.simulationobjects.ParameterSet(ptype, param_file=None)

Class to hold a collection of parameters

Attributes:
file_name : string

the name of the file where the parameters originated

params : list of Parameter objects

the parameters in this collection

Methods

get_par_by_index(ind) For an atm line index find matching the par line
get_params(ats) Find parameter sets for specific atoms.
read  
formats = {'angle': {'par': (<type 'str'>, <type 'int'>, <type 'float'>, <type 'float'>), 'atm': (<type 'str'>, <type 'str'>, <type 'str'>, <type 'str'>, <type 'int'>)}, 'bond': {'par': (<type 'str'>, <type 'int'>, <type 'float'>, <type 'float'>), 'atm': (<type 'str'>, <type 'str'>, <type 'str'>, <type 'int'>)}, 'clj': {'par': (<type 'str'>, <type 'int'>, <type 'str'>, <type 'int'>, <type 'float'>, <type 'float'>, <type 'float'>)}, 'dihedral': {'par': (<type 'str'>, <type 'int'>, <type 'int'>, <type 'int'>, <type 'int'>, <type 'int'>, <type 'int'>), 'atm': (<type 'str'>, <type 'str'>, <type 'str'>, <type 'str'>, <type 'str'>, <type 'int'>), 'term': (<type 'str'>, <type 'int'>, <type 'float'>, <type 'float'>, <type 'float'>, <type 'float'>)}}
get_par_by_index(ind)

For an atm line index find matching the par line

get_params(ats)

Find parameter sets for specific atoms.

Parameters:
ats : sequence of strings

the query atoms

Returns:
Parameter object

the found parameter

read(param_file)
class protomslib.simulationobjects.Residue(name='???', index=0)

Bases: object

Class for holding a set of PDB atoms

Attributes:
atoms : list of Atom objects

the atom of this residue

index : integer

the number of the residue in the pdb file

name : string

the name of the residue

center : NumPy array

the center of coordinates of all atoms

Methods

addAtom(atom) Adds an atom to this residue
getCenter() Sets the center of coordinates for this residue
isAminoacid() Checks is the residue name is an aminoacid name
addAtom(atom)

Adds an atom to this residue

Parameters:
atom : Atom object

the atom to add

getCenter()

Sets the center of coordinates for this residue

isAminoacid()

Checks is the residue name is an aminoacid name

class protomslib.simulationobjects.RestartFile(filename=None)

Class to hold the information on a restart file

Attributes:
filename : string

the name of the restart file

nsolutes : int

the number of solute molecules

ngcsolutes : int

the number of gcsolute molecules

solutesdic : dict

the correlation of solute ids with residue names

gcsolutesdic : dict

the correlation of gcsolutes ids with residue names

Methods

read(filename) Read the restart file from disc
read(filename)

Read the restart file from disc

Parameters:
filename : string

the name of the file to read

class protomslib.simulationobjects.ResultsFile

Class to store a collection of results

Attributes:
filename : string

the name of the file

snapshots : list of SnapshotResult

all the results

series : SnapshotResult

all the results in a single object

Methods

make_series() Make a snapshot with all the results
read(filename[, skip, readmax]) Read results from disc
make_series()

Make a snapshot with all the results

This routine creates a special SnapshotResults object, which has the same attributes as the objects in the snapshots list, but rather than storing single floats and integers, this object stores NumpyArrays of all the results

Returns:
SnapshotResults

the created object, also stored in self.series

read(filename, skip=0, readmax=None)

Read results from disc

The filename parameter can be either a single filename, in case it is assumed that it is a ProtoMS3.0 results file, i.e. it contains several snapshots if it is a list of strings, two behaviors are supported:

  1. if the length of the list is 2, and the first filename

is a dictionary, a glob is used to list all filenames in that directory that starts with the second filename in the list

  1. if the length is not two or you give two filenames,

each of them are processed individually and treated as individual snapshots

Parameters:
filename : string or list of strings

the name of the file to read

skip : int, optional

number of snapshot at the beginning to skip

readmax : int, optional

maximum number of snapshots to read

exception protomslib.simulationobjects.SetupError

Bases: exceptions.Exception

A general exception class to be raised by setup code

class protomslib.simulationobjects.SnapshotResults

Class to store a single snapshot of results

Not all attributes might not be set as they might not exists in the result file

internal_energies and interaction_energies are dictionary where the key is is the molecules involved in the energy, e.g. protein1, or protein1-solute1. The value is a list of EnergyResults objects, for the various energy components.

Attributes:
lam : float

the current lambda value

lamb : float

the previous lambda value

lamf : float

the next lambda value

datastep : int

the number of steps the averages are calculate over

lambdareplica : int

the index of the lambda replica

temperaturereplica : int

the index of the temperature replica

global replica : int

the index of the global replica

temperature : float

the temperature of the simulation

efftemperature : float

the effective temperature in case of REST

ngcsolutes : int

the number of GC solutes

nthetasolutes : int

the number of theta solutes

bvalue : float

the Adams value

solventson : int

the number of GC solutes that are turned on

pressure : float

the pressure

volume : float

the volume

seed : int

the random seed

backfe : float

the free energy to the previous lambda

forwfe : float

the free energy to the next lambda

total : EnergyResults object

the total energy

internal_energies : dictionary of lists of EnergyResults object

the internal energies

interaction_energies : dictionary of lists of EnergyResults objects

the interaction energies

capenergy : EnergyResults object

the energy of the cap

extraenergy : EnergyResults object

all extra energies

feenergies : dictionary of float

the total energy at various lambda values

gradient : float

the numerical gradient of lambda

agradient : float

the analytical gradient of lambda

thetavals : list of float

the theta values of all GC solutes

thetasolvals : list of float

the theta values of all theta solutes

Methods

parse(fileobj) Parse a file for results
parse(fileobj)

Parse a file for results

Parameters:
fileobj : file object

the file to read from

Returns:
string :

the line last read

class protomslib.simulationobjects.TemplateAtom(record)

Class to hold a template atom

Attributes:
name : string

the name of the atom

residue : string

the residue of the atom

param0 : int or ForceFieldParameter

the parameter associated with this atom at v0

param1 : int or ForceFieldParameter

the parameter associated with this atom at v1

bondedto : string or TemplateAtom

the atom this atom is bonded to

angleto : string or TemplateAtom

the atom this atom forms an angle with

dihedto : string or TemplateAtom

the atom this atom forms a dihedral with

Methods

parse(record) Parse a line from a ProtoMS file
zmat() Produces a simple z-matrix representation of this atom
parse(record)

Parse a line from a ProtoMS file

Parameters:
record : string

the line to parse from

zmat()

Produces a simple z-matrix representation of this atom

Returns:
returns

the z-matrix representation

class protomslib.simulationobjects.TemplateConnectivity(record=None)

Class to hold a solute bond, angle or dihedral

Attributes:
type : string

the kind of connectivity, i.e. bond, angle or dihedral

atoms : list of strings or TemplateAtom objects

the atoms involved in this connectivity

flex : float

the flexibility

param0 : int or ForceFieldParameter

the parameter associated with this connectivity at v0

param1 : int or ForceFieldParameter

the parameter associated with this connectivity at v1

dummy : bool

flag to indicate if this should be treated as a dummy

Methods

parse(record) Parse a line from a ProtoMS file
parse(record)

Parse a line from a ProtoMS file

Parameters:
record : string

the line to parse from

class protomslib.simulationobjects.TemplateFile(filename=None)

Class to hold a ProtoMS template file

Attributes:
bondparams : list of ForceFieldParameter objects

all bond parameters

bondatoms : list of AtomSet objects

all atoms associated with bond parameters

angleparams : list of ForceFieldParameter objects

all angle parameters

angleatoms : list of AtomSet objects

all atoms associated with angle parameters

dihedralterms : list of ForceFieldParameter objects

all dihedral term parameters

dihedralparams : list of ForceFieldParameter objects

all dihedral parameters

dihedralatoms : list of AtomSet objects

all atoms associated with dihedral parameters

cljparams : list of ForceFieldParameter objects

all clj parameters

templates : list of MolTemplate objects

all templates in this file

Methods

append(other) Adds another template file to this one
assign_paramobj() Replaces integers and strings with objects
read(filename) Read a template file from disc
write(filename) Write a template file to disc
append(other)

Adds another template file to this one

Will renumber the parameters of the other file

Parameters:
other : TemplateFile

the file to add

assign_paramobj()

Replaces integers and strings with objects

It replaces all indices to force field parameters with references to ForceFieldParameter objects

It replaces all atom in templates with references to TemplateAtom objects

read(filename)

Read a template file from disc

Parameters:
filename : string

the name of the file to read

write(filename)

Write a template file to disc

Parameters:
filename : string

the name of the file to write to

class protomslib.simulationobjects.TemplateSoluteAtom(record=None)

Bases: protomslib.simulationobjects.TemplateAtom

Class to hold a solute atom

This is a sub-class that implements functions that are specific for solute templates

Has the same attributes as TemplateAtom

Methods

parse(record) Parse a line from a ProtoMS file
zmat() Produces a simple z-matrix representation of this atom
parse(record)

Parse a line from a ProtoMS file

Parameters:
record : string

the line to parse from

zmat()

Produces a simple z-matrix representation of this atom

Returns:
returns

the z-matrix representation, viz. ATOM BONDEDTO ANGLETO DIHEDRALTO

protomslib.simulationobjects.angle(v1, v2)

Calculates the angle between two vectors

Parameters:
v1 : Numpy array

the first vecor

v2 : Numpy array

the second vecor

Returns:
float

the angle in radians

protomslib.simulationobjects.angle_atms(a1, a2, a3)

Calculates the angle between three atoms

Parameters:
a1 : Numpy array

the first atom

a2 : Numpy array

the second atom

a3 : Numpy array

the third atom

Returns:
float

the angle in radians

protomslib.simulationobjects.color(idx)

Returns a color of index

For instances when the index is larger than the number of defined colors, this routine takes care of this by periodicity, i.e. color at idx=0 is the same color as idx=n

Parameters:
idx : int

the index of the color

Returns:
list of floats

the color

protomslib.simulationobjects.find_box(pdbobj)

Parse box information from a pdb file

Looks for CENTER, DIMENSIONS and box keywords in the header of the pdb file

Parameters:
pdbobj : PDBFile object

the structure to parse

Returns:
dictionary of Numpy arrays

the center, length and origin of the box

protomslib.simulationobjects.is_solvent(name)

Check whether a given residue name is a solvent residue name

Parameters:
name : string

the residue name

protomslib.simulationobjects.merge_pdbs(pdbobjs)

Merge pdb files

Create a new PDBFile instance and add all residues and solvents from all pdb files given, renumbering residues

Parameters:
pdbobjs : list of PDBFile objects
Returns:
PDBFile object

the newly created and merged pdb structure

protomslib.simulationobjects.setup_logger(filename=None)

Setup ProtoMS logging system

Uses the standard logging module with two handlers, one that prints everything above DEBUG to standard output and one that prints everything, even DEBUG to a file

If filename is None no file handler is created

This should be called by protoms.py and all other stand-alone programs that uses the tools

Parameters:
filename : string, optional

the filename of the string

Returns:
a reference to the created logger
protomslib.simulationobjects.standard_filename(filename, folder)

Generates the filename of file in the standard ProtoMS file hierarchy

If $PROTOMSHOME is set, it uses it as the base-path, otherwise, it uses the location of this module as to create the base-path

Parameters:
filename : string

the name of the file

folder : string

a folder in the standard ProtoMS file hierarchy

Returns:
the full path to the file
protomslib.simulationobjects.write_box(filename, box)

Write a box in PDB file format

Parameters:
filename : string

the name of the file to write

box : dictionary of Numpy array

the box specification

Solvate

protomslib.solvate.solvate(box, ligand=None, protein=None, geometry='box', padding=10.0, radius=30.0, center='cent', namescheme='ProtoMS', offset=0.89)

Function to solvate ligand and/or protein structures using a pre-equilibrated box of waters.

Parameters:
box : String

String giving the name of a pdb file containing a box of pre- equilibrated water molecules for solvation

ligand : PDBFile or string, optional

Either pdb instance for ligand to be solvated, or string giving the name of a pdb file of the ligand. NOTE - either ligand or protein or both must be supplied

protein : PDBFile or string, optional

Either pdb instance for protein to be solvated, or string giving the name of a pdb file of the protein. NOTE - either ligand or protein or both must be supplied

geometry : string, optional

Shape of solvent box, options “box”, “droplet”, “flood”. Default box. Droplet defines a sphere of radius ‘radius’ around center ‘center’. Box defines a box that extends at least distance ‘padding’ from solute Flood defines the equivalent, but waters overlapping with the solute are not removed

padding : float, optional

Minimum distance between the solute and the box edge. Default 10.0A

radius : float, optional

The radius of the added droplet sphere. Default 30.0A

center : string, optional

Defines the center of the added droplet sphere. Options “cent”, one, two or three numbers. Default “cent”. Cent defines the center based on all solute atom coordinates One number defines coordinates for the center as identical x,y,z Two numbers defines coordinates for the center over an atom range Three numbers defines coordinates for the center as separate x,y,z

namescheme : string, optional

Naming scheme for printout, options “ProtoMS”, “Amber”. Default ProtoMS

Returns:
PDBFile

the created box of waters, including crystallographic water

Raises:
SetupError
neither protein nor ligand was given, or x-ray waters is not

same format as pre-equilibrated box

Template

Utils

protomslib.utils.rotmat_x(alpha)

3D rotation matrix for rotations about x-axis.

Parameters:
alpha : float or integer

the angle (in radians) by which rotation is performed about the x-axis.

Returns:
numpy array

the rotation matrix for the desired angle.

protomslib.utils.rotmat_y(beta)

3D rotation matrix for rotations about x-axis.

Parameters:
beta : float or integer

the angle (in radians) by which rotation is performed about the y-axis.

Returns:
numpy array

the rotation matrix for the desired angle.

protomslib.utils.rotmat_z(gamma)

3D rotation matrix for rotations about x-axis.

Parameters:
gamma : float or integer

the angle (in radians) by which rotation is performed about the z-axis.

Returns:
numpy array

the rotation matrix for the desired angle.