PyMOOSE: Interoperable Scripting in Python for MOOSE

Python is emerging as a common scripting language for simulators. This opens up many possibilities for interoperability in the form of analysis, interfaces, and communications between simulators. We report the integration of Python scripting with the Multi-scale Object Oriented Simulation Environment (MOOSE). MOOSE is a general-purpose simulation system for compartmental neuronal models and for models of signaling pathways based on chemical kinetics. We show how the Python-scripting version of MOOSE, PyMOOSE, combines the power of a compiled simulator with the versatility and ease of use of Python. We illustrate this by using Python numerical libraries to analyze MOOSE output online, and by developing a GUI in Python/Qt for a MOOSE simulation. Finally, we build and run a composite neuronal/signaling model that uses both the NEURON and MOOSE numerical engines, and Python as a bridge between the two. Thus PyMOOSE has a high degree of interoperability with analysis routines, with graphical toolkits, and with other simulators.

included optimized custom code that would allow the simulation to be run in affordable time and memory. This process of building domainspecifi c general simulators has continued with several simulators devoted to different aspects of computational and systems biology (e.g., VCell, Smoldyn, COPASI). This proliferation of simulators brings back the problems of model exchange and interoperability, albeit at a higher-level than raw Fortran or C code. While these simulators now have a common set of shared higher-level concepts (e.g., compartments, channels, synapses), they use entirely different vocabularies and languages for set up and control.
MOOSE is a new simulator project that supports simulations across a wide range of scales in computational biology, including computational neuroscience and systems biology. In order to improve interoperability, MOOSE uses two existing languages: the GENESIS scripting language, and Python. The Neurospaces (Cornelis and De Schutter, 2003; http://neurospaces.sourceforge. net/) project takes a distinct approach to supporting some GENESIS capabilities using backward-compatible scripting, and it too can utilize Python.
Most established simulators have their own scripting languages. For example, NEURON uses hoc along with modl fi les to set up simulations. GENESIS has its own custom scripting language. MOOSE avoids introducing a new language, and instead inherits the GENESIS parser. To increase compatibility, MOOSE has equivalents for most objects in GENESIS, and many old scripts can be run on MOOSE with little or no modifi cation. Given these existing capabilities, why add Python scripting? Despite its fl exibility, the GENESIS scripting language has several limitations:

INTRODUCTION
In computational biology there are two approaches to developing a simulation. First, write your custom program to do a specifi c simulation, and second, write a model and run it in a generalpurpose simulator. While the fi rst approach is very common, it requires the scientist to be a good programmer (or have one at her/his disposal) and moves the focus towards programming rather than science. Furthermore, it is very diffi cult for others to read such a program and understand how it relates to the targeted biological system. In this context, a model is a well-defi ned set of equations and parameters that is meant to represent and predict the behavior of a biological system. Ideally, a general-purpose simulator allows the model to be separated from the low-level data-structures and control. The scientist is no longer concerned with minutiae of software engineering and can concentrate on the biological system of interest. The model can be shared by other people and understood relatively easily using intermediate-level descriptions of the model with a more obvious mapping to the real biological system. General simulators also lend themselves to declarative, high-level model descriptions that have now become important part of scientifi c interchange in the computational neuroscience and systems biology communities (Beeman and Bower, 2004;Cannon et al., 2007;Hucka et al., 2002; http://www.morphml.org/; http://neuroml.org, http://sbml.org). The goal of this paper is to show how the simulator Multi-scale Object Oriented Simulation Environment (MOOSE; http://moose.ncbs.res.in/, mirrored at http://moose.sourceforge.net/) uses Python to address these issues of interoperability with analysis software, graphical interfaces, and other simulators.
General-purpose simulators have been in use since the venerable circuit simulator SPICE was utilized to solve compartmental models Segev et al., 1985). While this level of generality ran into limitations of computing power, more specialized neuronal simulators such as GENESIS and NEURON (Bower and interpreted functions are much slower than compiled code. The GENESIS scripting language itself provides for some degree of extensibility, but this is diffi cult to implement. Adding a single command requires implementation in C, as well as definition of the command in a confi guration fi le that must be pre-processed to include into the interpreter. The addition of a new class is still more involved. 3. Lack of existing libraries: The GENESIS scripting language is a special-purpose language and has no additional features other than those written into the language. 4. Syntax: The syntax is complex and inconsistent as a result of accretion of features by many developers and users. For example, arrays are implemented in three inconsistent ways in the GENESIS scripting language: as arrays of elements, entries within tables and extended fi elds.
To harness the capabilities provided by a modern widely used scripting language, we chose a Python interface. Among the plethora of programming languages, Python has some special advantages:

MATERIALS AND METHODS
There are two common approaches to create a Python interface to a C/C++ library: (1) statically link it with the Python interpreterwhich involves compiling the Python interpreter source-code, (2) create a dynamic link library and provide it as a Python module. We took the second approach as it provides more fl exibility on the choice of the Python interpreter and reduces the burden on the maintainer.

MAPPING MOOSE CLASSES INTO PYTHON
MOOSE has a set of built-in classes for representing simulation entities. These classes provide a mapping from the concept space to the computational space. Physical or chemical properties and other relevant parameters are accessible as member fi elds of the classes and the time-evolution of these parameters is calculated by a special process method of each class. These classes add another layer over ordinary C++ classes to provide messaging and scheduling as well as customized access to the member fi elds. MOOSE provides introspection (Maes, 1987;Smith, 1982), so that full fi eld information for each class is accessible to the programmer. This class information is statically initialized for each class at startup time. We utilized this class information and SWIG (Beazley, 1996; http://www.swig.org) to build the Python interface.
SWIG is a mature software with good support for Python and C/C++ interfacing as well as many other languages. While it is rather simple to create an interface for ordinary C++ class using SWIG, our task was complicated because MOOSE classes have another layer over ordinary C++ classes. For this reason we created a framework for Python interface with additional C++ classes to wrap MOOSE classes and a few classes to manage the system.

SIMULATOR CONTROL THROUGH PYTHON
All operations on MOOSE objects are carried out via a special class, Shell, of which there is a single instance on each processor node that is running MOOSE. In PyMOOSE we implemented a singleton context object to communicate with the Shell. The context object provides a set of functions that can be called to pass appropriate messages to the Shell. The user can call global MOOSE functions by calling the corresponding methods of the context object. Operations like creation of objects, setting integration time step, running the simulation are all done through the context object.
We created a one-to-one mapping of MOOSE classes to Python classes by means of light-weight C++ wrapper classes. All the wrapper classes were derived from one common base class. Each MOOSE object is identifi ed by an Identifi er (ID) fi eld. The main data content of a wrapper class instance is the ID of the corresponding object in MOOSE. Additionally, the wrapper classes have a static pointer to the single instance of the context object. Wrapper classes provide accessor methods that can be used to access the fi elds in the corresponding MOOSE object.
These C++ wrapper classes were input to SWIG to create the Python module. After translation to Python, the user sees the member fi elds in the Python classes in place of the accessor methods in the C++ wrapper classes. Behind the scene the Python interpreter calls these accessor methods whenever the user script tries to access MOOSE object fi elds ( Figure 1A).
Manually developing C++ wrapper classes for all MOOSE classes was a tedious but repetitive task. We therefore embedded stub code in the MOOSE initialization code to generate most of the wrapper code programmatically using Run-Time Type Information in C++. This auto-generated code was used with a few modifi cations to generate a Python module using SWIG. SWIG takes an interface fi le with SWIG-specifi c directives and generates a single C++ fi le for the library and a Python source-code fi le that contains support code. We completed the PyMOOSE code generation by compiling and linking the SWIG-generated C++ source-code as a dynamic library. This dynamic library can be imported in any Python program.

LEGACY MODELS AND PyMOOSE
The PyMOOSE context object keeps a single instance of the GenesisParser class in order to run legacy GENESIS scripts. Whenever the user asks for executing a GENESIS statement, the context object disconnects itself from the Shell and connects the GenesisParser object instead. The GENESIS statement string is passed to the GenesisParser object, which executes it as if the user typed it in at the MOOSE command prompt. After execution of the statement (or script) the GenesisParser object is disconnected from the Shell and the context object is reconnected ( Figure 1B).
While it is valuable to run GENESIS scripts within PyMOOSE, this feature is intended only to support legacy code and is better avoided in new model development. The use of GENESIS scripting language inside Python defeats the whole purpose of moving to a general-purpose programming language. It reduces readability and the user needs to know both languages in order to understand the code.

RESULTS
We used the Python interface of MOOSE to achieve three key targets: (1) Interfacing with standard libraries in a mature scientifi c computing language, (2) giving access to a portable GUI library for developing user interface and (3) enabling MOOSE to work together with other simulators.

INTERFACING SIMULATIONS WITH PYTHON LIBRARIES
We used Python scientifi c and graphing libraries to analyze and display the output of a PyMOOSE simulation. The interface with Python gives the user freedom to choose from a wide variety of scientifi c and numerical libraries available from third parties. We demonstrate the use of two libraries along with PyMOOSE for developing simulations with plotting and data analysis within Python. The fi rst of these, NumPy, is a library that provides data structures and algorithms for fast matrix manipulation (http:// numpy.scipy.org/). Even though Python is interpreted, with attendant slow execution, NumPy library provides access to compiled code and hence the functions from the library are as fast as compiled code. The second library, matplotlib, provides a rich set of functions for plotting 2D data both in hardcopy formats and interactively (http://matplotlib.sourceforge.net/). It can use NumPy for fast matrix operations in Python and several portable GUI toolkits (GTK/Qt/Tk/wxWidgets) as graphical back-end.
We implemented a simulation of the squid giant axon using Hodgkin-Huxley Na + and K + channels and parameters (script attached in Appendix). We applied an injection current with random amplitude uniformly distributed between 0 and 100 nA. We recorded the time-series for the membrane potential during the simulation in a MOOSE table object, which can accumulate a time-series of simulation output (Figure 2A). The interface to Python was done using the MOOSE table class. This class is exposed to Python with methods to emulate iterable type (Martelli et al., 2005). The array constructor in NumPy accepts an iterable object and creates a NumPy array with a copy of the contents of the object. Thus the user is relieved of explicitly iterating over the table entries and copying them to a NumPy array. This completes the interface from the MOOSE simulation output to NumPy ( Figure 2B). We used the fast Fourier transform operation available in NumPy to compute the discrete Fourier transform of the time-series of the simulated membrane potential. We used matplotlib to plot the original time-series, as well as the output of the FFT ( Figure 2C).
Overall, this example simulation illustrates how PyMOOSE facilitates interoperability of Python numerical and graphing libraries with MOOSE.

PORTABLE GUI THROUGH PYTHON
The use of Python separates the problem of GUI development from simulator development. Moreover, it gives one the freedom to choose from a number of free GUI toolkits. The major platform independent GUI toolkits with Python interfaces are Qt(TM) available as PyQt, wxWidgets (wxPython), Tk and GTK (http://wiki. python.org/moin/GuiProgramming; http://www.python.org/doc/ faq/gui/). We used PyQt4 to develop a simple user interface for a clone of the GENESIS squid tutorial in MOOSE. We selected Qt4 as it is a mature and clean toolkit that is freely distributed and runs well on all the major operating systems.
The program was divided into three modules -(1) the squid axon compartment with Hodgkin-Huxley channels, (2) a model object which combined a few tables with the squid compartment to record various parameters through the time of the simulation, and The Shell object is usually controlled through the PyMooseContext. When loading a GENESIS script, control is temporarily passed to the legacy GENESIS script language parser, and then returned to the PyMooseContext.
(3) the GUI to take user inputs and to plot data. We implemented the squid axon model as described in the previous section, using PyMOOSE to set up and parameterize the model. As before, the model was interfaced with table objects to monitor time-series output of the simulation. Finally, we implemented the GUI by loading in the PyQt4 libraries, and using Python calls to set up the interface (Figure 3). While there are Qt IDEs available (http://trolltech. com/products/qt/), we constructed the interface through explicit Python calls to create widgets, assign actions, and manage output data. Qt uses a signal-slot mechanism for passing event information. PyQt allows the use of arbitrary Python methods to be used as slots. Hence we could connect the GUI widgets to methods in the PyMOOSE model class and thus provided simulation control through the GUI in a clean manner. We used PyQwt, a Python interface of the Qt-based plotting library Qwt, for creating output graphs. Since PyQwt can take NumPy arrays as data, we converted the tables in MOOSE to NumPy arrays and used PyQwt plotting widgets to display them.
We based the layout of the simulation on the widely used GENESIS Squid tutorial program. To confi rm portability of the system, we ran the model on Linux as well as the Windows operating system. This exercise demonstrated the capability of PyMOOSE to draw upon existing graphical libraries for its graphical requirements. This is an important departure from GENESIS. The GENESIS graphical libraries (XODUS) were an integral part of the C code-base and XODUS objects were visible as, and manipulated in the same way as other GENESIS objects. In contrast, PyMOOSE did not need to implement any graphical objects within the MOOSE C++ code, but instead reused extant third-party graphical libraries available for Python. Furthermore the existing libraries are professionally designed and have a much more consistent look-and-feel than did the original GENESIS graphical library, XODUS (Bhalla, 1998).

SIMULATOR INTEROPERABILITY
With Python becoming a popular language for developing platform independent scripts, several neuronal simulators have implemented Python interfaces. This raises the possibility of using Python as a glue language to run simulations that span different simulators. As a fi nal demonstration of interoperability, we used PyMOOSE with PyNEURON to build a multi-scale, multi-simulator model that incorporates neuronal electrical activity as well as biochemical signaling ( Figure 4A).
We used NEURON to model a multicompartmental electrical model of a Type A neuron from the CA3b region of the rat hippocampus (Migliore et al., 1995; http://senselab.med.yale. edu/ModelDB/ShowModel.asp?model=3263). This is a morphologically detailed model with experimentally constrained distribution of membrane ion channels. It reproduces experimental observations of fi ring behavior and intracellular Ca 2+ dynamics. We modifi ed the hoc script for the model, to run it for arbitrary time intervals. We directed the output data to Vector objects in NEURON. The Python wrapper class for this model provided a handle for the simulation parameters and functions defi ned in the hoc script. As described in the PyNEURON documentation (http://www.neuron.yale.edu/neuron/docs/help/neuron/ neuron/classes/python.html), Python commands were directed to the NEURON engine by constructing hoc statement strings and executing them through the hoc interpreter instance provided by the neuron module. Moreover, hoc object references are directly available in Python as attributes of the hoc interpreter object. Thus accessing hoc objects was quite clean in Python ( Figure 4A).
We used MOOSE to model calcium-triggered biochemical signaling events at the synapse. We used a model of a bistable MAPK-PKC-PLA2 feedback loop that was originally implemented in GENESIS/Kinetikit (Ajay and Bhalla, 2004;Bhalla and Iyengar, 1999;Bhalla et al., 2002) and uploaded to the DOQCS database (http://doqcs.ncbs.res.in/template.php?&y=accessiondetails&an= 79). The model was defi ned in the GENESIS scripting language. We used the legacy scripting mode of PyMOOSE to load the GENESIS/ kinetikit model. The simulation objects thus instantiated were standard MOOSE objects, and were accessible using Unix-like path strings. The PyMOOSE interface exposed these objects as regular Python objects. Thus access to the MOOSE objects, representing GENESIS data concepts, was also straightforward in Python ( Figure 4A).
We used the Python interface to accomplish three critical operations to combine the two simulations: (1) Initialization, (2) runtime control and synchronization, and (3) variable communication and rescaling.
1. To initialize the models, we used PyNEURON command load_ fi le to load the hoc script. Once the script is loaded, variables and functions defi ned in the script become available as members of the hoc interpreter instance inside Python. In this case we defi ned a setup function to initialize the NEURON simulation. This function is called in the constructor (__init__) of the Python wrapper class over the NEURON simulation. At this  We wrote another higher-level function run to advance the coupled simulations using the two wrapper classes (not to be confused with the member method run of these classes). This function (1) creates instances of both wrappers, which involves initializing the models, (2) runs the NEURON simulation for 1 s, (3) reads out the calcium level, performs rescaling and updates the kinetic model with this value, (4) advances the kinetic simulation for 1 s to catch up with the electrical model, (5) reads out the activity level of MAPK from the GENESIS/Kinetikit model and modifi es the [Ca 2+ ] dependent K + channel conductances in the NEURON model in inverse proportion to this (Figures 4E,F).
Our simulated experiment is illustrated in Figure 4F. We loaded the models and allowed them to settle. We measured baseline neuronal responses at this stage using a 250-ms, 0.15 nA current pulse. Following this we used the run function for the further timeevolution of the system. We applied a strong LTP-inducing stimulus to the neuronal model for 7 s, and then allowed the simulation to continue for 183 s. Finally we repeated the 250 ms, 0.15 nA test for neuronal responses.
The time-evolution of membrane potential, Ca 2+ levels, and MAPK activity are shown in Figures 4D,E. The initial and fi nal burst waveforms of the neuron are shown in Figures 4C,G. We observe that the coupled model shows how electrical stimulation can lead to signaling events, with feedback effects on the electrical properties of the neuron. We should point out that this simulation is only a demonstration and the relationship between the chemical system and the biophysical properties of the neuron is over-simplifi ed, although the two component models we used are realistic within their respective domains.
This example also illustrates the effi ciency of using Python for data transfer when traffi c volumes are small compared to the computational times. The neuronal calculations in NEURON took about 91% of the simulation run-time, the signaling calculations in MOOSE took ∼8.5%, and the data transfer through Python accounted for only around 0.5%. As we discuss below, there may be other interface contexts where more effi cient, low-level data transfer protocols may be needed, and the relatively facile Python interface may not be appropriate.

DISCUSSION
We have used PyMOOSE, the Python interface to MOOSE, to achieve interoperability at three levels. First, we used standard mathematical packages in Python to analyze MOOSE output. Second, we used the QT graphical toolkit from within Python to build a GUI for a MOOSE simulation. Third, we used Python as a glue language to run a cross-simulator model combining an electrophysiological model set up in NEURON with a biochemical signaling model set up in GENESIS/Kinetikit.

ISSUES WITH PYTHON INTEROPERABILITY
The strengths of the Python language make it perhaps too easy to repeat well-known mistakes in simulation development. We consider two such issues. First, Python is an interpreted language in most implementations. In the context of simulations, it is not meant for number crunching. Well-designed libraries like NumPy can hide some of these limitations from the user, and fast hardware can conceal other ineffi ciencies. However, given the same specialized algorithms, a compiled language will perform better than an interpreted one. Therefore, for large simulations, we need to combine the best possible algorithms with optimized and compiled languages. MOOSE has as one of its goals the capability of managing the low-level, high-traffi c fl ow of data between different numerical engines incorporated into MOOSE. We do not consider Python appropriate for such operations. Second, many aspects of model specifi cation should be done using declarative rather than procedural approaches (Cannon et al., 2007;Crook et al., 2005Crook et al., , 2007. However, Python makes procedural model defi nition very easy, and may even provide a certain level of interoperability if several simulators provide equivalent calls for model setup. For example, there are some impressive recent efforts to develop a standard vocabulary for network defi nitions across simulators (http://neuralensemble. org/trac/PyNN/; this issue). While the presence of Python as a common link language may temporarily address the interoperability issues of this approach, we feel that it would be a cleaner design to use a separate, declarative defi nition for networks such as NeuroML (http://neuroml.org). Nevertheless, we completely agree that a standard vocabulary for model defi nitions is an important fi rst step toward this goal.

MODEL SPECIFICATION VS. SIMULATOR CONTROL
Model specifi cation and exchange issues have been ably addressed by the communities developing model specifi cation languages (Le Novère et al., 2005;Qi and Crook, 2004; http://neuroml.org; http://sbml.org). The current paper focuses on the second problem, that of making it easier for researchers to control and set up these diverse simulation tools. We have shown how this can be done with the simulator MOOSE, using Python as a glue language. Run-time communication between simulators has previously been achieved using the NEOSIM framework, which uses Java Howell et al., 2002). More recently, the MUSIC framework specifi es an API for simulators to use to communicate with each other (Ekeberg and Djurfeldt, 2008). Our study is novel in two respects. First, we use the built-in Python capabilities of two simulators to achieve run-time communication, without the need to modify either simulator or to build an additional framework for communication. Second, we carry out bidirectional communications across scales (biophysical to biochemical models) and involving continuous data types (channel conductance and calcium concentrations) rather than spike events.
The evolution of neuronal simulator technology has seen a gradual separation of different aspects of modeling, with a corresponding improvement in interoperability. The fi rst step was to develop higher-level simulation tools (e.g., NEURON and GENESIS) to separate the numerical and housekeeping code from the modelspecifi c code. This let people share models, provided they were written for the same simulator. The second was the development of declarative model specifi cations that were separate from the simulator. This initially took the form of semi-declarative cell morphology fi les (NEURON '.geom' fi les and GENESIS '.p' fi les), which required additional fi les for channel specifi cation. This process of separation of model defi nition from simulator control has continued. The Neuroconstruct suite refi nes the declarative definition of models, with NeuroML and ChannelML as declarative defi nitions suffi cient for most single-neuron models. Importantly, at this level quite different simulators can use the same original model defi nition to run simulations. A third stage is the convergence of different simulators to use the same link language, in this case Python. This makes it possible to explicitly separate model defi nition from simulator control. In the current paper, we have illustrated this with a composite signaling-neuronal model drawing on NEURON and MOOSE. We have utilized two legacy models, one written for NEURON, and one written for GENESIS. Even though the legacy models themselves were not entirely set up in a declarative manner, we used the original model defi nitions only to load in the model specifi cations. We used Python as the procedural language to control these operations, and to mediate communication between the models at run-time.

SUSTAINABILITY OF PYTHON INTEROPERABILITY
Simulator interoperability has long been regarded as important (Crook et al., 2005(Crook et al., , 2007. Such projects have been diffi cult to execute, and still harder to maintain, because they depend on multiple underlying simulator projects, each with different APIs, directions and life-cycles. Python is a potential way out of this problem. First, Python itself is a well-established language with a strong community and support. Second, the issues of interfacing to Python are now being undertaken by individual simulator development teams. Interoperability emerges from these independent efforts rather than requiring a separate project to achieve coordination. Third, PyMOOSE itself will be maintained for the long-term, since Python will be the default scripting language for MOOSE. We suggest that long-term improvements in interoperability will be driven both by widespread simulator support for declarative model specifi cations, and by a richer ecosystem of simulators fl uent in Python.