Current proposal for the function syntax and semantics

[Daniel Lee]

Moving email into its own thread.

On Thu, Feb 20, 2014 at 7:33 PM, Bob Carpenter wrote:

Here’s the current proposal for the function syntax and semantics:

https://github.com/stan-dev/stan/wiki/Function-Syntax-and-Semantics-Design

I was going to send it to the user list, but I thought that might just
frustrate people more than make it feel like an open process if I didn’t
take the suggestions. I do not want to add named arguments, nor do I want
to have any Stan variables brought into scope.

And of course, I need to balance ease of implementation with ease of use.

- Bob

[Michael Betancourt]

Looks great.  Having a language guy on the team certainly puts

Stan far ahead of the competition.  Some comments:

- I disagree with Marcus — having (transformed) data in scope may lead

to confusion down the road.  It’s easy enough to pass any necessary

objects into the function explicitly and it makes the users intent much clearer.

- Explicit container sizes, please!  They will make functions and models

that use them infinitely easier to understand and debug.  

- I’m not sure I understand the bit about not allowing constrained inputs/outputs.

Does this mean, for example, that functions can’t take covariance matrices 

as arguments? That seems like it would be pretty limiting.

- Given the special scope of “subroutines” can we give them a more

informative name?  Anything that makes it clear that they are the only

functions that have access to lp.  Maybe “submodel” or “model component”

or the like.  As Marcus noted, it may be best to avoid the use of “void”

if only to limit confusion.

--
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/groups/opt_out.

[Marcus Brubaker]

On Fri, Feb 21, 2014 at 5:22 AM, Michael Betancourt <betan...@gmail.com> wrote:

Looks great.  Having a language guy on the team certainly puts

Stan far ahead of the competition.  Some comments:

- I disagree with Marcus — having (transformed) data in scope may lead

to confusion down the road.  It’s easy enough to pass any necessary

objects into the function explicitly and it makes the users intent much clearer.

I'm not sure what confusion you're foreseeing here.  At any rate, I'm not wed to the idea, but it seems like it could get tedious to write functions in some cases.

 

- Explicit container sizes, please!  They will make functions and models

that use them infinitely easier to understand and debug.  

- I’m not sure I understand the bit about not allowing constrained inputs/outputs.

Does this mean, for example, that functions can’t take covariance matrices 

as arguments? That seems like it would be pretty limiting.

- Given the special scope of “subroutines” can we give them a more

informative name?  Anything that makes it clear that they are the only

functions that have access to lp.  Maybe “submodel” or “model component”

or the like.

+1.  The function vs subroutine distinction is a fairly obscure one that triggers painful flashbacks to Fortran with the occasional Pascal crossover nightmare.  I was going to mention it but couldn't come up with a better name but I like "submodel", makes it clear not only what can be done but also what it's meant for.

Also...for some comic relief: http://thequickword.wordpress.com/2014/02/16/james-irys-history-of-programming-languages-illustrated-with-pictures-and-large-fonts/ Now, back to Olympic hockey...

Cheers,

Marcus

[Michael Betancourt]

- I disagree with Marcus — having (transformed) data in scope may lead

to confusion down the road.  It’s easy enough to pass any necessary

objects into the function explicitly and it makes the users intent much clearer.

I'm not sure what confusion you're foreseeing here.  At any rate, I'm not wed to the idea, but it seems like it could get tedious to write functions in some cases.

I think it’s just easier to have the scope limited, especially when we deal with huge models.  That said, I would not argue against

having data, transformed data, parameters, and transformed parameters in scope for the subroutines.

Also...for some comic relief: http://thequickword.wordpress.com/2014/02/16/james-irys-history-of-programming-languages-illustrated-with-pictures-and-large-fonts/ Now, back to Olympic hockey…

Ha, solid.  The link, not the hockey.  Not the hockey.  -________________-

[Bob Carpenter]

On Feb 21, 2014, at 6:11 PM, Marcus Brubaker <marcus.a...@gmail.com> wrote:

> On Fri, Feb 21, 2014 at 5:22 AM, Michael Betancourt <betan...@gmail.com> wrote:
> Looks great. Having a language guy on the team certainly puts
> Stan far ahead of the competition. Some comments:
>
> - I disagree with Marcus — having (transformed) data in scope may lead
> to confusion down the road. It’s easy enough to pass any necessary
> objects into the function explicitly and it makes the users intent much clearer.
>
> I'm not sure what confusion you're foreseeing here. At any rate, I'm not wed to the idea, but it seems like it could get tedious to write functions in some cases.
>
>
> - Explicit container sizes, please! They will make functions and models
> that use them infinitely easier to understand and debug.
>
> - I’m not sure I understand the bit about not allowing constrained inputs/outputs.
> Does this mean, for example, that functions can’t take covariance matrices
> as arguments? That seems like it would be pretty limiting.

No, it just means it can't specify them as arguments.

Note that none of the Stan functions that operate on
triangular or correlation matrices or positive-constrained
variables actually require arguments of that type.

The motivation is that only variables have those primitive
types and if you have any kind of expression, it's just too hard
to do the type inference.

So it's up to the function to deal with enforcing constraints.

>
> - Given the special scope of “subroutines” can we give them a more
> informative name? Anything that makes it clear that they are the only
> functions that have access to lp. Maybe “submodel” or “model component”
> or the like.
>
> +1. The function vs subroutine distinction is a fairly obscure one that triggers painful flashbacks to Fortran with the occasional Pascal crossover nightmare. I was going to mention it but couldn't come up with a better name but I like "submodel", makes it clear not only what can be done but also what it's meant for.

"submodel" sounds good.

The basic idea in my mind is that "subroutines" applied to arguments
create statements, whereas functions applied to arguments create
expressions.

It only makes sense to me to call the subroutines in transformed parameters
or model blocks. Or maybe even just in the model block.

The other question is whether to allow functions to be used as statements,
and we'd want to do it in general. It's allowed in languages like C/C++/Java.
The only use case we have on the table is if the function threw an exception,
because there are no other side effects (other than log prob increments).

> Also...for some comic relief: http://thequickword.wordpress.com/2014/02/16/james-irys-history-of-programming-languages-illustrated-with-pictures-and-large-fonts/

Mitzi and I just woke the neighbors laughing through that one.
Especially the Wall and McCarthy bits.

And the Ministry of Silly Walks image is awesome!

> Now, back to Olympic hockey...
>

Go Canada!

- Bob

[Ben Goodrich]

On Friday, February 21, 2014 8:24:08 PM UTC-5, Bob Carpenter wrote:

The basic idea in my mind is that "subroutines" applied to arguments
create statements, whereas functions applied to arguments create
expressions.  

Is there a plan on how to grow this into something where the intermediate parameters can be the solution to an optimization problem? So if the user does something like

functions {
  real loss(real theta, x, y, z)
}
...
model {
  real theta;
  theta = 0.0;
  theta <- minimize(theta, loss, x, y, z);
}

But for this we need some way to specify that the user's loss function is to be minimized. Most libraries do this with a functor, but I guess there are other ways to do it.

Ben

[Bob Carpenter]

On Feb 22, 2014, at 3:45 AM, Ben Goodrich <goodri...@gmail.com> wrote:

> On Friday, February 21, 2014 8:24:08 PM UTC-5, Bob Carpenter wrote:
> The basic idea in my mind is that "subroutines" applied to arguments
> create statements, whereas functions applied to arguments create
> expressions.
>
> Is there a plan on how to grow this into something where the intermediate parameters can be the solution to an optimization problem? So if the user does something like

Absolutely --- the main motivation right now is fitting the diff eq
models we've been playing with, where we'll integrate rather than minimize.

> functions {
> real loss(real theta, x, y, z)
> }
> ...
> model {
> real theta;
> theta = 0.0;
> theta <- minimize(theta, loss, x, y, z);
> }
>
> But for this we need some way to specify that the user's loss function is to be minimized. Most libraries do this with a functor, but I guess there are other ways to do it.

That's what's going to have to go on internally and I'll either write
the code genrator so that loss is a functor or wrap it in something like
that in calls to solve, minimize, integrate, etc.

- Bob

[Michael Betancourt]

But for this we need some way to specify that the user's loss function is to be minimized. Most libraries do this with a functor, but I guess there are other ways to do it.

That's what's going to have to go on internally and I'll either write
the code genrator so that loss is a functor or wrap it in something like
that in calls to solve, minimize, integrate, etc.

Any chance the sub models would be separately accessible or is that going to get

too nasty?  This might be a sneaky way into getting separate prior/likelihood

specifications.

[Ben Goodrich]

On Thu, Feb 20, 2014 at 7:33 PM, Bob Carpenter wrote:

And of course, I need to balance ease of implementation with ease of use.

 
Are we envisioning that these things can do PRNG? As in do-it-yourself Gibbs sampling by drawing discrete parameters from their full-conditional distribution?

Ben

[Bob Carpenter]

I'm not sure what you mean by separately accessible.

One thing we could do is have them return a log probability.
So the idea would be that they had their own local lp__ variable
that gets incremented rather than passed in and then that gets
returned.

I'm not sure what the syntax for that would look like:

sub_model foo(vector y, matrix x, vector beta, real sigma) {
y ~ normal(x * beta, sigma);
}

(I'm just using the simplified arg structure to save typing here.)

And then it would get called by:

increment_log_prob(foo(y,x,beta,sigma));

or you could set it:

real temp;
temp <- foo(y,x,beta,sigma);

And then yes, you coudl define one submodel for the likelihood and
another for the prior. Or one for each bit of data.

- Bob

[Bob Carpenter]

Good point.

The PRNG functions are currently set up to only be accessible
in the generated quantities block. But in that block, there's
no way to read in previous values --- the generated quantities
aren't accessible from the previous iteration.

I don't know how we'd do Gibbs in the parameter block because
the parameters are controlled on the outside, not set on the inside
of a model.

Anybody have any ideas?

- Bob

On Feb 22, 2014, at 12:35 PM, Michael Betancourt <betan...@gmail.com> wrote:

>

[Daniel Lee]

Thanks, Bob!

Function Block

* Looks great.
* Can we consider including other files? I'm envisioning a set of functions that someone might want to reuse / distribute. It'll essentially be allowing users to extend the language as a user base without hacking into Stan itself.

Syntax of Function Declarations and Definitions
* Looks great.
* I'm thinking that down the line we'd want to provide some way to debug functions. I know there are print statements and exceptions, but one of these days, we should really give a way for people to call their functions stand-alone and in a testing framework. I don't think this affects the design of the functions, but maybe the generated code?

Container Sizing in Argument and Return Types
* +1 for declared sizes (first option)

Variable Scope and Local Variables
* Looks good.

Return Statements
* Thanks for working it out. This is trickier than I imagined.
* Last bullet: what if we threw a fatal exception instead? Perhaps something like "User function is not specified properly. foo(x,y,z) with x=...,y=...,z=... does not return a value."

Input and Output Validation
* My guess is that not allowing constrained inputs / outputs is easier. Is that true?
* If that's true, my suggestion is to start by implementing that first. Allowing constrained variables looks like an extension to the language and would be backwards compatiable. We could get it out with no constrained variable declarations, and if that ends up hindering us(ers), we could work on constrained inputs / outputs.

Function application is an expression
* Looks good.

Subroutines
* Are these going to be declared in the same block as functions? From Marcus's and Michael's comments, it sounds like it should before the model block. Seems like a weird place to insert some code, though...
* I think if we had these in a separate block, we could get around the necessity of having to specify the return type of "void," but maybe I'm not seeing things clearly.
* Would it make sense to have a named subroutine (or "submodel") declared inside the model block? I'm looking at your example and think it really does belong within the model block, but I might be looking at this really shortsided.

Call by Constant Reference
* Looks good.

[Bob Carpenter]

We should talk about this at the next Stan meeting I'm at
(not this Thu but next Thu).

More inline.

On Feb 23, 2014, at 2:37 PM, Daniel Lee <bea...@alum.mit.edu> wrote:

> Thanks, Bob!
>
> Function Block
> * Looks great.
> * Can we consider including other files? I'm envisioning a set of functions that someone might want to reuse / distribute. It'll essentially be allowing users to extend the language as a user base without hacking into Stan itself.
>

The easiest route would be to write a pre-processor
that read the files and concatenated them all together
to be passed into the parser.

The question is how we specify paths. Right now, we have
the inconvenience that everything has to be run from
the home directory of Stan, so I suppose a relative path to
that would make sense? Or we might do something clever with
stripping apart command line arguments to make it relative
to the file being processed?

> Syntax of Function Declarations and Definitions
> * Looks great.
> * I'm thinking that down the line we'd want to provide some way to debug functions. I know there are print statements and exceptions, but one of these days, we should really give a way for people to call their functions stand-alone and in a testing framework. I don't think this affects the design of the functions, but maybe the generated code?

I'm not sure most of our users will want to interact with
the C++ side alone.

Any idea how we'd do something like standalone testing of functions?
Ideally we'd have a nice unit testing framework, but that would require
a way to call the functions. I could imagine writing something in the
simple Stan language that compiles to C++ and including some test
functions. So users would have something like a test block after functions
that could optionally be called.

If we want to make data in scope, then the functions won't be
standalone, but will have to be coded as member functions on
a model. It's one argument for

>
> Container Sizing in Argument and Return Types
> * +1 for declared sizes (first option)
>

I'm really surprised everyone wants that option. Certainly
easiest to implement. But it's not like the way either Eigen
or C++ writes functions to deal with containers, which contain
their own sizing information. (And our variables all have their
own sizing built in because they're either std::vector or Eigen::Matrix.)

> Variable Scope and Local Variables
> * Looks good.
>
> Return Statements
> * Thanks for working it out. This is trickier than I imagined.
> * Last bullet: what if we threw a fatal exception instead? Perhaps something like "User function is not specified properly. foo(x,y,z) with x=...,y=...,z=... does not return a value."

We could do that by appending a final return statement and deal with
the "function not reachable" warning messages we'll get under compile for when they
do specify returns. Or we could do the full analysis of functions to only add
it if there's not a final return somewhere (that is, look at the end of each
statement and add an exception throw if there's not a return). I'd think rather
than doing that, we might just want to do a test at compile time and have them
rewrite the function.

> Input and Output Validation
> * My guess is that not allowing constrained inputs / outputs is easier. Is that true?

It's not only easier, it goes along with how we've done functions
already. You'll note that things like bernoulli_log don't take
constrained inputs.

It's partly because at runtime, the C++ translation of a Stan type
is unconstrained. That is, a real<lower=0,upper=1> just translates
to a double type.

We could have constrained inputs or outputs in the function specs and
then automatically enforce those constraints. But you'd still be able
to pass in any old expression, not just a variable declared as that type.
That is, if I had a function declared

real foo(real<lower=0,upper=1> theta);

then I could do something like this:

real x;
real y;
...
y <- foo(x);

and then foo() would check that x is between 0 and 1.

> * If that's true, my suggestion is to start by implementing that first. Allowing constrained variables looks like an extension to the language and would be backwards compatiable. We could get it out with no constrained variable declarations, and if that ends up hindering us(ers), we could work on constrained inputs / outputs.
>

It's not that big a deal one way or the other. The code
to do all the constraint checks is there already.

> Function application is an expression
> * Looks good.
>
> Subroutines
> * Are these going to be declared in the same block as functions? From Marcus's and Michael's comments, it sounds like it should before the model block. Seems like a weird place to insert some code, though...

I was planning on having everything up top and not depending on
data. So yes, that was my proposal.

One problem with allowing dependence on data variables (not hard
to implement with members) is that we'd probably want to allow functions
to be called in transformed data and allow them to be defined in terms
of transformed data. So you'd want to allow functions to be declared before
trnasformed data AND after transformed data.

I have a pretty strong preference to just put them all at the top and not
let them depend on data. That would also make standalone debugging a lot
easier.

> * I think if we had these in a separate block, we could get around the necessity of having to specify the return type of "void," but maybe I'm not seeing things clearly.

Let's just forget I ever mentioned "void". It's not necessary
no matter where they go. It was just there in my mind for symmetry
with the C++ that gets generated. But I can generate it in C++ w/o
having it be there in Stan.

> * Would it make sense to have a named subroutine (or "submodel") declared inside the model block? I'm looking at your example and think it really does belong within the model block, but I might be looking at this really shortsided.

I think of the functions as more outside the whole
model --- as you say, we may want to include them.

No reason to put them in the model block. The model
block will use them, but I don't see any reason to declare
or define the functions or submodels there.

- Bob

[Michael Betancourt]

I mean accessible in the API, so we would have the potential to code
up algorithms that require a prior/likelihood factorization.

[Bob Carpenter]

It'd be defined in the model, not in the API per se.

- Bob

[Bob Carpenter]

I revised the proposal w.r.t. everyone's suggestions and removed
alternatives that are now off the table. Some other changes:

* allow constrained types in declarations and return that get
checked on input and output

* changed subroutines to "submodels", removed "void" altogether
and put them in their own block to clarify

* added a proposal for functions ending in "_rng" to have access
to the PRNG base in order to allow meaningful generated quantities
functions

* added some other suggestions (solvers for roots, non-constant args,
Boost tuples [better idea than structs!]) and tried to clearly mark
things I won't implement in the first pass as FUTURE

Let me know if there are any more comments, because I'd like to
start building this as soon as the parser fixes get merged.

- Bob

[Michael Betancourt]

Everything looks great, although I am slightly concerned about the _rng
functions (unless I’m misreading your proposal). Will _rng functions be
callable anywhere or will they be limited to the generated quantities?

Also, any reason we don’t allow them in the transformed data for
random generation of pseudo data?

[Bob Carpenter]

_rng() will only be callable from generated quantities, matching
our current restrictions on the PRNGs.

We could allow RNG in transformed data, but the problem we
discussed with that before (or maybe I just discussed with
myself before we had a bigger team) was that there's no way to
get the transformed data back out.

We could probably generalize the commands somehow to allow transformed
data to have access. We'd have to push the seed through to the
model constructor, as right now it's only being passed in
later.

But that's an independent issue. The "_rng" functions will follow
the use of the other PRNGs.

I really don't like these naming conventions, by the way. Right
now, the compiler looks to see if a function ends in "_rng" and if
so, implicitly passes the PRNG base function to it.

- Bob

[Bob Carpenter]

One thing I just realized is that if we prohibit access to increment_log_prob
in functions and restrict it to submodels, there'll be no way to write
natural transform functions, as in:

real logit_inverse_transform(real logit_a) {
real a;
a <- inv_logit(logit_a);
increment_log_prob(log(a) + log1m(a));
return a;
}

Instead, we'd need to write the function and then a Jacobian
submodel and make sure to call both.

I could add lp__ as an implicit argument to functions, and
set it to a dummy double value for blocks other than transformed
parameters and model. Would that be too confusing?

- Bob

[Ben Goodrich]

Can the user-defined function end with _log and get lp__ passed to it automatically?

Ben

[Bob Carpenter]

It could, but I really think that was a very poor choice on
my part to begin with given the popularity of log! So I'd
rather choose something else if we do it explicitly, like maybe
_lp?

- Bob

[Michael Betancourt]

I like the idea and agree that we need a better name.  _dist?

[Ben Goodrich]

Sure, but I would probably leave normal_log etc. as is instead of changing it to normal_lp.

Ben

[Bob Carpenter]

I think "_dist" is misleading because the function isn't necessarily
going to be a distribution --- it might just be a transform with
a Jacobian adjustment, for example.

- Bob

[Bob Carpenter]

Yes, I plan to leave "normal_log" as is. The only thing that happens
if the user adds "_log" to one of their functions,

real foo_log(T y, ...)

is that

y ~ foo(...)

becomes synonymous with

increment_log_prob(foo_log(y,...))

give or take some constant terms.

- Bob

[Michael Betancourt]

I know it's too late, but overloading log like that is confusing.  "prob" or "dist" would have been better.

Anyways, this seems like an expert-level feature and I can't think of a way that you could code it

where the function would work on both the LHS and the RHS of a sampling statement.

[Bob Carpenter]

On Mar 4, 2014, at 6:57 PM, Michael Betancourt <betan...@gmail.com> wrote:

> I know it's too late, but overloading log like that is confusing. "prob" or "dist" would have been better.

Agreed --- very confusing.

I added "_log" at the end to indicate they were on the log scale
and then didn't feel like adding yet another decoration on the functions.
Given that we don't do anything on the linear scale, I probably could've
gotten rid of the "_log" and just had that be the default and then gone
with "_dist".

> Anyways, this seems like an expert-level feature and I can't think of a way that you could code it
> where the function would work on both the LHS and the RHS of a sampling statement.

Right.

The topic I wanted feedback on, though, wasn't this one, but how or even
whether to deal with functions that want to manipulate the log probability
accumulator.

- Bob

[Michael Betancourt]

>
> The topic I wanted feedback on, though, wasn't this one, but how or even
> whether to deal with functions that want to manipulate the log probability
> accumulator.

I’m just not sure what you would need such a function for that you couldn’t
accomplish with increment_log_prob(function(inputs));

[Bob Carpenter]

That's why I included the transform example that returns the transformed
parameter and increments the log probability with the Jacobian adjustment.
Then you could write:

transform(x,...) ~ ...;

instead of

increment_log_prob(log_jacobian_transform(x,...));
transform(x,...) ~ ...

Right now, that would be useful in cases of vectors with varying
lower and upper bound constraints, which we don't currently support.

- Bob

[Michael Betancourt]

I dunno, that seems like a really error prone approach —
it’s really easy to mess it up and there would be little diagnostics
to identify the problem. I would personally prefer to not encourage
user to do such things by hand.

[Ben Goodrich]

The limitation that lower and upper bounds can only be scalars is pretty limiting in my opinion.

Ben

[Bob Carpenter]

It could, but I really think that was a very poor choice on
my part to begin with given the popularity of log! So I'd
rather choose something else if we do it explicitly, like maybe
_lp?

- Bob

On Mar 4, 2014, at 4:09 PM, Ben Goodrich <goodri...@gmail.com> wrote:

[Daniel Lee]

What if you had to explicitly list it as an argument?

[Bob Carpenter]

Won't work if we pass by const reference. Also, direct
access of lp__ is deprecated. We decided we'd rather not
give users that much rope with which to hang themselves. (They
[and we] were making mistakes like setting rather than incrementing,
thereby bypassing the Jacobians on the constraints.)

- Bob

[Ben Goodrich]

On Wednesday, March 5, 2014 12:12:27 PM UTC-5, Bob Carpenter wrote:

Won't work if we pass by const reference.  Also, direct
access of lp__ is deprecated.  We decided we'd rather not
give users that much rope with which to hang themselves.  (They
[and we] were making mistakes like setting rather than incrementing,
thereby bypassing the Jacobians on the constraints.)

 
On that point, the parser warning is still phrased in terms of the old way

DIAGNOSTIC(S) FROM PARSER:

Warning (non-fatal):     sampling statement (~) contains a transformed parameter or local variable.
     You must increment lp__ with the log absolute determinant of the Jacobian of the transform.
     Sampling Statement left-hand-side expression:
          sigma2 ~ inv_gamma_log(...)

So, it should say something to the effect of: You need to do increment_log_prob(LOG_ABS_JACOBIAN_DETERMINANT);

Ben