Yaps: Python Frontend to Stan

TL;DR

Yaps is a new Python frontend for Stan that allows users to write Stan models using Python syntax, improving usability and integration while maintaining Stan's semantics.

Contribution

This paper introduces Yaps, a novel reinterpretation of Python syntax for Stan, enabling seamless Python-style modeling in probabilistic programming.

Findings

Successfully tested on over a thousand Stan models

Open-source implementation available for community use

Improves user experience by integrating Python syntax with Stan semantics

Abstract

Stan is a popular probabilistic programming language with a self-contained syntax and semantics that is close to graphical models. Unfortunately, existing embeddings of Stan in Python use multi-line strings. That approach forces users to switch between two different language styles, with no support for syntax highlighting or simple error reporting within the Stan code. This paper tackles the question of whether Stan could use Python syntax while retaining its self-contained semantics. The answer is yes, that can be accomplished by reinterpreting the Python syntax. This paper introduces Yaps, a new frontend to Stan based on reinterpreted Python. We tested Yaps on over a thousand Stan models and made it available open-source.

Full text

††conference: ; ††journalyear: ;††doi: ††copyright: none

1. Introduction

A probabilistic model is a mathematical model for explaining real-world observations as being generated from latent distributions (Ghahramani, 2015). Probabilistic models can be used for machine learning, and compared to alternative approaches, have the potential to make uncertainty more overt, require less labeled training data, generate synthetic data, and improve interpretability (Baudart et al., 2018). The key abstractions for writing probabilistic models are sampling of latent variables and observations and inference of latent variables (Gordon et al., 2014). Probabilistic modeling is supported by several stand-alone domain-specific programming languages, e.g., Stan (Carpenter et al., 2017). On the other hand, machine learning is supported by many Python-based packages. To capitalize on Python’s packages and familiarity, PyStan (Carpenter et al., 2017), PyCmdStan (Woodman, 2017), and other efforts (Uber, 2017; Salvatier et al., 2015; Tran et al., 2017) embed probabilistic abstractions into Python.

In programming, a watertight abstraction provides a basis for coding or debugging without leaking information about lower-level abstractions that it builds upon. Unfortunately, probabilistic abstractions in Python offered by existing efforts (Uber, 2017; Salvatier et al., 2015; Tran et al., 2017) are not watertight (Baudart et al., 2018). To code with those packages, one must also use lower-level packages such as NumPy, PyTorch, or TensorFlow. Furthermore, bugs such as tensor dimension mismatches often manifest at those lower levels and cannot be reasoned about at the probabilistic level alone. This paper introduces Yaps (Yet Another Pythonic Stan), a watertight embedding of Stan into Python. Yaps is available as open-source code (https://github.com/ibm/yaps).

Figure 1 shows a Yaps example in a Jupyter notebook (Pérez, 2014). Cell 1 imports Yaps and other Python packages and makes PyCmdStan less verbose. Cell 2 Line 1 uses the @yaps.model decorator to indicate that the following function, while being syntactically Python, should be semantically reinterpreted as Stan. Since the code is reinterpreted, its original Python interpretation is no longer available, and thus, does not leak abstractions. Line 2:2 declares coin as a probabilistic model with one observed variable x, representing ten coin tosses, each of which is either tails (0) or heads (1). The type int(lower=0, upper=1)[10] comes from Stan, where it is used for informative error messages and compiler optimizations. Line 2:3 declares a latent variable theta for the unknown bias of the coin. The initialization <****~ uniform(0, 1) samples theta with the prior belief that any bias is equally likely. Lines 2:4–2:5 indicate that each of the coin tosses x[i] is sampled from a Bernoulli distribution with the same latent theta parameter.

While watertight abstractions are essential for modeling, interaction with the host language and other libraries is essential for probabilistic inference and for exploring the results of inference. Cell 3 Line 1 represents concrete observed coin flips using the popular NumPy package (van der Walt et al., 2011). Line 3:2 passes the observed data as an actual argument to the model. And Line 3:3 calls inference to infer a joint posterior distribution using the sample function of PyCmdStan, with a random seed for reproducibility (Woodman, 2017). Cell 4 shows how to interactively explore the results of inference. Yaps creates a posterior object with fields for latent variables such as theta. Cell 4 Line 1 retrieves the posterior distribution of theta. Line 4:2 uses the popular MatPlotLib package (Hunter, 2007) to plot its histogram, and Line 4:3 prints its mean. The inferred posterior belief is that the coin is biased towards tails, with a mean theta of 0.254. That makes sense, since the observed data from Cell 3 Line 1 contained more tails than heads.

To reinterpret Python syntax, we must first parse it, and once parsed, we can easily visualize its dependencies. Cell 5 yields a visual rendering of the graphical model. Yaps compiles to Stan code, and Cell 6 prints the result, showing Stan’s explicit code blocks for observations (data), latent variables (parameters), and the actual model. In Python, we opted for a more concise (but equally watertight) syntax with implicit blocks inspired by SlicStan (Gorinova et al., 2018).

This paper argues that even when a domain-specific language such as Stan is embedded into a host language such as Python, it should avoid leaky abstractions and gratuitous use of host-language strings. This paper shows how to accomplish these objectives via reinterpreted Python, while also providing high-quality error messages. But the main contribution is a practical one: this paper introduces Yaps, a new Python-embedded frontend to Stan.

2. Related Work

PyStan (Carpenter et al., 2017) and PyCmdStan (Woodman, 2017), the existing Python interfaces for Stan, focus on communication with the inference engine. They read Stan models from a file or from a multi-line Python string and send them to the Stan compiler. Building on this bridge, Yaps adds an additional layer, supporting Python syntax for writing models. Figure 2 shows the PyCmdStan version of Cells 1–4 of the Yaps example from Figure 1 (Cells 5–6 illustrate Yaps features not present in PyCmdStan). The Yaps version is shorter than the PyCmdStan version. Both versions produce the same results because the actual Stan model in Figure 2 Cell 2 is equivalent to the model generated by the Yaps compiler in Figure 1 Cell 6.

There are various other Python-embedded probabilistic programming languages, such as PyMC3 (Salvatier et al., 2015), Edward (Tran et al., 2017), and Pyro (Uber, 2017). Compared to these, Yaps is more watertight. Those efforts use lazy evaluation: they overload operators that appear to do eager computation to instead generate a computational graph to be evaluated during inference. Lazy evaluation is a popular approach for embedding domain-specific languages into host languages (Hudak, 1998; Rompf and Odersky, 2012). Unfortunately, lazy evaluation in Python does not track local variable names, does not cleanly isolate the embedded language, and leads to verbose syntax that is less similar to stand-alone probabilistic languages. In contrast, the reinterpretation approach of Yaps avoids those disadvantages. Both approaches (lazy evaluation and reinterpretation) have the advantage of working in pure Python without a separate preprocessor.

Like our work, several other recent efforts also recognize the potential to reinterpret Python code by parsing it. Both Tangent and Myia build a computational graph from which they derive derivatives that are crucial for gradient-descent based machine learning (Wiltschko et al., 2017; Breuleux and van Merriënboer, 2017). And Relay builds a computational graph both for automatic differentiation and for mapping to heterogeneous hardware (Roesch et al., 2018). In contrast to these papers, our paper focuses on a watertight embedding of a probabilistic programming language in Python.

3. Design and Implementation

Language design and rationale.

The design of Yaps follows the motto “Stan-like for probabilistic features, Python-like for everything else”. Yaps uses familiar Python syntax for non-probabilistic features such as for loops, type declarations, or function declarations. Stan-specific features are expressed with syntactically-valid Python syntax that resembles the original syntax: <****~ for sampling or $x\!$.T[a,b] for truncated distribution. Yaps uses the syntax <****~ for sampling because it resembles ~ from stand-alone probabilistic languages such as Stan, is concise, and is syntactically valid Python (infix operator < followed by prefix operator ~).

Since observed variables are free during modeling but bound during inference, Yaps make them formal arguments of the model (Figure 1 Cell 2 Line 2) and actual arguments when preparing the model for the inference call (Figure 1 Cell 3 Line 2). Here, reinterpretation enabled more intuitive syntax. Whereas Stan requires users to place functions, parameters, transformed parameters, model, transformed data, and generated quantities into separate code blocks, Yaps uses SlicStan-style program analysis (Gorinova et al., 2018) to let users place these at the top-level in the function. This syntax is more concise and flexible and SlicStan shows it can even improve modularity. For users who prefer explicit blocks, Yaps offers a syntax based on Python with statements.

Implementation difficulties and solutions.

The top portion of Figure 3 shows the stages and representations of the Yaps compiler. Writing our own embedded Python parser would have been cumbersome, but fortunately, the Python standard library modules inspect and ast solved that for us. Our @yaps.model decorator uses those modules to replace the Python function by an intermediate representation (IR) suitable for visualization, compilation to Stan, and inference. One difficulty was to implement a suitable syntax for sampling. We first tried =****~, where = is Python’s assignment, but that was insufficient, because the left-hand side of Stan’s sampling is more expressive than that of Python’s assignment (Burroni et al., 2018). Hence we settled on <****~ instead, which we substitute with is before parsing, since that does not occur in Stan and has a low precedence in Python.

Another difficulty was name resolution. Identifiers in Yaps code refer to Stan types, functions, and distributions; for watertight abstractions, they should not be resolved to Python entities. The Python interpreter does not attempt to resolve names in Yaps function bodies. The Python interpreter does report errors for unknown names in function signatures, e.g., def model2(N: int, y: real[N]). To avoid that, Yaps provides stubs (e.g., int, real) and enables users to declare other required identifiers (e.g. N = yaps.dependent_type_var()).

Some tokens that are identifiers in Stan are keywords in Python, such as lambda. Yaps users cannot use them as identifiers in their models.

Finally, building a robust interface between Python and the Stan compiler would also have been cumbersome, but fortunately, PyCmdStan (Woodman, 2017) solved that for us. The only difficulty was that PyCmdStan error messages refer to locations in generated Stan code. We implemented a reverse source location mapping and used that to make error messages refer to source locations in Yaps code instead.

4. Evaluation

We conducted three kinds of tests to evaluate Yaps. First, we tested our compiler on a large number of programs; second, we performed runtime tests on Yaps code; and third, we tested whether Yaps yields good error messages.

To obtain many realistic programs to test our compiler on, we built a second compiler that goes from Stan to Yaps, shown at the bottom of Figure 3. This stan2yaps compiler can be used to import existing Stan models to Yaps to ease the transition. In this paper, we used the round-trip of the two compilers to enable the following experiments:

• •

We extracted 61 examples from the official Stan manual. The round-trip test succeeded for all of them (100%).

• •

We tried the round-trip test on 721 programs from the Stan dev repository. It succeeded for 700 of them ($97\%$). The failed tests used deprecated syntax.

• •

We tried the round-trip test on 500 programs from the Stan examples repository. It succeeded for 411 of them ($82\%$). Again, the failed tests used deprecated syntax.

To test whether the output of Yaps code matches the output of equivalent Stan code, we performed runtime tests as follows: we picked 13 Stan models with corresponding datasets and used our round-trip setup as above to generate Stan code through Yaps. We then compared the output of the sample function of PyCmdStan (Woodman, 2017) on the original Stan model code and the generated Stan model code. The output matched for all of the 13 models.

To test whether Yaps improves error messages, we conducted two experiments. First, we intentionally misspelled a variable name and looked at the buggy program with popular Python tools. Figure 4 shows that the VSCode editor puts a red squiggly line under variable thetap and PyLint explains that it is undefined. Similarly, a green squiggly line warns that variable theta is unused. If the code had used a string instead of a variable, VSCode and PyLint would not have been able to detect and explain these mistakes.

For the second error-message experiment, we again injected an undefined variable, then used Yaps and Stan to compile the code, see Figure 5. Stan finds the error in the compiled code, and Yaps intercepts the error message from the Stan compiler, mapping it back to the Python code. As is common in Jupyter notebooks, the error message gets displayed as part of the cell output.

5. Conclusion

This paper introduces Yaps, a new Python embedding of the Stan probabilistic programming language that embraces Python syntax. Yaps interfaces seamlessly with popular Python packages such as NumPy and MatPlotLib, while at the same time remaining well-separated from them. For testing purposes, we have done round-trip translations of over a thousand models from Stan to Yaps and back. Future work includes adding local variable type inference to make Yaps models leaner without compromising their static typing guarantees; adding meta-programming features so users can experiment with different variants of their models; and possibly adding an interactive debugger. Yaps is available at https://github.com/ibm/yaps.

Acknowledgements

We thank Javier Burroni for insightful discussions and the Stan team for their encouragement.