From: Adelman, Joshua Lev (jla65_at_pitt.edu)
Date: Thu Mar 19 2015 - 08:29:51 CDT

Hi Nicolas,

I'm still a little unclear about your counting method for the description you provided, but I think you might want to be a little careful. You can't really get accurate transitions from just looking at frames n and n - 1 in isolation. It is more typical to look at the entire history of each water molecule and split your simulation box into 3 sections using 2 planes (lower bulk, channel, upper bulk). You only want to count transitions that go from the upper bulk to the lower bulk (or lower to upper), that go through the channel. It's not sufficient to look only at the crossings between pairs of regions, due to either periodic boundary crossings or movement between the channel and bulk that don't represent full transitions. We and others have used this method to look at permeation counts of waters and other molecules through channels and pores. An example of a detailed description of this is found in:

http://www.ncbi.nlm.nih.gov/pmc/articles/PMC3984995/

and there is a discussion of the errors involved in the calculation due to finite time resolution in the trajectory in the supporting materials.

While it's possible to do this using TCL scripting, I have always preferred to do this analysis outside of VMD, using a combination of Python and low-level c code. I can make this available to you, but I haven't had the time to package it in a clean way.

There is also a gromacs tool that you might be able to use by Oliver Beckstein that already has this method coded:
https://github.com/orbeckst/g_count

Best wishes,
Josh

On Mar 19, 2015, at 4:04 AM, nicolas martin wrote:

Dear all,

I come to you with a problem which is as much related to scripting than it is to chemistry.

I'm have a long trajectory of a protein channel in explicit solvent. Since I'm watching the transitions from opened to closed configuration of this channel I'd like to measure the quantity of water going through along time.

I've figured out that counting the number of molecules that are up to the constriction point at frame n-1 and down to this constriction point at frame n or the other way around would give me the answer. But it's impossible (and pretty useless) to do for the water molecules of the whole system. so I came up with the idea to consider waters inside the channel only.

It sounds really complicated to me to create an array which would contain the idea and z coordinated or the waters selected at frame n-1 and frame n so I can compare and see if a transition happens.

So here are my questions. Are you aware of a script which could do that for me? Do you think my approach is meaningful ? Do you think there are better/faster/easier way to do it?

Thanks a lot for your help,

NM

--
Nicolas Martin
PhD candidate at ISIS
67000 Strasbourg
06.29.73.12.73