ase python chemistry html 28 May 2017

Visualizing ASE structures in Jupyter notebook

A few options worth considering.

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.

html

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.

nglview

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

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

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.