search

giorgiomaccari / mdanalysis

Automatically exported from code.google.com/p/mdanalysis
0 stars 0 forks source link

AtomSelect using pre-existing AtomGroup properties #42

Closed GoogleCodeExporter closed 9 years ago

GoogleCodeExporter commented 9 years ago 
from MDAnalysis import Universe
u = Universe(pdb, xtc) 
PO4 = u.AtomSelect("name PO4")
resPO4 = u.Atomselect( "byres  PO4")

## However there is no way of selecting properties from another predefined ## 
AtomGroup to use in the Atomselect command.

## Normally this wouldn't be a problem by combining the two - i.e.

resPO4 = u.AtomSelect("byres name PO4")

## However I am using LeafletFinder which returns two AtomGroups
## containing the PO4 atoms in each leaflet. I then want to select 
## the residues from the universe corresponding to the residues from each 
## AtomGroup...

L = MDAnalysis.analysis.leaflet.LeafletFinder(u,"name PO4") 
L.atoms(0)  ### Atomgroup containing PO4 atoms in leaflet 0 
L.atoms(1)  ### Atomgroup containing PO4 atoms in leaflet 1
res0 = u.Atomselect( "byres  L.atoms(0)" )

### Is there a way that Atomgroups can be defined and referenced for use
### within the AtomSelect command. This would also be useful for selecting
### between universes - i.e. selections in one universe based on an
### AtomGroup in another.

Original issue reported on code.google.com by joseph.g...@gmail.com on 20 Sep 2010 at 1:51

GoogleCodeExporter commented 9 years ago

Original comment by joseph.g...@gmail.com on 20 Sep 2010 at 1:51

GoogleCodeExporter commented 9 years ago 
This is related to Issue 10.

Maybe we can introduce "named selections"?

  # ... as above
  L0 = L.atoms(0)        # produces an AtomGroup
  L0.name = "leaflet_0"  # assign a name (or use the CHARMM-like "L0.define"?)
  res0 = u.selectAtoms("byres [leaflet_0]")   # new syntax for a named selection

The Universe/AtomGroup.selectAtoms() method should also get a new keyword 
argument name='selection_name'. If this is set (or a non-None value is assigned 
to the AtomGroup.name) then we should store the selection in a *global* 
(weak-ref) dictionary. In this way it will also become possible to do 
selections between different universes.

Original comment by orbeckst on 21 Sep 2010 at 10:41

GoogleCodeExporter commented 9 years ago 
I'm a little unclear on how you would like to do selections between different 
universes. Currently AtomGroups are simply collections of the atom indices in 
the structure file. So you could only consistently apply an AtomGroup from one 
Universe to another if the underlying structure file is the same (ie the same 
number and ordering of atoms between groups - basically between separate 
trajectory files). Is this what you mean?

Original comment by naveen.m...@gmail.com on 28 Sep 2010 at 2:37

GoogleCodeExporter commented 9 years ago 
rom a user's perspective it would make sense to load, say, all solved 
structures of a enzyme, generating a list of universes, and then asking for the 
coordinates of the active site residues in order to compute an 'average 
picture' of the active site:

multiverse = [Universe(pdb, name=pdb) for pdb in pdbfiles]         # name each 
universe
selection = " or ".join(["universe %s" % pdb for pdb in pdbfiles]) # select any 
universe
selA = "resid 31 and ( " + selection + " )"
selB = "resid 111 and ( " + selection + " )"
selC = "resid 391 and ( " + selection + " )"

u = multiverse[0]
resA = u.selecAtoms(selA)    # just use any universe to access multiverse; resA 
= mixed AtomGroup
avgA = resA.coordinates().mean(axis=0)  # average coordinates of active site 
residue A across all pdbs
...

Or: I run simulation A of a protein, and B is a homologous protein with a 
ligand. I fit trajectory B so that the protein superimposes, so everything is 
happening in the same region of space. Now I want to find all residues in A 
that are within a certain distance of a cutoff in order to define the binding 
site in A based on my knowledge of B.

The example is a bit contrived but that's essentially what I mean. I suppose 
the main use would be distance-based selections. But I remember using VMD and 
occasionally being annoyed that I cannot use properties of one molecule to look 
for something in the other.

Example for doing a selection from an arbitrary ATomGroup:

Artificial example:

opep = u.selectAtoms("name OW and around 4.0 name N", name="opep")  # water 
oxygens near peptide bond
wpep = u.selectAtoms("byres selection opep")  # introduce new 'selection' 
keyword to refer to named AtomGroup

and/or

opep = u.selectAtoms("name OW and around 4.0 name N")  # water oxygens near 
peptide bond
wpep = opep.selectAtoms("byres self")   # maybe use special 'self' keyword to 
refer to the group itself

wpep should be an AtomGroup containing all atoms in the water molecules having 
an oxygen with 4.0 Å to any N atom.

This should work for any given AtomGroup, which can be generated without ever 
going through the selection parser.

Original comment by orbeckst on 28 Sep 2010 at 2:48

GoogleCodeExporter commented 9 years ago 
Maybe we leave this issue concerned with selections from arbitrary atom groups 
and open another one for the "multiverse" selections.

Original comment by orbeckst on 28 Sep 2010 at 2:49

GoogleCodeExporter commented 9 years ago 
A workaround to get all residues for a give group of Atoms is to use the new 
AtomGroup.residues attribute, which lists all residues that contain at least 
one of the atoms in the group (available since r488).

from MDAnalysis import Universe
u = Universe(pdb, xtc) 
PO4 = u.AtomSelect("name PO4")
resPO4 = PO4.residues

Original comment by orbeckst on 15 Oct 2010 at 12:56

GoogleCodeExporter commented 9 years ago 
Using the OO hierarchy for groups goes some way to alleviate the problem; 
anything else requires a major investment of coding so we leave it as won't fix 
– if anyone really wants any of the features discussed here then please open 
a new issue with the desired specifications.

Original comment by orbeckst on 13 Sep 2011 at 3:28