Augment model with python script

[ajkluber]

I want to add a feature that allows the user to modify the model object (e.g. by adding interactions) using an external python script. I am proposing a solution but I want to get feedback (mainly from @TensorDuck ) about how to implement it without disrupting recent changes to inputs.py.

I have found this to be necessary when creating the model of SNARE that @gsc4 has been using, which required a harmonic interaction between the endpoints to model the effects of the pulling potential. However, I think this solution will solve all these kinds of 'tinkering' problems we may have.

Specifically, I want to add an argument to the .ini file that defines the path to an external python script that the user has created (e.x. path_to_py = /home/ajk8/scratch/modify.py). The external script has a function that takes in a model object and returns a model object.

I basically just want to add the following code to inputs.py,

if not (modelopts["path_to_py"] is None):
    if not os.path.exists(path_to_py):
            raise IOError(path_to_py + " does not exist!")
    else:
        import imp
        modify_py = imp.load_source("DUMMY", path_to_py)
        model = modify_py.augment_model(model)

[TensorDuck]

In short, I think it's an okay modification ,but I have some reservations on the implementation that I think we should consider.

When we first set out to do the refactor, we were aiming to make the package more similar to the pyemma package. I think the philosophy of that package has been that the built in functions can do all the heavy lifting, but specific little tinkering things like adding a potential here, or there, is left to the user to add in some way. My understanding of this "some way" meant the user would have their python script, and inside the script they would load the model with the inputs.py, and then do any modifications to the model followed by using the writer objects to output the gromacs files.

So what I imagine would your current python script is something like this:

import model_builder

model = model_builder.inputs.load_model("bpti.ini")
model.add_umbrella_potential(minima, K, *args)

And essentially the augment_model object would just carry out the model.add_umbrella_potential function. A more meaningful example of what I mean would be something like:

import model_builder

model = model_builder.inputs.load_model("bpti.ini")
elastic_network_model_contacts = [ ... ]
for enm in elastic_network_model_contacts:
    model.add_harmonic_potential(*args)

Then this would make more sense to want to handle it automatically, as how you construct the enm might change a lot depending on your specific model (Ca vs CaCb, cutoff distances, etc.). We don't need to add a helper function in the model objects just to add the harmonic potential, but for example it could be, model.Hamiltonian.add_potential("harmonic", *args).

Some thoughts on this:

1. There needs to be a check in place on the model object. With your proposed example, the model object that comes out my not even be the same object that went in. Or, more likely, it might use a potential that doesn't exist, or use features that won't be interpretable by the writer objects later.
2. I think to stay in line with the philosophy of the package, the potentials you use should be implemented in the model_builder package if it's a reasonably common potential. Some wacky polynomial doesn't need to be implemented, but the harmonic potential should, and I think a generic pulling potential wouldn't be bad.
3. The inputs.py scripts is basically our automated way to make a model, so it is appropriate to include the option to supply the augment_model function. Especially in a case of an ENM, you might want to automatically determine which ENM based on some criteria based on the mapping, so you include the option to do so an call the model to add the necessary potentials.

[ajkluber]

I should provide a little more context that led me to this idea. I was using my inherent_structure package that loads in a model in the usual way in the source which doesn't allow me to add these extra interactions before the calculation runs (as I normally would in the manner you suggested if I was running simulation). In a similar way, this would be a problem for any external program that wants to create a model from just an .ini file.

In response to your thoughts:

1. Yeah, I think the fact that model = augment_model(model) could break later functionality in an unforeseeable way makes this especially tricky. There are some checks for supported functions in the output writers, but this could yield some bugs that are hard to track down.
2. Yeah, I agree that harmonic potentials are probably a common case but I am thinking along the lines of allowing for a user-defined pariwise potential. For example, I see a great benefit from having a custom pairwise potential class so the user can define functions without having to explicitly take the derivative. This would allow us to play around with more functional forms. Greg and I have been thinking of an interpolation formula for our desolvation potential which uses an unknown number of terms.
3. Yea, agreed. I think the best solution would allow an .ini file to be a self-contained way to create a model. I liked the idea of listing a script as an argument in the .ini so that the model could be changed in different contexts if necessary. These modifications could either be "extra" pieces to the Hamiltonian, like an umbrella potential, or just some custom interaction as I mention above.

I will wait to implement this for the time being.