Computational Statistics in Python (re: PyMC3, PyStan)

331 views
Skip to first unread message

Allen B. Riddell

unread,
Jun 16, 2015, 2:42:04 PM6/16/15
to stan...@googlegroups.com
Very impressive (draft) introduction to using Python in statistics which
happens to feature a side-by-side comparison of PyMC3 and PyStan:

http://people.duke.edu/~ccc14/sta-663/PyMC3.html

Bob Carpenter

unread,
Jun 16, 2015, 3:07:00 PM6/16/15
to stan...@googlegroups.com
That's awesome.

But I didn't see any commentary, just running through a few models
in both systems.

- Bob
> --
> You received this message because you are subscribed to the Google Groups "stan development mailing list" group.
> To unsubscribe from this group and stop receiving emails from it, send an email to stan-dev+u...@googlegroups.com.
> For more options, visit https://groups.google.com/d/optout.
>

Allen B. Riddell

unread,
Jun 16, 2015, 8:40:45 PM6/16/15
to stan...@googlegroups.com
On 06/16, Bob Carpenter wrote:
> That's awesome.
>
> But I didn't see any commentary, just running through a few models
> in both systems.
>

yes, I suppose side-by-side was a bit of an exaggeration. There is some
prose comparing the two:

> The nice thing about PyMC is that everything is in Python. With PyStan,
> however, you need to use a domain specific language based on C++ synteax
> to specify the model and the data, which is less flexible and more work.
> However, in exchange you get an extremely powerful HMC package (only
> does HMC) that can be used in R and Python.

Andrew Gelman

unread,
Jun 16, 2015, 9:10:35 PM6/16/15
to stan...@googlegroups.com
Is it really more work to write a statistical model in Stan than in PyMC? I doubt it.

Bob Carpenter

unread,
Jun 16, 2015, 9:48:00 PM6/16/15
to stan...@googlegroups.com
Here are some more thoughts on various topics. Bottom line
is we should stay on our toes, but I think we'll be OK for
the foreseeable future.

This isn't so much just PyMC3, of course --- there are lots of
other packages out there with functionality very much overlapping
Stan's (BUGS, JAGS, Laplace's Demon, Church and its progeny,
Figaro, emcee, NONMEM, GPstuff, lme4, NIMBLE, INLA, etc. etc.).


METAPROGRAMMING

The key drawback for Stan is that we can't do any
metaprogramming relative to the Stan object language because
the Stan program itself isn't an object.

I'm genuinely curious what people would use it for.

I think this is an area that's really worth exploring.
It really struck me looking at Jeromy Anglin's workflow
for building and comparing multiple models in JAGS. He
was reduced to cutting and pasting strings to swap in
and out priors, etc., which is pretty painful.

One thing we could do is define a graphical language within
Python and translate to Stan. Perhaps even PyMC's if it's
decent --- I don't quite understand its boundaries from
the examples. I was hoping we could piggyback on NIMBLE for
that in R.


FLEXIBILITY

Exactly what I'm trying to understand is how flexible
the PyMC modeling language is given that it's backed by
Theano. Rob was digging around in the Python source today
trying to figure out which probability functions they
support.

I have the same question about Church, etc.


HARD WORK

Learning the Stan language involves work --- a lot of work to
learn it well. I'd compare that to learning the PyMC structures
you need in Python. In the end, assuming Python competence in
users and reasonable design choices, PyMC3 should be
less work to learn.

PyMC3 is also completely in-workflow --- you don't need to call
another language and have source files for it. This is good, but
to the extent that we could define a model by a simple Python
structure, we might be able to swing something along these lines
ourselves. In which case, it'd just come down to ease of installation
and who designed the better Python interface.

From what I've seen, Stan's language is more direct for statistical
modeling, but is obviously way more limited in what it can do than
Python. We can't stop and open up the Twitter API for more data in
the middle of a run, for instance --- I'm not sure if PyMC3 would
give you that kind of control (defining a variable in the model
as doing external work or getting updated asynchronously).


SPEED

Someone should measure this. We've put a ton of effort into
our derivative speed and multivariates. It'd be interesting
to see how PyMC3 compares.


FEATURES

Also consider things like constrained variables, covariance priors,
Cholesky factors, ODE solvers, etc. These should all be possible
in something like PyMC3, but I have no idea to tell what they've
already implmented.

They certainly already have optimization, but I don't know if the
do the lm-like uncertainty estimation. Of course, that'd be pretty
easy.


PORTABILITY

Stan models are portable across R, Python, command-line, etc.
That's good from a model-building community point of view, but
sub-optimal from an infrastructure perspective and it means we
have to do a ton of work as a team to field all these interfaces (the
point of the ongoing Stan 3 design and Allen's visit was to figure out
simpler interfaces without losing any power).

The same goes for platforms. Windows is an ongoing pain, but there
are a lot of users there. R's also been restricting what we can
do, as is Python to some extent in terms of moving to C++11, updating
dependent packages, structuring interface data types, etc.

For all we know, Julia's the future and Python and R are soon-to-be
zombies! I'd rather hedge our bets on the platform choices.

- Bob

Allen B. Riddell

unread,
Jun 17, 2015, 1:35:15 PM6/17/15
to stan...@googlegroups.com
On 06/16, Bob Carpenter wrote:
> PyMC3 is also completely in-workflow --- you don't need to call
> another language and have source files for it. This is good, but
> to the extent that we could define a model by a simple Python
> structure, we might be able to swing something along these lines
> ourselves. In which case, it'd just come down to ease of installation
> and who designed the better Python interface.
>
> From what I've seen, Stan's language is more direct for statistical
> modeling, but is obviously way more limited in what it can do than
> Python. We can't stop and open up the Twitter API for more data in
> the middle of a run, for instance --- I'm not sure if PyMC3 would
> give you that kind of control (defining a variable in the model
> as doing external work or getting updated asynchronously).
>

This is certainly where someone might see the advantage to "staying in"
a language like R or Python.

My sense is that we haven't fully explored what is possible here. There
aren't many concrete limitations to how Python and R can interact with
the Stan C++ library. If we identify something we'd like to achieve that
PyMC3 or a similar library can do (e.g., stopping and fetching some data
from somewhere), I suspect there's some way to accomplish it. Both
Python and R can transparently wrap and manipulate just about anything
that one can create in C++.

Best,

Allen

Bob Carpenter

unread,
Jun 17, 2015, 4:09:54 PM6/17/15
to stan...@googlegroups.com

> On Jun 17, 2015, at 1:35 PM, Allen B. Riddell <a...@ariddell.org> wrote:
>
> On 06/16, Bob Carpenter wrote:
>> PyMC3 is also completely in-workflow --- you don't need to call
>> another language and have source files for it. This is good, but
>> to the extent that we could define a model by a simple Python
>> structure, we might be able to swing something along these lines
>> ourselves. In which case, it'd just come down to ease of installation
>> and who designed the better Python interface.
>>
>> From what I've seen, Stan's language is more direct for statistical
>> modeling, but is obviously way more limited in what it can do than
>> Python. We can't stop and open up the Twitter API for more data in
>> the middle of a run, for instance --- I'm not sure if PyMC3 would
>> give you that kind of control (defining a variable in the model
>> as doing external work or getting updated asynchronously).
>>
>
> This is certainly where someone might see the advantage to "staying in"
> a language like R or Python.
>
> My sense is that we haven't fully explored what is possible here.
...

I'm very much interested in a language for metaprogramming
models and it'd make sense to have that be native Python.

What I really want is modularity --- some way to take reasonable
programming practice to designing a model where you can test
a component as a module before plugging it into something larger.
Right now, our programs are big blobs.

But it'd be fun to brainstorm ideas in general as to what we
could do with metaprograms.

Does anyone know what people do with this capability in other
environments?

- Bob

Joshua N Pritikin

unread,
Jun 17, 2015, 8:33:28 PM6/17/15
to stan...@googlegroups.com
On Wed, Jun 17, 2015 at 04:08:26PM -0400, Bob Carpenter wrote:
> > My sense is that we haven't fully explored what is possible here.
> ...
>
> I'm very much interested in a language for metaprogramming
> models and it'd make sense to have that be native Python.
>
> What I really want is modularity --- some way to take reasonable
> programming practice to designing a model where you can test
> a component as a module before plugging it into something larger.
> Right now, our programs are big blobs.
>
> But it'd be fun to brainstorm ideas in general as to what we
> could do with metaprograms.
>
> Does anyone know what people do with this capability in other
> environments?

That's one of the primary design ideas of OpenMx,
http://openmx.psyc.virginia.edu/

--
Joshua N. Pritikin
Department of Psychology
University of Virginia
485 McCormick Rd, Gilmer Hall Room 102
Charlottesville, VA 22904
http://people.virginia.edu/~jnp3bc

Bob Carpenter

unread,
Jun 18, 2015, 1:10:55 AM6/18/15
to stan...@googlegroups.com

> On Jun 17, 2015, at 8:33 PM, Joshua N Pritikin <jpri...@pobox.com> wrote:
>
> On Wed, Jun 17, 2015 at 04:08:26PM -0400, Bob Carpenter wrote:
>>> My sense is that we haven't fully explored what is possible here.
>> ...
>>
>> I'm very much interested in a language for metaprogramming
>> models and it'd make sense to have that be native Python.
>>
>> What I really want is modularity --- some way to take reasonable
>> programming practice to designing a model where you can test
>> a component as a module before plugging it into something larger.
>> Right now, our programs are big blobs.
>>
>> But it'd be fun to brainstorm ideas in general as to what we
>> could do with metaprograms.
>>
>> Does anyone know what people do with this capability in other
>> environments?
>
> That's one of the primary design ideas of OpenMx,
> http://openmx.psyc.virginia.edu/

So what kind of things do you do with the model as object?

- Bob

P.S. Congrats for getting onto CRAN! We're still working on it for Stan.

Joshua N Pritikin

unread,
Jun 18, 2015, 1:23:04 AM6/18/15
to stan...@googlegroups.com
On Thu, Jun 18, 2015 at 01:09:15AM -0400, Bob Carpenter wrote:
> > On Jun 17, 2015, at 8:33 PM, Joshua N Pritikin <jpri...@pobox.com> wrote:
> > On Wed, Jun 17, 2015 at 04:08:26PM -0400, Bob Carpenter wrote:
> >>> My sense is that we haven't fully explored what is possible here.
> >> ...
> >>
> >> I'm very much interested in a language for metaprogramming
> >> models and it'd make sense to have that be native Python.
> >>
> >> What I really want is modularity --- some way to take reasonable
> >> programming practice to designing a model where you can test
> >> a component as a module before plugging it into something larger.
> >> Right now, our programs are big blobs.
> >
> > That's one of the primary design ideas of OpenMx,
> > http://openmx.psyc.virginia.edu/
>
> So what kind of things do you do with the model as object?

Here's a relevant excerpt from Pritikin, Hunter, & Boker (2015):

One of the main contributions of \OpenMx has been to encourage
treatment of statistical models and their components as things that
can be generated and manipulated within a programming environment.
The \OpenMx IFA module also uses this style of model specification.
Not all users will prefer to exploit such a low-level user interface.
Our goal is to extend the range of utility in both directions,
toward both expert and novice users.
Sophisticated users are empowered to build non-standard models
or to simplify application of routine analyses by
developing higher level interfaces.
For example, the \texttt{metaSEM} package \parencite{metaSEM} simplifies
the specification of \OpenMx models for meta-analysis.
A similar package could simplify the specification of
nominal testlets to make the technique simpler for non-specialists
to put into practice.

> P.S. Congrats for getting onto CRAN! We're still working on it for
> Stan.

Yeah, it was a major effort. Whew! ;-)

Matt Hoffman

unread,
Jun 18, 2015, 11:33:19 AM6/18/15
to stan...@googlegroups.com
Regarding speed, Theano has a couple of potential advantages:
• Symbolic differentiation, which allows for some automatic pre-compile-time optimizations. It also makes it easyish to compute higher-order derivatives, diagonals of Hessians, etc., since it's just running things through the gradient operator again.
• GPU/BLAS support.

Matt

Bob Carpenter

unread,
Jun 18, 2015, 1:21:30 PM6/18/15
to stan...@googlegroups.com

> On Jun 18, 2015, at 11:33 AM, Matt Hoffman <mdho...@cs.princeton.edu> wrote:
>
> Regarding speed, Theano has a couple of potential advantages:

> • Symbolic differentiation, which allows for some automatic pre-compile-time optimizations. It also makes it easyish to compute higher-order derivatives, diagonals of Hessians, etc., since it's just running things through the gradient operator again.

What does it do with matrices? I'd be interested
in finding a symbolic math package that could deal with
matrices as matrices.

Also, any idea if it can handle loops, conditionals, etc.?

We could do similar pre-compile-time optimizations at the
level of the generated Stan code --- it's harder with autodiff
unless we went the expression template route of Adept.

> • GPU/BLAS support.

That's a biggy. Hopefully we'll be able to plug this into
Eigen as we go --- they have some CUDA code ready to go.

- Bob

Matt Hoffman

unread,
Jun 18, 2015, 2:16:09 PM6/18/15
to stan...@googlegroups.com
Theano deals well with matrices, at least for relatively basic stuff. (It knows things like ∂log|X|=X^{-1}, but I doubt that it knows obscure things like the differential of an eigenvalue or whatever.) It's the only package I know of that does. It is not, however, very focused on display, so it's better used for computation than derivation. (Theano is really optimized for neural networks, so dealing well with matrices isn't optional for them.)

Loops are a bit of a pain—they can handle them using their scan() function, but it's not very intuitive or well documented. But if you spent a day or so trying to understand how scan() works you could probably do things. Support for conditionals is limited, since everything needs to be done at compile time.

Matt

Reply all
Reply to author
Forward
0 new messages