ase python chemistry html 28 May 2017
Since I first heard about the Jupyter notebook (around that time it was called ipython notebook) I instantly adopted it into my work flow for rapid testing, developing and most of all experimenting and playing around with Python code. I realized quickly how comfortable it is as a working environment and started creating notebooks using more and more of the capabilities. I began integrating plots, images, markdown notes and even creating slide presentation in jupyter.
Since I'm a computational chemist I manipulate the chemical structures a lot and the package I've been using a lot lately is the Atomic Simulation Environment (ASE). Mostly because it's written in Python and has a great assortment of methods for running and analyzing DFT calculations. It also provides a set of convenient methods for handling chemical structures including periodic ones. The package comes with it's on GUI, that does a decent job when it comes to displaying and handling structures however when working in a notebook it would be great to have a tool that allows embedding the viewer inside the notebook. Some time ago there was even a thread on the ase-users mailing list where a question about embedding ASE structures in a Jupyter notebook was raised.
I found a few ways that enable the interactive visualization of chemical
structures given as the Atoms
objects (internal representation in ASE) in
the notebook, that I think are interesting to try out:
The list above is by no means complete and there are probably some other great tools that provide similar functionality, so if you think I skipped some important alternative - let me know.
ASE has a builtin format
converter (or
writer) to html that uses the x3dom library to
create an interactive view of the molecular structure once you open or embed
the generated html in a browser. You can interact with your molecule by
rotating, translating zooming and panning the view. To embed the view in one
of the Jupyter output cells we can use the native HTML
function
from the IPython.display
module
from IPython.display import HTML
A small hurdle is that we would like to have the html as a string that can be
passed to the HTML
functions but the ASE html writer needs to write a physical
file. One of the ways of fixing this behavior is a to use a
named temporary file
from the tempfile package and
a custom function that takes the Atoms
objects and returns the html string.
def atoms_to_html(atoms):
'Return the html representation the atoms object as string'
from tempfile import NamedTemporaryFile
with NamedTemporaryFile('r+', suffix='.html') as ntf:
atoms.write(ntf.name, format='html')
ntf.seek(0)
html = ntf.read()
return html
Then to display any structure available as an Atoms
instance, we could do
from ase.build import molecule
tbut = molecule('trans-butane')
tbut_html = atoms_to_html(tbut)
HTML(tbut_html)
If everything went well you should see a view widget with a trans-butene molecule displayed, similar to the one below. If you want to try it yourself on a couple more examples see, the short notebook I used to test this code.
The most feature rich option is offered by the nglview package that provides a Jupyter widget for interactive visualizations of chemical structures and trajectories. It is built on top of the NGL Viewer and it supports some of the more popular formats through a number of convenience functions. The supported formats include:
Most importantly it also supports ase.Atoms
though nglview.showase
function, which
displays a single structure and nglview.showasetraj
that is capable of showing an
animation based on the images in the trajectory file.
The installation is pretty straightforward with conda:
$ conda install nglview -c bioconda
or pip:
$ pip install nglview
but I found I also need to enable it with
$ jupyter-nbextension enable nglview --py --sys-prefix
The feature that I really like, that is not available in other viewers, is displaying trajectories and sequences of structures, which makes it possible to visualize molecular vibrations, structure relaxation steps or reaction paths. As an example consider viewing a vibrational mode from the trajectory file:
import ase.io
mode30 = ase.io.read('vib.30.traj', index=':')
nglview.show_asetraj(mode30)
Moreover nglview let you save the movie as gif file.
If you need more control for tweaking the display you can activate the gui
mode by
passing an additional parameter to the viewer
nglview.show_asetraj(mode30, gui=True)
that will show a menu with a lot of options to adjust.
imolecule is an embeddable WebGL viewer that comes with a convenient chemical format converting features that largely rely on the capabilities on openbabel.
The embedding capabilities are nicely illustrated on the page provided by the author and the compatibility with Jupyter is demonstrated in this example notebook.
The installation is as simple as running:
pip install imolecule
The imolecule package can handle multiple file formats, and can take string
representation of structures as well a read files, however it cannot handle
the default ASE format. It uses it's own molecule representation internally that is
based on json so we can write a small
function for converting ase.Atoms
to json. The are probably a few alternative ways
of getting it done but I chose to convert ase.Atoms
to OBMol
first and then
use the converter from imolecule to do the rest. I chose that solution since
by default ase.Atoms
has no information about chemical bonds, and we can
infer that information by using openbabel. The bonding data is then used to
correctly render the stick part of the ball and stick representation of structures.
The atoms_to_json
function does the above conversion, it takes the ase.Atoms
objects and returns a json representation compatible with imolecule.
def atoms_to_json(aseatoms, infer_bonds=True):
'''
Convert ASE Atoms isntance into the json format compatible with
Args:
aseatoms : ase.Atoms
Instance of Atoms from ase package
infer_bonds : bool
If `True` bonds will be inferred using openbabel
Returns:
mol : dist
A dictionary with the json format of the molecule
'''
import pybel
ob = pybel.ob
obmol = ob.OBMol()
obmol.BeginModify()
for atom in aseatoms:
obatom = obmol.NewAtom()
obatom.SetAtomicNum(int(atom.number))
obatom.SetVector(*atom.position.tolist())
# If there is no bond data, try to infer them
if infer_bonds:
obmol.ConnectTheDots()
obmol.PerceiveBondOrders()
# Check for unit cell data
if any(aseatoms.pbc):
uc = ob.OBUnitCell()
uc.SetData(*(ob.vector3(*v) for v in aseatoms.get_cell()))
uc.SetSpaceGroup('P1')
obmol.CloneData(uc)
obmol.EndModify()
mol = pybel.Molecule(obmol)
return imolecule.format_converter.pybel_to_json(mol)
chemview is another option which utilizes the WebGL and three.js libraries. It is pretty easy to install since it can be done using conda:
conda install -c http://conda.binstar.org/gabrielelanaro
however I had to clone latest repository and install it through pip
git clone https://github.com/gabrielelanaro/chemview
cd chemview
pip install .
due to some dependency conflicts for Python 3.6.1. After installing I had to manually enable the widget with jupyter:
jupyter nbextension enable widgetsnbextension --user --py
jupyter nbextension install --user --py --symlink chemview
jupyter nbextension enable --user --py chemview
The MolecularViewer
class that is responsible for creating the visualization
accepts the (x, y, z) coordinates and additional topology
dictionary that specifies
a list of chemical symbols and chemical bonds as a list of connected atom indices.
As in the case of imolecule we can convert the ase.Atoms
instance into the required
format using the code below.
def atoms_to_chemview(atoms, infer_bonds=True):
'''
Convert ase.Atoms instance into a dict of kwargs for the chemview.MolecularViewer
Args:
atoms : ase.Atoms
infer_bonds : bool
Get the list of indices of connected atoms, (requires pybel)
Returns:
data : dict
A dict with kwargs that can be passed to chemview.MolecularViewer
'''
data = dict()
# convert the coordinates to nanometers
data['coordinates'] = atoms.get_positions() / 10.0
data['topology'] = dict()
data['topology']['atom_types'] = atoms.get_chemical_symbols()
if infer_bonds:
import pybel
ob = pybel.ob
obmol = ob.OBMol()
obmol.BeginModify()
for atom in atoms:
obatom = obmol.NewAtom()
obatom.SetAtomicNum(int(atom.number))
obatom.SetVector(*atom.position.tolist())
obmol.ConnectTheDots()
obmol.PerceiveBondOrders()
bonds = [[b.GetBeginAtom().GetIndex(), b.GetEndAtom().GetIndex()]
for b in ob.OBMolBondIter(obmol)]
data['topology']['bonds'] = bonds
return data
Here I used pybel again to get the information about chemical bonds. Notice also that I converted the coordinates to nanometers since it looks a lot nicer, although I didn't found any information in the chemview's docs that this is the default unit used.
If you want to see a working example see this notebook where the above function is used to display a couple of ASE structures.