Implementation of the Martini3 Go-model

[fgrunewald]

work in progress but almost finished;

To Do

• Write the atomtypes file for the Go Potential
• Check why MAD and Poma server + Martinize2 yield different Go models
• Check IDRs are correctly biased
• Try implement Chris's dihedral fix for pure IDPs
• align ouput file names with CreateGoVirt script
• fix the CLI arguments to pass all necessary information for the IDP/IDR functionality
• write tests
• write integration tests
• make gmx topfile writer a thing
• have the domain argument like for En

@pckroon I've added the missing Go functionality so far that it can now read a contact map, but not do the contact analysis. Because the functionality missing from the Go code concerns interactions written to the atomtypes/nonbondparams directives it is not a processor in the original sense. Hence I've dumped the files in an extra toplevel directory. Does that make sense or should we rather have it as a processor?

@Lp0lp Can you generate a reference contact map using the CreateGoVirt script for lysozyme? The canonical AA structure can be found here. Unfortunately, it doesn't run on Mac (maybe just my Mac). Do you know if the format of the contact map in terms of columns is always the same?

Also this project has a hard deadline October 10.

[pckroon]

First of all, great! I'm very happy to have this functionality in vermouth/martinize2.

I do think it should be a Processor (which can live in its own folder, see dssp), the question is what to do with the non-bonded parameter info. I'm tempted to overengineer this, but YAGNI. So maybe just dump them in a format that's convenient for you in either System or Forcefield. I think my preference would be in System (analogous to Molecule.interactions). Add a comment to that to indicate we'll engineer something proper when we need to touch non-bonded parameters and atom types for a second time.

align ouput file names with CreateGoVirt script

This is super optional to me. I'm more than happy to make backward breaking changes wrt that script.

[Lp0lp]

Heyo! I've generated the contact map for lysozyme using the old, new, and locally executable contact map generator. Annoyingly, the format is not mega-consistent across the 3 versions. They have pretty different headers (not a big deal, I think) but more annoying is their substantially different column layouts. I think the create_goVirt script can only use the map from the old webserver currently.

On top of that... besides the different file layout the output is also quite different..... but that is a whole other deal...

I've attached below the maps for you to have a look at:
exec.txt
webservernew.txt
webserverold.txt

EDIT: Also note that this is probably why MAD gives different go-models. Since mad uses the locally executable version of the contact map generator, which besides having a different layout seems to also spit out a different map altogether.

[Lp0lp]

align ouput file names with CreateGoVirt script

This is super optional to me. I'm more than happy to make backward breaking changes wrt that script.

The file names out of creategovirt are a mess, if we can clean them up now then all the better I think.

[Lp0lp]

They have pretty different headers (not a big deal, I think) but more annoying is their substantially different column layouts. I think the create_goVirt script can only use the map from the old webserver currently.

I just revised the contact maps again, the executable and old webservers are actually consistent. I was wrongly looking at the first table when it is the second one (line 5454 and onwards) that the script reads. This residue-residue contact table is in the same format as the executable it just has more stuff printed out before it. IMO, we drop support for the output out of the new webserver for now, and I can send them an email later asking them to make it more consistent perhaps.

[fgrunewald]

@pckroon why do we want the atomtypes to be dumped into system or molecules? They can be written to file directly because they are coupled to molecule i.e. two molecules cannot have the same set of atomtypes

[pckroon]

Because it fits better into the general workflow imho. It also makes the include statements easier in the top file. I do understand your point though.

[fgrunewald]

@Lp0lp Last but not least I'm refactoring the water biasing a little, such that it can be applied without the Go model. I think that's also the original idea of the IDPs. One question I have is whether you think someone would only want to bias IDRs in a structured protein and NOT the rest of the protein (i.e. applying idp fix only). My idea was a CLI like so:

martinize2 -water-bias -> applies the auto water bias
martinize2 -water-bias-eps -> you can set the epsilon for all sec-struc types including an "idr" key for IDPs/IDRs
martinize2 -id-regions -> define which regions are IDRs

In this scenario, any IDR biasing overwrites the default auto biasing. If we have default setting for the sec-struc biases the only tricky part would be to only apply biases to the IDP regions because that would require setting the specific biases to 0.

Any thoughts?

[pckroon]

all sec-struc types including an "idr" key for IDPs/IDRs

What do you mean by this? Would it make sense to allow users to set which sec-struct types to apply the bias to as well (with a sane default)?

[fgrunewald]

like as follows:

-water-bias-eps H:2.1 C:2.3 idr:2.1

This means applying a bias of 2.1 kJ/mol to helix but 2.3 to coil. All intrinsically disordered regions as defined by the region argument are biased with 2.1. Regions always superseded the other assignments.

[pckroon]

Alright, that looks cool.

[Lp0lp]

@Lp0lp Last but not least I'm refactoring the water biasing a little, such that it can be applied without the Go model. I think that's also the original idea of the IDPs.

Perfect! Would be helpful to also be able to combine the water biasing with elastic networks for peptides, etc.

One question I have is whether you think someone would only want to bias IDRs in a structured protein and NOT the rest of the protein (i.e. applying idp fix only). My idea was a CLI like so:

martinize2 -water-bias -> applies the auto water bias
martinize2 -water-bias-eps -> you can set the epsilon for all sec-struc types including an "idr" key for IDPs/IDRs
martinize2 -id-regions -> define which regions are IDRs

In this scenario, any IDR biasing overwrites the default auto biasing. If we have default setting for the sec-struc biases the only tricky part would be to only apply biases to the IDP regions because that would require setting the specific biases to 0.

Could be a possibility for some I guess... If I understood correctly you define 4 assignments (helix, coil, strand and idr) with specific eps' and where idr overrides the other 3. What is there to stop you from setting helix coil and strand to 0 with -water-bias-eps and leaving only the IDR with a eps != 0? That would let users only apply the bias to the IDP region, no? Bit ugly CLI wise and writes a few useless BB-W interactions with a bunch of (harmless, I think) zeros though.

[pckroon]

I'll continue reviewing later today

[pckroon]

I went over the code again and still found some things. Most are very small, some are really just nitpicks.
I also noticed some of my old comments didn't get addressed or resolved, could you also go back and have a look at those.

As for the signature of the GoPipeline, leave it as is for now.

[pckroon]

Did we ever reach a conclusion on this? My original comment mentioned generating the _path s based on moltypes. IIRC the outcome was that we can't do that since things like nbparams are system-wide attributes.

[pckroon]

Nice work!