Given the 3D structure of a molecule (e.g. as *.mol) and a probe radius r, I want to calculate the solvent excluded surface (SES) and the excluded volume enclosed by it. Note the difference between SAS and SES:

enter image description here

I'm looking for a command line program, C/C++ or Python library that can accomplish such a task.



MSMS is a rather old (1996) command line tool that computes the Solvent Excluded Surface (SES). It is wrapped in a python module in mgtool or in openstructure which may be more convenient to use.

  • 1
    $\begingroup$ Includes everything I was looking for and the Python bindings is a neat plus. You set the bar very high for potential competing answers. $\endgroup$ – Kexanone Mar 15 at 19:57
  • 1
    $\begingroup$ Nice to hear you got an answer that you like! The other answer is substantially longer, was that suggestion not useful for you? $\endgroup$ – Nike Dattani Mar 19 at 16:24
  • 1
    $\begingroup$ I personally went for MSMS, but I will give him the bounty for his effort. $\endgroup$ – Kexanone Mar 21 at 20:21


EDTSurf[1] is a common and reasonably well regarded approach. The source code is available at https://zhanglab.ccmb.med.umich.edu/EDTSurf/

I tested the mac version and it worked with the following:

./EDTSurf_mac -i kras.pdb

which outputs kras.ply. A slight hiccup is the website describes different command line arguments - but the default is the SES. You'll see they have windows and linux versions, too.

Assuming you're comfortable with python, you could then visualize the output with Open3D:

import open3d as o3d
mesh = o3d.io.read_triangle_mesh("./kras.ply")

Now, for volume. Volume of the SES is difficult to calculate exactly, but the volume of a mesh is trivial to calculate. Open3D has a get_volume method:


unfortunately that froze my jupyter notebook every time - you may have better luck with a fresh install.

But the following will also work:

import numpy as np

verts = np.array(mesh.vertices)
faces = np.array(mesh.triangles)

def signed_vol_of_triangle(p1, p2, p3):
    v321 = p3[0]*p2[1]*p1[2]
    v231 = p2[0]*p3[1]*p1[2]
    v312 = p3[0]*p1[1]*p2[2]
    v132 = p1[0]*p3[1]*p2[2]
    v213 = p2[0]*p1[1]*p3[2]
    v123 = p1[0]*p2[1]*p3[2]
    return (1 / 6)*(-v321 + v231 + v312 - v132 - v213 + v123)

def make_vol(pts):
    return signed_vol_of_triangle(pts[0], pts[1], pts[2])

print(sum([make_vol(v) for v in verts[faces]]))

In short, this loops through each triangle in the mesh and calculates the signed volume of the triangular prism with the fourth point at the origin. The theory behind this is in [2], and I used the implementation in this SO answer which you might want to also consult: https://stackoverflow.com/a/1568551/3089865

While I haven't rigorously tested that volume-calculating approach for correctness, it returns 21058.959617488006 angstrom cubed, which seems reasonable to me (~21 nm^3) given KRAS fits in a ~4x3x3nm box.

[1] Xu, Dong, Hua Li, and Yang Zhang. "Protein depth calculation and the use for improving accuracy of protein fold recognition." Journal of Computational Biology 20.10 (2013): 805-816.

[2]Zhang, Cha, and Tsuhan Chen. "Efficient feature extraction for 2D/3D objects in mesh representation." Proceedings 2001 International Conference on Image Processing (Cat. No. 01CH37205). Vol. 3. IEEE, 2001.


Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.