niedziela, czerwca 29, 2008

More "binding site apples"

It took me more time then I've expected, but finally I've updated my previous script. And here it is.

bp_surface.py(Toggle Plain Text) 



# -*- coding: utf-8 -*-

# -------------------------------------------------------------------------------
# bp_surface.py - Binding pocket surface generator for Pymol version 2.0
# -------------------------------------------------------------------------------

# Copyright (C) 2008 by Lightnir - lightnir@gmail.com

# This script is free software; you can redistribute it and#or modify
# it under the terms of the GNU Library General Public License as
# published by the Free Software Foundation; either version 2 of the
# License, or (at your option) any later version.

# This script is distributed in the hope that it will be useful,
# but WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
# GNU General Public License for more details.

# You should have received a copy of the GNU Library General Public
# License along with this program; if not, write to the
# Free Software Foundation, Inc.,
# 59 Temple Place - Suite 330, Boston, MA  02111-1307, USA.

from pymol import cmd
from Tkinter import *
from pymol.cgo import *
from math import *
import Pmw
import tkColorChooser

BPSFrame = None

def __init__(self):
    self.menuBar.addcascademenu('Plugin', 'LightnirsPlugins', 'Lightnirs Plugins',
                                    label='Lightnirs Plugins')
    self.menuBar.addmenuitem('LightnirsPlugins', 'command',
                                 'BPsurface',
                                 label='BP surface',
                                 command = lambda s=self: BPSurface(s))
class BPSurface:
    # Class variables
    selections = []
    cv = IntVar()
    pp = IntVar()
    cp = IntVar()
    rv = IntVar()
    colorize = IntVar()
    Default_Color = None
    
    res_states = [['Hydrophobic','Gray',0], ['Aromatic','White',0], ['Polar','Green',0], ['Positive','Red',0], ['Negative','Blue',0], ['Cysteine','Yellow',0], ['Proline','Magenta',0]]

    def __init__(self,app):
        global BPSFrame

        if BPSFrame==None:
            BPSFrame = Toplevel(width=400, height=300)              # create a root window
            BPSFrame.title("BP surface editor")                     # Set window title
            BPSFrame.resizable(0, 0)                                # Disable window resize
            BPSFrame.geometry('-40+40')                             # Set window placement
            self.get_selections()
            self.Default_Color = BPSFrame.cget("bg")
            # Create Ligand selection Combobox
            BPSFrame.cb1 = Pmw.ComboBox(BPSFrame,
                label_text = 'Ligand:',
                labelpos = 'nw',
                scrolledlist_items = BPSurface.selections
            )
            BPSFrame.cb1.grid(row=0, column=0,columnspan=2)

            # Create binding pocket Combobox
            BPSFrame.cb2 = Pmw.ComboBox(BPSFrame,
                label_text = 'Binding pocket:',
                labelpos = 'nw',
                scrolledlist_items = BPSurface.selections
            )
            BPSFrame.cb2.grid(row=1, column=0,columnspan=2)

            # Create Z axis Slider
            BPSFrame.slider1 = Scale(BPSFrame, from_=-20, to=20, command=self.scaling, tickinterval=10, resolution=0.5, label='Z axis cutoff', orient=HORIZONTAL)
            BPSFrame.slider1.grid(row=2, column=0, columnspan=2, sticky=N+W+S+E,)

            # Create X axis Slider
            BPSFrame.slider2 = Scale(BPSFrame, from_=-90, to=90, command=self.scaling, tickinterval=30, resolution=1, label='X axis rotation',orient=HORIZONTAL)

            # Create Y axis Slider
            BPSFrame.slider3 = Scale(BPSFrame, from_=-90, to=90, command=self.scaling, tickinterval=30, resolution=1, label='Y axis rotation', orient=HORIZONTAL)

            # Add Checkboxes
            BPSurface.cv.set(1)
            cbx = Checkbutton(BPSFrame, text="Create new object", command=None, variable=BPSurface.cv, onvalue="1", offvalue="0")
            cbx.grid(row=5,column=0, columnspan=2, sticky=W)
            self.CreatePlane()
            
            BPSurface.pp.set(1)
            cbx2 = Checkbutton(BPSFrame, text="Plane preview", command=self.TogglePlane, variable=BPSurface.pp, onvalue="1", offvalue="0")
            cbx2.grid(row=6,column=0, columnspan=2, sticky=W)

            BPSurface.cp.set(0)
            cbx3 = Checkbutton(BPSFrame, text="Use custom plane", command=self.ToggleCustomPlane, variable=BPSurface.cp, onvalue="1", offvalue="0")
            cbx3.grid(row=7,column=0, columnspan=2, sticky=W)

            BPSurface.rv.set(0)
            cbx4 = Checkbutton(BPSFrame, text="Make residue selections", command=self.ToggleResidueSelect, variable=BPSurface.rv, onvalue="1", offvalue="0")
            cbx4.grid(row=8,column=0, columnspan=2, sticky=W)

            # Add some buttons
            btn1=Button(BPSFrame, text='Refresh', padx=0, pady=0, command = self.refresh_list)
            btn1.grid(row=9, column=0, sticky=N+W+S+E, padx=0, pady=0,ipadx=0,ipady=0)
            btn2=Button(BPSFrame, text='Clear', padx=0, pady=0, command = lambda r=1,d=1: self.clear_flags(r,d))
            btn2.grid(row=9, column=1, sticky=N+W+S+E, padx=0, pady=0,ipadx=0,ipady=0)
            btn3=Button(BPSFrame, text='Apply', padx=0, pady=0, command = self.ok)
            btn3.grid(row=10, column=0, sticky=N+W+S+E,)
            btn4=Button(BPSFrame, text='Close', padx=0, pady=0, command = self.myHide)
            btn4.grid(row=10, column=1, sticky=N+W+S+E)

            # Residue selection column
            # Create frame for residue content
            BPSFrame.Res = Frame(BPSFrame)

            # Create Colorize button
            self.colorize.set(1)
            col_cbx = Checkbutton(BPSFrame.Res, text="Colorize residues", command=self.ToggleColor, variable=self.colorize, onvalue=1, offvalue=0)
            col_cbx.grid(row=0,column=0, sticky=N+W)

            # Create group widget
            BPSFrame.gr = Pmw.Group(BPSFrame.Res, tag_text='Residues')
            BPSFrame.gr.grid(row=1, column=0, padx=3, pady=3, sticky=N+W+E)

            # Create residues buttons and checkboxes
            BPSFrame.ccb = []
            for i in range(0,len(self.res_states)):
                var = IntVar()
                var.set(1)
                BPSFrame.ccb.append(Button(BPSFrame.gr.interior(), text='', bg=self.res_states[i][1], activebackground=self.res_states[i][1], relief=RAISED, overrelief=RAISED, height=1, width=2, padx=0, pady=0, command=lambda nr=i: self.changecolor(nr)))
                BPSFrame.ccb[i].grid(row=i,column=0, sticky=N+W)
                rs_cbx = Checkbutton(BPSFrame.gr.interior(), text=self.res_states[i][0], variable=var, onvalue="1", offvalue="0")
                rs_cbx.grid(row=i,column=1, sticky=N+W)
                self.res_states[i][2]=var
            undo1 = Button(BPSFrame.Res, text='Undo all', padx=0, pady=0, command = lambda u="all":self.undo(u))
            undo1.grid(row=2,column=0, sticky=N+W+E)
            undo2 = Button(BPSFrame.Res, text='Undo last', padx=0, pady=0, command = lambda u="last":self.undo(u))
            undo2.grid(row=3,column=0, sticky=N+W+E)
            # Define own colors
            self.create_bp_colors()
            # Make a backup session
            cmd.save('bp_undoall.pse', '(all)')

            # create callback to prevent window kill
            BPSFrame.protocol("WM_DELETE_WINDOW", self.myHide)
            BPSFrame.mainloop()
        else:
            if BPSFrame.state() == "normal":
                self.myHide()
            elif BPSFrame.state() == "withdrawn":
                self.myShow()


    def ToggleColor(self):
        if self.colorize.get()==1:
            for i in range(0,len(self.res_states)):
                BPSFrame.ccb[i].config(state=ACTIVE, bg=self.res_states[i][1], activebackground=self.res_states[i][1], relief=RAISED, overrelief=RAISED)
        else:
            for i in range(0,len(self.res_states)):
                BPSFrame.ccb[i].config(state=DISABLED, bg=self.Default_Color , activebackground=self.Default_Color, relief=RAISED, overrelief=RAISED)

    def changecolor(self,nr):
        new_color = tkColorChooser.askcolor(self.res_states[nr][1])
        rgb_color = map(lambda color:  '%.4f'%float(color/float(65536)) ,new_color[0])
        self.res_states[nr][1]=new_color[1]
        print rgb_color,new_color[1]
        cmd.set_color('bp_plugin_color_'+str(nr),new_color[0],0)
        BPSFrame.ccb[nr].config(bg=new_color[1], activebackground=new_color[1])

    def create_bp_colors(self):
        for i in range(0,len(self.res_states)):
            rgb_color = map(lambda color:  '%.4f'%float(color/float(65536)) ,BPSFrame.winfo_rgb(self.res_states[i][1]))
            cmd.set_color('bp_plugin_color_'+str(i),rgb_color)

    def state(self):
        return map(lambda var: var[2].get(), self.res_states)

    def undo(self,u):
        if u=="all":
            cmd.load('bp_undoall.pse')
        if u=="last":
            cmd.load('bp_undolast.pse')
        self.TogglePlane()

    def scaling(self,v):
        self.TogglePlane()

    def myShow(self):
        self.refresh_list()
        BPSFrame.deiconify()

    def myHide(self):
        cmd.delete("plane_preview")
        BPSFrame.withdraw()

    def get_selections(self):
        BPSurface.selections = []
        for item in cmd.get_names("all"):
            if cmd.get_type(item)=="object:molecule":
                BPSurface.selections.append(item)
            if cmd.get_type(item)=="selection":
                if item[0]<>"_":
                    BPSurface.selections.append(item)

    def refresh_list(self):
        if BPSurface.pp.get()==1:
            self.CreatePlane()
        # Save old values
        old = BPSurface.selections
        try:
            old_cb1=BPSFrame.cb1.getvalue()[0]
        except:
            old_cb1=''
        try:
            old_cb2=BPSFrame.cb2.getvalue()[0]
        except:
            old_cb2=''
        self.get_selections()

        # Update selectionlist if needed
        if not(old==BPSurface.selections):
            BPSFrame.cb1.setlist(BPSurface.selections)
            BPSFrame.cb2.setlist(BPSurface.selections)

        # Is ligand selection still available?
        if not(old_cb1 in BPSurface.selections):
            # No - clear the combobox entry
            BPSFrame.cb1.component('entryfield').clear()
        else:
            # Yes - set the old value
            BPSFrame.cb1.selectitem(old_cb1,setentry=1)

        # Is binding pocket selection still available?
        if not(old_cb2 in BPSurface.selections):
            # No - clear the combobox entry
            BPSFrame.cb2.component('entryfield').clear()
        else:
            # Yes - set the old value
            BPSFrame.cb2.selectitem(old_cb2,setentry=1)

    def clear_flags(self,r,d):
        cmd.flag('ignore','all','clear')
        if r==1:
            cmd.rebuild('all')
        if d==1:
            cmd.delete('bp_cutoff')
            cmd.delete('bp_surface')
            cmd.delete('bp_hydrophobic')
            cmd.delete('bp_aromatic')
            cmd.delete('bp_negative')
            cmd.delete('bp_positive')
            cmd.delete('bp_proline')
            cmd.delete('bp_cysteine')
            cmd.delete('bp_polar')

    def CreatePlane(self):
            view = cmd.get_view()
            # choose the size of the plane
            plane_size = abs(view[11])/4.0
            obj = []

            # Create a plane in camera space
            plane = [
                [ -plane_size,  plane_size, BPSFrame.slider1.get() ],
                [  plane_size,  plane_size, BPSFrame.slider1.get() ],
                [ -plane_size, -plane_size, BPSFrame.slider1.get() ],
                [  plane_size, -plane_size, BPSFrame.slider1.get() ]]

            # Using custom plane?
            if BPSurface.cp.get()==1:
                sinx = sin(pi*BPSFrame.slider2.get()/180)
                siny = sin(pi*BPSFrame.slider3.get()/180)
                cosx = cos(pi*BPSFrame.slider2.get()/180)
                cosy = cos(pi*BPSFrame.slider3.get()/180)
                # Map XY axis rotation
                # cosY  sinX*sinY   -cosX*sinY
                # 0     cosX        sinX
                # sinY  -sinX*cosY  cosX*cosY
                plane = map( lambda p: [
                cosy*p[0] + sinx*siny*p[1] + (-cosx*siny*p[2]),
                cosx * p[1] + sinx* p[2],
                siny * p[0] + (-sinx*cosy*p[1]) + cosx*cosy*p[2]
                ], plane )

            normal = [ 0.0, 0.0, 1.0 ]

            # Transform plane coordinates into model space
            plane = map( lambda p,v=view: [
                v[0] * p[0] + v[1] * p[1] + v[2]* p[2],
                v[3] * p[0] + v[4] * p[1] + v[5]* p[2],
                v[6] * p[0] + v[7] * p[1] + v[8]* p[2]
            ], plane )

            normal = apply( lambda p,v=view:[
                v[0] * p[0] + v[1] * p[1] + v[2]* p[2],
                v[3] * p[0] + v[4] * p[1] + v[5]* p[2],
                v[6] * p[0] + v[7] * p[1] + v[8]* p[2]
            ], (normal,) )

            # Transform position relative to the camera
            plane = map( lambda p,v=view: [
                p[0] + v[9 ] + v[12],
                p[1] + v[10] + v[13],
                p[2] +       + v[14],
            ], plane )
            obj.extend( [ COLOR, 0.5, 0.5, 1.0 ] ) # blueish
            obj.extend( [ BEGIN, TRIANGLE_STRIP ] )
            obj.append( NORMAL )
            obj.extend( normal )
            # draw the plane
            for a in plane:
                obj.append( VERTEX)
                obj.extend(a)
            obj.append( END )
            # delete existing object (if any)
            cmd.delete("plane_preview")
            # now load the new object without zooming
            cmd.set('auto_zoom', 0, quiet=1)
            auto_zoom = cmd.get('auto_zoom')
            cmd.load_cgo(obj,'plane_preview')
            cmd.set('auto_zoom', auto_zoom, quiet=1)

    def TogglePlane(self):
        if BPSurface.pp.get()==0:
            cmd.delete("plane_preview")
        if BPSurface.pp.get()==1:
            self.CreatePlane()

    def ToggleCustomPlane(self):
        self.TogglePlane()
        if BPSurface.cp.get()==0:
            BPSFrame.slider2.grid_remove()
            BPSFrame.slider3.grid_remove()
        if BPSurface.cp.get()==1:
            BPSFrame.slider2.grid(row=3, column=0, columnspan=2, sticky=N+W+S+E,)
            BPSFrame.slider3.grid(row=4, column=0, columnspan=2, sticky=N+W+S+E,)

    def ToggleResidueSelect(self):
        if BPSurface.rv.get()==1:
            BPSFrame.Res.grid(row=0, column=2, rowspan=10, sticky=N+W+E)
        if BPSurface.rv.get()==0:
            BPSFrame.Res.grid_remove()
    
    def ok(self):
        # Has cb1 a value?
        sel1 = BPSFrame.cb1.getvalue()
        if len(sel1) == 0:
            dialog = Pmw.MessageDialog(BPSFrame,
                title = 'Error',
                defaultbutton = 0,
                buttons = ('OK',),
                message_text = 'ERROR: No ligand selection. Please\n choose one from the dropdownlist.'
            )
            dialog.iconname('Selection Error')
            dialog.bell()
            dialog.activate()
            return 0

        # Has cb2 a value?
        sel2 = BPSFrame.cb2.getvalue()
        if len(sel2) == 0:
            dialog = Pmw.MessageDialog(BPSFrame,
                title = 'Error',
                defaultbutton = 0,
                buttons = ('OK',),
                message_text = 'ERROR: No binding pocket selection. Please\n choose one from the dropdownlist.'
            )
            dialog.iconname('Selection Error')
            dialog.bell()
            dialog.activate()
            return 0

        # Are cb1 and cb2 the same?
        if sel1[0]==sel2[0]:
            dialog = Pmw.MessageDialog(BPSFrame,
                title = 'Error',
                defaultbutton = 0,
                buttons = ('OK',),
                message_text = 'ERROR: Ligand and binding pocket must have different selections'
            )
            dialog.iconname('Selection Error')
            dialog.bell()
            dialog.activate()
            return 0

        # Set orgin and normal vector
        cmd.origin(sel1[0])

        vector = [[  0,  0,  1]]

        # Using custom plane?
        if BPSurface.cp.get()==1:
            sinx = sin(pi*BPSFrame.slider2.get()/180)
            siny = sin(pi*BPSFrame.slider3.get()/180)
            cosx = cos(pi*BPSFrame.slider2.get()/180)
            cosy = cos(pi*BPSFrame.slider3.get()/180)
            # Remap the normal vector to custom
            vector = map( lambda p: [
            cosy*p[0] + sinx*siny*p[1] + (-cosx*siny*p[2]),
            cosx * p[1] + sinx* p[2],
            siny * p[0] + (-sinx*cosy*p[1]) + cosx*cosy*p[2]
            ], vector )

        # Remap the normal vector to the model space
        d=cmd.get_view(0)
        vector = map( lambda p,v=d: [
            v[0] * p[0] + v[1] * p[1] + v[2]* p[2],
            v[3] * p[0] + v[4] * p[1] + v[5]* p[2],
            v[6] * p[0] + v[7] * p[1] + v[8]* p[2]
            ], vector )

        # Calculate plane coefficients
        # Ax + By + Cz + D = 0
        # D = -Ax -By -Cz
        A = vector[0][0]
        B = vector[0][1]
        C = vector[0][2]
        D = -vector[0][0]*d[12] - vector[0][1]*d[13] - vector[0][2]*d[14]

        atoms = cmd.get_model(sel2[0])

        # Surface stuff
        self.clear_flags(0,1)
        for at in atoms.atom:
            if (A*(at.coord[0])+B*at.coord[1]+C*at.coord[2]+D-BPSFrame.slider1.get()>0):
                if not('bp_cutoff' in cmd.get_names("selections")):
                    cmd.select('bp_cutoff','index '+str(at.index))
                else:
                    cmd.select('bp_cutoff','bp_cutoff or index '+str(at.index))
        if ('bp_cutoff' in cmd.get_names("selections")):
            if BPSurface.cv.get()==0:
                cmd.flag('ignore','bp_cutoff','reset')
                cmd.delete('indicate')
                cmd.deselect()
                cmd.show('surface',sel2[0])
            else:
                self.clear_flags(0,0)
                cmd.hide('surface',sel2[0])
                cmd.set('auto_zoom','off')
                cmd.create('bp_surface',sel2[0]+' and not bp_cutoff')
                cmd.set('auto_zoom','on')
                cmd.show('surface','bp_surface')
            cmd.rebuild(sel2[0])
        # Selections stuff
        if BPSurface.rv.get()==1:           # make residue selections
            if self.state()[0]==1:          # select hydrophobic
                cmd.select('bp_hydrophobic', 'byres ('+sel2[0]+' and resn ala+gly+ile+leu+met+val)')
            if self.state()[1]==1:          # select aromatic
                cmd.select('bp_aromatic', 'byres ('+sel2[0]+' and resn phe+trp+tyr)')
            if self.state()[2]==1:          # select polar
                cmd.select('bp_polar', 'byres ('+sel2[0]+' and resn his+asn+gln+ser+thr)')
            if self.state()[3]==1:          # select positive
                cmd.select('bp_positive', 'byres ('+sel2[0]+' and resn lys+arg)')
            if self.state()[4]==1:          # select negative
                cmd.select('bp_negative', 'byres ('+sel2[0]+' and resn asp+glu)')
            if self.state()[5]==1:          # select cysteine
                cmd.select('bp_cysteine', 'byres ('+sel2[0]+' and resn cys)')
            if self.state()[6]==1:          # select proline
                cmd.select('bp_proline', 'byres ('+sel2[0]+' and resn pro)')
        cmd.deselect()
        cmd.save('bp_undolast.pse', '(all)')
        if (BPSurface.rv.get()==1 and BPSurface.colorize.get()==1): # colorize?
            if self.state()[0]==1:          # hydrophobic
                cmd.color('bp_plugin_color_0','bp_hydrophobic')
            if self.state()[1]==1:          # select aromatic
                cmd.color('bp_plugin_color_1','bp_aromatic')
            if self.state()[2]==1:          # select polar
                cmd.color('bp_plugin_color_2','bp_polar')
            if self.state()[3]==1:          # select positive
                cmd.color('bp_plugin_color_3','bp_positive')
            if self.state()[4]==1:          # select negative
                cmd.color('bp_plugin_color_4','bp_negative')
            if self.state()[5]==1:          # select cysteine
                cmd.color('bp_plugin_color_5','bp_cysteine')
            if self.state()[6]==1:          # select proline
                cmd.color('bp_plugin_color_6','bp_proline')


The major updates of bp_surface.py version 2.0 are:
  • Ligand and binding site can now be either a selection or an object
  • The orientation of the cutting plane now can be customized
  • Binding site residues can be selected by type: hydrophobic, aromatic, polar, charged positive or negative, proline, cysteine
  • Selected residues can be colorized
  • If you screw up the color settings you can undo the last action or recover the state of the binding site when the script was started the first time (use it with caution)
Here's how the new GUI looks in full view. I like it simple so you only see this window if you check the "Use custom plane" and "Make residues selections". Otherwise parts of the window are collapsed. It's also highly customizable and self-explaining. Lets say I would like to color all polar residues in the binding pocket of 1P93 in red, then I have to check the "Make residues selections" and "Colorize residues" checkboxes, uncheck all types of residues except "polar" and set the color to red by clicking the color button to the left of the "polar" checkbox. In addition I found polar contacts with the ligand and binding pocket residues and hided the generated surface. This looks like this:

Same molecule (Dexamethasone) but this time hydrophobic residues were colored orange. Notice that large parts of the binding site surface are hydrophobic, which isn't a surprise since the bounded ligand is a steroid.

I hope this tool will be helpful ]:)
Autor: o
Etykiety: , , ,

4 komentarze:

Felix pisze...

hi, very cool and great pictures! i am looking forward to using it when I have some time. I already thought of the system

1:38 PM
Lightnir pisze...

Thanks for the feedback. You point it out - time flies by.

Wunderliches Wort "Die Zeit vertreiben"!
Sie zu halten wäre das Problem.
~ R.M. Rilke

I would also like to play a bit with it and make some code optimization. Unfortunately I have no time for it right now since in one week I'll have to defend my Magister (Master of Science) degree. After that I will look up into it again hopefully with some fresh ideas.

6:24 PM
Felix pisze...

good luck for your defense.
for me it's actually going to be about another year until I can get my Magister.

3:27 PM
Vladimir Chupakhin pisze...

Thank's a lot!

10:49 PM

Prześlij komentarz

Subskrybuj: Komentarze do posta (Atom)