[BioPython] protein-ligand interactions
Peter Cock
p.j.a.cock at googlemail.com
Fri Mar 20 13:36:47 UTC 2009
On Fri, Mar 20, 2009 at 12:18 PM, mitlox <mitlox at op.pl> wrote:
> Thank you very much for your code, it works and the output is exactly for
> what I was looking for.
>
> I try to get a structureCA object to write out the results in a PDB file
> (outCA.pdb) like this:
> ATOM 5275 CA ILE A 881 17.242 57.141 22.062 1.00 38.49
> C
> ATOM 5283 CA VAL A 882 16.292 57.880 25.678 1.00 38.90
> C ....
>
> Unfortunately I get this error with ... Here is the code:
> ...
> structureCA = []
> ...
> io=Bio.PDB.PDBIO()
> io.set_structure(structureCA)
> io.save('outCA.pdb')
Your structureCA object is just a python list, containing Residue objects.
Instead you need to create a new object with the partial chain - which
can be done by creating structure, model and chain objects manually.
However, I suggest you re-read pages 5 and 6 of the Bio.PDB
documentation for the recommend approach:
http://biopython.org/DIST/docs/cookbook/biopdb_faq.pdf
In your case, you'll want to write your own selection class using the
residue distance to the ligand. I recognise this might seem rather
complicated for a python novice as you have to create your own
class - so here is my solution:
import Bio.PDB
import numpy
pdb_code = "1E8W"
pdb_filename = "1E8W.pdb"
structure = Bio.PDB.PDBParser().get_structure(pdb_code, pdb_filename)
model = structure[0]
chainA = model["A"]
def residue_dist_to_ligand(protein_residue, ligand_residue) :
"""Returns distance from the protein C-alpha to the closest ligand atom."""
distances = []
for atom in ligand_residue :
diff_vector = protein_residue["CA"].coord - atom.coord
distances.append(numpy.sqrt(numpy.sum(diff_vector * diff_vector)))
return min(distances)
class NearLigandSelect(Bio.PDB.Select):
def __init__(self, distance_threshold, ligand_residue) :
self.threshold = distance_threshold
self.ligand_res = ligand_residue
def accept_residue(self, residue):
if residue == self.ligand_res :
return True #change this to False if you don't want the ligand
else :
dist = residue_dist_to_ligand(residue, self.ligand_res)
return dist < self.threshold
io=Bio.PDB.PDBIO()
io.set_structure(structure)
#From looking at the PDB file, ligand is last residue in chain A
ligand_res = chainA.child_list[-1]
#Going to use a distance theshold of 4A
io.save("near_ligand.pdb", NearLigandSelect(4, ligand_res))
print "Done"
Peter
More information about the Biopython
mailing list