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
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
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
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
$ 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
that will show a menu with a lot of options to adjust.
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
OBMol first and then
use the converter from imolecule to do the rest. I chose that solution since
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.
atoms_to_json function does the above conversion, it takes the
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)
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
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.