CHARMM Development Project Forums User Discussion & Questions Setup, I/O, and Basic questions Selecting only lipids under protein

Hi Charmm users,

I tried to come up with a smart selection to only select the lipids under my protein. Since the protein moves on the X/Y-Plane I had the following plan:

Get the center of mass for my protein, extract the X/Y-Coordinates and select all lipids within 15A of these Coordinates. Is this possible within one selection?

I am interested in the dynamics, density, overall charge, etc.. of the lipids below the protein compared to the lipids outside the selection.

The selection would have to be updated every step. Is this possible with a selection or do I have to iterate every step and save which lipids are at which step under my protein and then calculate ${WHATEVER} separately for this specific residue at the specific steps?

I hope I could make myself clear.

Thanks for any suggestions.
Cheers
Bjoern

Something like this snippet should do what you are asking about (but I am not sure that your criterion for lipids being under the protein is optimal):

label loop
traj read
coor stat sele segid prot end
define lipunder sele segi lipi .and. point ?XAVE ?YAVE ?ZAVE CUT 15 end
.
.
.
goto loop

For time series analysis involving the notion of nearby molecules, the shell feature should be useful; see shell.doc and correl.doc for more info.

I agree, for many analyses of solvent behavior around a solute the existing tools (COOR ANALysis or SHELL) do a very good job. For the kinds of analyses indicated by the original poster a CHARMM loop is likely to be the better solution; in particular since the concept of "under" is not so clearly tied to "being close to".

True, but in many contexts, such as the above snippet, atom selection is not always "image aware" (shell code is), which may be an issue. Depending on the cell size, a 15 A shell based on lipid P and N atoms (assuming PC lipid) may be fairly reasonable.

For a loop over frames, the image cell boundary problem could also be handled by a coordinate translation in the xy plane followed by an image update, but that would slow things down a bit more.

I often try several different approaches in this situation (something to do with eggs and baskets ).

Looks promising, thanks a lot!
I would have to add an offset to ZAVE or increase the CUT since the max diameter in the current selection will be somewhat over the membrane interface.




This is sadly so true, since the protein has a "hole" in the middle an if I would use nearby those lipids under the center of protein would most likely be excluded. Maybe this ASCII-Art can clear things up:





Thanks for pointing out this problem.
Could you give an code example for"coordinate translation in the xy plane". I know the purpose behind the translation and in fact I probably would have run into this pitfall, since my protein is is moving around and reaches the end of the box here and there so I guess I have to deal with the boundary conditions.

How much would this slow things down? It's going to be a very long simulation (so far ~400ns but still crunching to reach some microseconds and I am dreaming of milliseconds)
The simulations are consecutive runs of 1ns slices. So for each unreduced trajectory (1ns, 500000 steps, write out freq. 500, integration time 2 femtoseconds) has the size of 1.5G.

Sadly, running multiple small simulations is not an option here

If the slowdown is way to big I maybe should go with a more inaccurate selection (COOR ANALysis or SHELL) in favor for a reasonable time.

I think I'll play a bit around with both options and see what happens.

Thanks a lot!
Bjoern

You ought to be able to use MERGE RECEnter to get the protein in the middle of the box in all frames.
To save diskspace you could increase the write out frequency, say to 2500. This still gives you 200 frames/ns which is plenty for such a long simulation.

I am not entirely familiar how recenter in detail works, but I only used it to eliminates translation of the entire system. But since the Protein is moving relative to the membrane this would have no effect on the protein moving around on the membrane. Please correct me if understood recentering wrong. Another concern I do have regarding the recentering on the protein. If I should choose to recenter on the protein will this not influence the lipids I see passing/resting/following below the protein, will it? I am confused and I couldn't dissipate my doubt by reading the RECEnter notes in dynam.doc.


Well, I am not sure about this, since I once ran into trouble with a too low resolution (1ns) while calculating relaxation times. Comparing the simulation to the same setup with 100ps time resolution yielded different results. So far we didn't investigate any further but the difference might be related to the time resolution. So I'd rather throw in some more discs then throw away time by recalculating with a higher time resolution.

Deleting trajectories and storing reduced ones is quick and always possible ;-)

Thanks for your help and suggestions.
I appreciate it very much!

Cheers
Bjoern

I hesitate to use MERGE COOR RECENTER, since the the protein will be shifted in z as well as x and y, which may not be desirable. Also, if (1) the bilayer remains centered around z=0 during the simulation, and (2) the distance between the protein an bilayer surface is variable, you may wish to use a fixed point in z instead of ?ZAVE, which is the center of the protein. If the protein is always above (+z) the bilayer, ?ZMIN is a better choice as well. For the x,y shift--

traj read
coor stat sele segid prot end
coor tran xdir -?XAVE ydir -?YAVE
update imgfrq 1 ...
define lipunder sele segi lipi .and. point 0. 0. ?ZAVE CUT 15 end

Well, this looks like a plan.

Thanks a lot for your suggestions!

Cheers
Bjoern