rendering zdock server results using python and UCSF Chimera

Problem:
I was interested in learning where the protein-protein interaction sites were on a particular protein (the receptor). I had pdb files for that protein, and for three proteins that it is known to bind to (the ligands). There are specific amino-acids of interest on the receptor protein, where there is variability among different species. For every combination of receptor and ligand, I want to perform a prediction of where the proteins interact, and then generate an image where the variable amino-acids are highlighted.



Strategy:
Use I-TASSER to generate pdb files for all of the proteins of interest (two isoforms of the receptor, three different ligands). Submit the pdb files to ZDOCK in receptor-ligand pairs (six runs total). Use Chimera scripts to render each combination, highlighting the proper residues in the rendering.

Procedure:
Submitting sequences to I-TASSER and ZDOCK is straightforward. For I-TASSER results, I usually download the tar.bz file, then use 7-zip to extract it. I then typically just use the file called model1.pdb for downstream analyses.

For the ZDOCK submission, I did not select any specific residues to be part of the binding pocket, or excluded from the binding pocket.

To download the ZDOCK results, I just click the four links in the "Download Files" box. An additional file that is of great interest to me is the representation of the the centers of mass of the top 500 docking positions as spheres in space. There isn't a direct link to that file. Instead, you have to enter the url of the file directly into the address bar of your browser. The file itself is called ligcms.pdb. To obtain it, just append that name to the end of the url for your results. For example, one of my results urls was: http://zdock.umassmed.edu/results/88cedc8055/
So to get the file of centers of mass, I used http://zdock.umassmed.edu/results/88cedc8055/ligcms.pdb

Instead of manually highlighting the same residues in every rendering, I'm using a Chimera Python script to do the work for me.

Chimera can be executed from the command line, and given a python script as an argument. For example: chimera myscript.py

Chimera scripting is a little quirky (or I'm just bad at it, probably both). Here are a few quirks I had trouble with:
To change the color of something, you must first get a color object representing the color you want. To get a color object, you must use the convertColor function from the Commands module. For example: red = Commands.convertColor('red')

Amino acid residues in a peptide are zero-indexed, rather than one-indexed as is the convention for biologists. The first amino acid is called amino acid 0, not 1. So when you use the model.findResidue(x) command, be aware that x should be 1 less than the conventional position.

Chimera has some trouble reading the ZDOCK pdb files. Fix that by removing everything after the x y z columns.

There are two distinct syntaxes for executing commands in Chimera Python scripts, one is through the Python API. The other is using command line functions. My impression is that the Python API is more powerful, but it is also exceedingly poorly documented (as far as I can tell). The command line functions are well documented, and can be accessed through Python via the chimera.runCommand function. In the code below, I use a mix of both kinds of commands.

Code:

import os
import chimera
from chimera import runCommand as rc # use 'rc' as shorthand for runCommand
from chimera import replyobj # for emitting status messages
import Commands
from glob import glob

#Note: the columns after XYZ in the pdb file for the receptor from ZDOCK must be removed before chimera can read the file
#This is done automatically in the script loop, later


zdock_path="E:/data/docking/zdock"

red = Commands.convertColor('red')
sky_blue = Commands.convertColor('sky blue')
res_to_highlight = [3,11,14,16,32,54,55,57,69,73,85,95] #1-indexed 
res_to_highlight = [x - 1 for x in res_to_highlight] #residues in the Chimera model are 

for dir in glob(zdock_path + "/*"):
  out_name = os.path.basename(dir)
  rc("background solid white")
  
  #remove columns from pdb files that cause an error.
  with open(os.path.join(dir,'Fd.pdb'), "r") as infile:
    with open(os.path.join(dir,'Fd_truncated.pdb'), "w") as outfile:
      for line in infile:
        outfile.write(line[0:58]+"\n")
  receptor = chimera.openModels.open(os.path.join(dir,'Fd_truncated.pdb'),type="PDB")
  cloud = chimera.openModels.open(os.path.join(dir,'ligcms.pdb'),type="PDB")

  m=receptor[0]

  m.color=sky_blue
  for x in res_to_highlight:
    m.findResidue(x).ribbonColor = red
  
  #rc("savepos original")
  rc("turn y -45")
  rc("turn x -100")
  if "minor" in dir:
    rc("turn z 60")
  #rc("~focus")
  rc("window")
  #rc("scale 0.8")
  rc("unset depthCue")

  rc("copy file " + out_name + "_zoom.png" + " width 3300 supersample 3")
  #rc("reset original")
  rc("close session") 

Result:
 

References:

Quick introduction to scripting Chimera. Shows how to open files, and save images.
https://www.cgl.ucsf.edu/chimera/docs/ProgrammersGuide/basicPrimer.html

In chimera the camera does not rotate, it always points in the -z direction. It can pan and zoom and focus, but it cannot rotate. To get the effect of rotating the camera, you must instead rotate all of the models. This can be done using the python interface, or using the command line "turn" function.
https://www.cgl.ucsf.edu/chimera/docs/ProgrammersGuide/faq.html#q4


An index of all of the commands that can be used through the Chimera command line, or through the chimera.runCommand function in python.
https://www.cgl.ucsf.edu/chimera/current/docs/UsersGuide/framecommand.html