CaviTracer tool for detecting channels, tunnels, pores

[karolamik13]

CaviTracer, created by Eryk Trzcinski, as his MSc work under my supervision.
Additional installation of open3D package is required:
pip install open3d

Test functions:
from prody import *

1. Single PDB:
   p = parsePDB('1tqn.pdb')
   atoms = p.select("protein")
   channels, surface = calcChannels(atoms, output_path='channels_output.pdb')
   showCavities(surface)

To create 3D model to display it in ProDy
vmd_path = '/usr/local/bin/vmd'
model = getVmdModel(vmd_path, atoms)

To display channels (all or one by one):
showChannels(channels, surface=surface, model=model)
showChannels(channels, model=model)
getChannelParameters(channels)
showCavities(surface, show_surface=True)

Save all channels independently as PDB:
output = 'chanels_test'
channels, surface = calcChannels(atoms, output, separate=True)

Display channel 1:

showChannels(channels[1], model)

selected_channels = [channels[1], channels[8]]
showChannels(selected_channels, model)

selected_channels = channels[1:4]
lengths, bottlenecks, volumes = getChannelParameters(selected_channels)
selected_channels_atoms = getChannelAtoms(selected_channels)

atoms_protein = getChannelAtoms(channels, atoms)

Extracting residues that are forming channel[1]:
atoms_protein = getChannelAtoms(channels[1], atoms)
len(atoms_protein.getResnames())

distance = 4
residues = atoms_protein.select('same residue as exwithin '+str(distance)+' of resname FIL')
residues.select('name CA').getResnames(), residues.select('name CA').getResnums()
list(zip(residues.select('name CA').getResnames(), residues.select('name CA').getResnums()))

2. Multiple frames:
   p2 = parsePDB('LXmulti.pdb')
   channels2, surfaces2 = calcChannelsMultipleFrames(p2)

To display particular channel:
channels2[4][1]

model2 = getVmdModel(vmd_path, p2)
showChannels(channels2[4][1], model2)
showChannels(channels2[4][2], model2)

3. Trajectories:
   PDBfile = 'L_traj.pdb'
   DCDfile = 'L_traj.dcd'
   atoms = parsePDB(PDBfile)
   dcd = Trajectory(DCDfile)
   dcd.link(atoms)
   dcd.setCoords(atoms)

channels3, surfaces3 = calcChannelsMultipleFrames(atoms, dcd)
all_channels, all_surfaces=calcChannelsMultipleFrames(atoms, dcd, output_path = 'channels_dcd', separate=True)

[karolamik13]

Now option 'separate' will return (1) separate channels as PDBs and (2) all channels in one PDB file. When I was analyzing systems, I wanted to keep both for further analysis and I had to compute everything twice. Now we can generate all the useful data at once.

[jamesmkrieger]

I updated recent changes in ProDy to be sure that everything is working fine. Additionally, I added option to save the results as 'pqr' file. channels, surface = calcChannels(atoms, output_path='channels_pLoxA_ALL.pqr') If name is not with '.pdb' or '.pqr', it will be save as 'pdb'.

I suggest having the default as pqr

[jamesmkrieger]

Otherwise, this all sounds good. Hopefully I can get a chance to try it in a few weeks

[karolamik13]

I updated recent changes in ProDy to be sure that everything is working fine. Additionally, I added option to save the results as 'pqr' file. channels, surface = calcChannels(atoms, output_path='channels_pLoxA_ALL.pqr') If name is not with '.pdb' or '.pqr', it will be save as 'pdb'.

I suggest having the default as pqr

Ok. I can do that.

[karolamik13]

All the channels will be saved as PQR file. Now we will have:
(1) calcOverlappingSurfaces() - to analyze all PQR files with channels in the current directory
(2) calcOverlappingSurfaces(pqr_files='./test_1') - to analyze only test_1 folder
(3) calcOverlappingSurfaces(pqr_files=['All_channels_for_2sblA.pqr', 'All_channels_for_4wfoA.pqr', 'All_channels_for_2iukB.pqr']) - to analyze the list of selected PQR files with channels generated by CaviTracer.

PQR and PDB files have the same format but when we upload PQR file into VMD it will give us better visualization with radius per predicted bead of the channel. If someone wants PDB just change the name "All_channels_for_2sblA.pqr" to "All_channels_for_2sblA.pdb" and it will work, but then I would recommend to use QuickSurf in VMD and radius 0.5 or 0.4 to see the prediction correctly.

[karolamik13]

After changes getChannelParameters() can analyze automatically channels (single PDB) and channels2 (mult-model PDB) or channels3 (trajectory results, see in the previous comments)

[karolamik13]

One important thing left: getChannelResidueNames() that will analyze DCD and multi-model PDB. Currently it is working with single PDB file (so you can write a loop and use it anyway for DCD and PDBEnsemble), but I will get there in the upcoming days and provide an automatic function.

[jamesmkrieger]

This all sounds very good.

I am just thinking that calcOverlappingSurfaces is too general a name and sounds related to what we have in insty for hydrophobics. Maybe it can be called calcChannelSurfaceOverlaps or something like that.

[karolamik13]

Good suggestion James. Thank you. I replaced it.
Now we have:
All the channels will be saved as PQR file. Now we will have:
(1) calcChannelSurfaceOverlaps() - to analyze all PQR files with channels in the current directory
(2) calcChannelSurfaceOverlaps(pqr_files='./test_1') - to analyze only test_1 folder
(3) calcChannelSurfaceOverlaps(pqr_files=['All_channels_for_2sblA.pqr', 'All_channels_for_4wfoA.pqr', 'All_channels_for_2iukB.pqr']) - to analyze the list of selected PQR files with channels generated by CaviTracer.

[AnthonyBogetti]

@jamesmkrieger did you get a chance to test this?

[jamesmkrieger]

@jamesmkrieger did you get a chance to test this?

I don't think so really

[karolamik13]

At the current stage, a few more changes are needed. CaviTracer can be omitted in the upcoming ProDy release.

[jamesmkrieger]

Someone still needs to go through and check this in detail, which I haven't managed to do yet

At any rate, it would be good to have unit tests please

[karolamik13]

Hi James,
I plan to do it in this week. I will prepare tutorial in Jupyter Notebook file for other people who would be checking/testing CaviTracer functions.

[jamesmkrieger]

Sounds good!