pushed bracket changes

[Bob Carpenter]

I did the global search and replace on both the
doc and all the models. The unit tests now all work
again. I skimmed an output report and did some other
investigations and I think the models should now all be
OK, too.

The new syntax for dimension declarations, interval
constraints, and truncated probability functions is as
follows:


real a[3,4,5];

vector[4] b;

matrix[4,5] c;

real[-1,1] rho;

vector[4] d[3,2];

a ~ normal(0,1) T[0,];

I didn't check the Java program in, but I have it in
case anyone else needs it for something.

- Bob

[Andrew Gelman]

B
As I sent in my last email msg to the group, I'm not so thrilled with notation such as real, vector, matrix (I assume we're using "array" for arrays of more than 2 dimension?). I'd prefer if we could just use "real" for real, "integer" for integer, etc. Perhaps we can discuss, this isn't something that is necessary perhaps, I'd be interested in hearing the rationale for doing it the way you're doing it now.
A

[Bob Carpenter]

On 8/8/12 3:29 PM, Andrew Gelman wrote:
> B
> As I sent in my last email msg to the group,

Could you resend -- I think I missed it.

> I'm not so thrilled with notation such as real, vector, matrix (I assume we're using "array"
> for arrays of more than 2 dimension?). I'd prefer if we could just use "real" for
> real, "integer" for integer, etc. Perhaps we can discuss, this isn't something
> that is necessary perhaps, I'd be interested in hearing the rationale for doing it the way
> you're doing it now.

The rationale for "int" vs. "integer" is that it's
conventional in every programming language. You
already had us change "double" to "real" and it's no
problem to change "int" to "integer".

NOW is the time to do it. I will not want to change
these kinds of superficial details after our first stable
release because it'd require our users to change all their
models.

There are full details of variable types and their
declarations in the variable declaration
chapter of the doc. Could you maybe take a look through
the doc and send us other comments you have about syntax?

Here's some examples of what things look like now
for simple integer, real, and array declarations:

real a; // unconstrained real
int b; // unconstrained integer

real c[5]; // 5-array of reals
real d[5,6]; // 5 x 6-array of reals

int e[3,4,5]; // 3 x 4 x 5 array of integers

Vectors and matrices have the following syntax:

vector[5] f; // 5-vector (of reals)
row_vector[4] g; // 4-row-vector (of reals)
matrix[2,3] h; // 2 x 3-vector (of reals)

In particular, note that a 2x3 matrix is NOT the same
thing as a 2x3 array. This distinction is built in at
a fundamental level.

The basic data types are:

int (with optional range constraints)
real (with optional range constraints)
vector (require # of rows)
simplex (require # of rows)
ordered (require # of rows)
row_vector (require # of cols)
matrix (require # of rows and cols)
corr_matrix (require # of rows/cols)
cov_matrix (require # of rows/cols)

The indentation is on purpose -- simplex and ordered types are
implemented as vector types and corr_matrix and cov_matrix types
are instances of matrix. So you can assign a vector to a simplex
or ordered type, but the constraint will be enforced at the end
of the block in which the variable is declared.

Arrays of any dimensionality can be made out of any of
the basic types. So, for instance,

simplex[4] theta[10,5];

declares theta to be a 10 x 5 array of 4-simplexes.
theta[7] is a 5-array of 4-simplexes, theta[7,2] is
a 4-simplex, and theta[7,2,1] is a real value between
0 and 1 (inclusive).

The reason we need to declare the constraints on parameters
is that we need to do the transforms. The reasons there
are constraints elsewhere is for error checking.

- Bob

[Andrew Gelman]

Hi, Bob. I don't care about "int" vs. "integer". If "int" is standard, let's do "int." I'm also happy with the array notation you have below ("int e[3,4,5]", etc). I thikn that sometimes it seems like I'm being picky but really I'm just not understanding what's going on. You can ignore my earlier email to the list regarding constraints, I'd misunderstood what was happening.

The thing that had been freaking me out were the declarations "vector" and "matrix". I just didn't understand why these were not subsumed under "real".

But now I understand that you can have an array of vectors or of matrices. That's great. It's so powerful, I didn't even realize it was possible. Cool!

This also implies, I suppose, that there are multiple ways of declaring variables. For example, in the 8 schools, we could say "real theta[8]" or we could say "vector[8] theta". Is that correct?

Anyway, this is all super-cool, and it's all new to me. I was so used to Fortran style that I just missed completely what you were doing here.

I'm assuming that everything is real by default, i.e."matrix" and "vector" are real, and if we ever want to have matrices or vectors of integers (e.g., for latent class models), we'd need to define a new type such as "integer_matrix" or "integer_vector". Is that correct? Also, I could imagine that it might make sense to define a "binary" class rather than call it integer[0,1]? But I don't really know.

Regarding constraints, again I defer to you but here are my thoughts:

- I worry that angle brackets will be confusing because it's one more thing for the user to remember: "use square brackets everywhere, except for bounds where you need angle brackets." I prefer the R-style convention of indicating arguments via names rather than punctuation. Thus, for example,

matrix[floor=0,ceiling=1] x[3,4,5];
or
matrix[bounds=[0,1]] x[3,4,5];
or
matrix[floor=0] x[3,4,5];
or
matrix[bounds=[0,]] x[3,4,5];

What do you think?

A

[Bob Carpenter]

On 8/8/12 4:46 PM, Andrew Gelman wrote:
> Hi, Bob. I don't care about "int" vs. "integer". If "int" is standard, let's do "int." I'm also happy with the array notation you have below ("int e[3,4,5]", etc). I thikn that sometimes it seems like I'm being picky but really I'm just not understanding what's going on. You can ignore my earlier email to the list regarding constraints, I'd misunderstood what was happening.
>
> The thing that had been freaking me out were the declarations "vector" and "matrix". I just didn't understand why these were not subsumed under "real".
>
> But now I understand that you can have an array of vectors or of matrices. That's great. It's so powerful, I didn't even realize it was possible. Cool!

I feel your pain. You and Michael M. and everyone else
was very patient with me getting up to speed on stats
notation and terminology.

> This also implies, I suppose, that there are multiple ways of declaring variables. For example, in the 8 schools, we could say "real theta[8]" or we could say "vector[8] theta". Is that correct?

Right, either would work. As would row_vector. But
they're not totally interchangeable. Some contexts require
vectors or matrices, such as vector/matrix multiplication
and parameters/outcomes for probabilty functions.

My feeling was that we'd use vectors/matrices whenever
we wanted to apply vector or matrix operations or whenever
we need parameters for multivariate distros.

The main reason for having these types is that we can
infer that for

row_vector * vector

the result type is real.

> Anyway, this is all super-cool, and it's all new to me. I was so used to Fortran style that I just missed completely what you were doing here.
>
> I'm assuming that everything is real by default, i.e."matrix" and "vector" are real, and if we ever want to have matrices or vectors of integers (e.g., for latent class models), we'd need to define a new type such as "integer_matrix" or "integer_vector". Is that correct? Also, I could imagine that it might make sense to define a "binary" class rather than call it integer[0,1]? But I don't really know.

Yes, all of the matrix types are real only. We could extend
that later as Eigen supports any type as long as the basic
scalar operations are supported.

Do you think we'll have cause to use matrix operations on
integer-based matrices?

> Regarding constraints, again I defer to you but here are my thoughts:
>
> - I worry that angle brackets will be confusing because it's one more thing for the user to remember: "use square brackets everywhere, except for bounds where you need angle brackets." I prefer the R-style convention of indicating arguments via names rather than punctuation. Thus, for example,
>
> matrix[floor=0,ceiling=1] x[3,4,5];
> or
> matrix[bounds=[0,1]] x[3,4,5];
> or
> matrix[floor=0] x[3,4,5];
> or
> matrix[bounds=[0,]] x[3,4,5];
>
> What do you think?

There aren't any angle brackets yet -- we're just making
sure that the current language is extendable to something
reasonable.

We could start naming arguments. The only issue is verbosity.
So we could also have an order-only usage like in R.

A bigger issue is moving the matrix sizes into
the array declaration. What we have now is:

matrix[4,5] x[3];

for a 3-array of 4x5 matrices to

matrix x[3,4,5]

This is possible, but I worry that it won't be clear where
the matrix is. An alternative would be to put the named arguments
on matrix:

matrix[rows=4, cols=5] x[3];

or

matrix<rows=4, cols=5> x[3];

The angle brackets are traditional in C++ and Java for type
modifiers. Using separate notation for different
things is easier to read, if not easier to remember. (Brackets
are very ugly typographically for some reason I don't understand.)

If we do go with names, I'd prefer "lb" and "ub" to
"floor" and "ceiling", because the latter have integer-rounding
interpretations. So that'd give us something like:

real[lb = -1, ub = 1] x[12,7];

for a 12x7 array of reasl bounded in [-1,1]. Then we'd be
able to write

real[-1,1] x[12,7];

as we do now. It'd also open the way for us to write:

real[lb = -1] x[12,7];

instead of

real[-1,] x[12,7];

And then maybe instead of writing:

simplex[5] theta[3,4];

for a 3x4 array of 5-simplexes, we could write:

matrix[simplex,rows=5] theta[3,4];

or

matrix[constraint=simplex, rows=5] theta[3,4];

and infer the number of columns matches the number of rows?

- Bob

[Marcus Brubaker]

I think ideally it would be good to do something more like you're expecting from Fortran but right now we don't have a general multidimensional array object which can map efficiently onto an Eigen::Matrix object.  The practical side of that is the following:

real theta1[8]; // theta1 is an array

and

vector[8] theta2; // theta2 is a vector

results in objects of two different types which are handled differently.  For instance, you can't write

real A[8,8];
b = A * theta1;

but you can write

matrix[8,8] A;
b = A * theta2;

This is the matrix-vs-array distinction that Bob mentioned in the previous email.

Basically, if you want to do linear algebra, you have to use a vector or a matrix type.  If you don't care about that and are only accessing elements directly then you can use vector, matrix or an array.

If at some point in the future we implement a multidimensional array object which can interface with Eigen::Matrix, then the vector and matrix types can become aliases to the array class and declaring bounds on elements becomes consistent across the board.  This is part of why I was thinking of using matrix<real[0,1]> because eventually this notation would no longer be needed so there's no point making the real[0,1] case ugly.

Also, the angle brackets are for declaring a type on a matrix, the bound would still be with square brackets.  In other words: use square brackets for applying bounds to variables, and use angle brackets to declare the element type of a matrix.  Also note that this latter case (the angle brackets) would only be necessary if someone wanted to include a bound on the elements of the matrix.

Cheers,
Marcus

[Marcus Brubaker]

On Wed, Aug 8, 2012 at 5:11 PM, Bob Carpenter <ca...@alias-i.com> wrote:

On 8/8/12 4:46 PM, Andrew Gelman wrote:

This also implies, I suppose, that there are multiple ways of declaring variables. For example, in the 8 schools, we could say "real theta[8]" or we could say "vector[8] theta".  Is that correct?

Right, either would work.  As would row_vector.  But
they're not totally interchangeable.  Some contexts require
vectors or matrices, such as vector/matrix multiplication
and parameters/outcomes for probabilty functions.

My feeling was that we'd use vectors/matrices whenever
we wanted to apply vector or matrix operations or whenever
we need parameters for multivariate distros.

The main reason for having these types is that we can
infer that for

    row_vector * vector

the result type is real.

I know we've talked about this before but I forget the results.  Couldn't we do that by making a 1x1 matrix implicitly convertible to a scalar and vice-versa?

A bigger issue is moving the matrix sizes into
the array declaration.  What we have now is:

  matrix[4,5] x[3];

for a 3-array of 4x5 matrices to

  matrix x[3,4,5]

This is possible, but I worry that it won't be clear where
the matrix is.

I would strongly argue against this.  I think it would be confusing and lead to a lot of problems.  For instance, with the above, x[1,:,:] would be a 4x5 matrix but what would x[:,:,1] be?  If it's a 3x4 matrix, then it would involve a lot of unpleasant copying.  So long as internally a matrix has to be 2D, that should be explicit to users, like it is now.

In a world with more hours in the day I would go and implement a multidimensional array with a Map for Eigen::Matrix to allow us to just do:

real x[3,4];

instead of

matrix[3,4] x;

or

real[0,1] x[7,8];

and be done with it :)

At any rate, since this is all started based on a feature that's not yet implemented (vectors and matrices of bounded variables) and is not likely to be widely used at first I think we can let it rest for the near term.  There are options which are consistent with the current notation (matrix<real[0,1]>[3,4] for instance) and don't cause parsing problems.  And, since it relates to the larger issue of multidimensional arrays, it may get solved in another way at some point in the future.

Cheers,
Marcus

[Bob Carpenter]

On 8/8/12 5:29 PM, Marcus Brubaker wrote:

> On Wed, Aug 8, 2012 at 5:11 PM, Bob Carpenter <ca...@alias-i.com <mailto:ca...@alias-i.com>> wrote:
> The main reason for having these types is that we can
> infer that for
>
> row_vector * vector
>
> the result type is real.
>
>
> I know we've talked about this before but I forget the results. Couldn't we do that by making a 1x1 matrix implicitly
> convertible to a scalar and vice-versa?

I'm not saying it's impossible, but I don't see how
we could make it work. It would require a richer
type-inference scheme than what we have, which only
casts an int up to a real. Given the way code's generated,
it'd be easiest if we could do it on the C++ side.
But there are several issues (see below).

Overall, though, I prefer stricter typing with less
coercion. I think it makes the programs easier to understand
and debug.

What we've done in the probability functions is make
a vector view wrapper so that we can treat a scalar
the same way as a vector. But if we started allowing

> A bigger issue is moving the matrix sizes into
> the array declaration. What we have now is:
>
> matrix[4,5] x[3];
>
> for a 3-array of 4x5 matrices to
>
> matrix x[3,4,5]
>
> This is possible, but I worry that it won't be clear where
> the matrix is.
>
>
> I would strongly argue against this. I think it would be confusing and lead to a lot of problems. For instance, with
> the above, x[1,:,:] would be a 4x5 matrix but what would x[:,:,1] be?

Right now, x[1] is well formed, but neither x[1,,] nor x[1,:,:] nor x[:,:,1] is.
I started going down that route. It would involve either more views or
copying. Either way, it would be nice to support.

> If it's a 3x4 matrix, then it would involve a lot
> of unpleasant copying. So long as internally a matrix has to be 2D, that should be explicit to users, like it is now.

That's my preference, too.

> In a world with more hours in the day I would go and implement a multidimensional array with a Map for Eigen::Matrix to
> allow us to just do:
>
> real x[3,4];
>
> instead of
>
> matrix[3,4] x;
>
> or
>
> real[0,1] x[7,8];
>
> and be done with it :)

This isn't quite going to yield the same result because
it doesn't distinguish row vectors from column vectors.
The problem I see is with declarations

real x[5];
real y[5];

the expression (x * y) is ambiguous. If x is a row
vector and y a column vector, then (x * y) is a scalar.
If x is a column vector and y a row vector, then
(x * y) is a matrix. Or maybe we'd only promote 1D arrays
to column vectors? Then (x * y) would be undefined, but
(x' * y) and (x * y') would be well-defined.

And this isn't just a type issue of distinguishing a 1x1
matrix from a scalar -- there are different result shapes!

At the underlying storage level it would make sense
to do our own implementation rather than using std::vector.
It'd be faster and less bulky.

Then we could use as_matrix(x), as_vector(x), or
as_row_vector(x) to create a matrix view of the underlying
data wherever necessary. I don't know how hard it would be
to implement a new Eigen-like matrix wrapper (I haven't studied
how Map is implemented -- looking into Eigen always makes my
head spin).

In Stan, we could also just do the copying. So we could have:

real[,] to_array(matrix);
real[] to_array(vector);
real[] to_array(row_vector);

matrix to_matrix(real[,]);
vector to_vector(real[]);
row_vector to_row_vector(real[]);

With the way we're using Eigen on an operation-by-operation basis
(I know, Ben, this isn't optimal!), this might not make that
much of a difference performance-wise.

> At any rate, since this is all started based on a feature that's not yet implemented (vectors and matrices of bounded
> variables) and is not likely to be widely used at first I think we can let it rest for the near term. There are options
> which are consistent with the current notation (matrix<real[0,1]>[3,4] for instance) and don't cause parsing problems.
> And, since it relates to the larger issue of multidimensional arrays, it may get solved in another way at some point in
> the future.

We just need to nail down what the basic syntax is.
I guess we've found our bike-shed issue

http://en.wikipedia.org/wiki/Parkinson's_Law_of_Triviality

But then getting basic syntax right is important for
understandability of the code.

- Bob

[Andrew Gelman]

B
See comments below.
A

On Aug 8, 2012, at 11:11 PM, Bob Carpenter wrote:

> The main reason for having these types is that we can
> infer that for
>
> row_vector * vector
>
> the result type is real.

If you think that's a good idea, we can go with it. My inclination (based on my background from R) is to _not_ allow just plain "*" to do this. I would prefer something more explicitly such as
inner (a, b);
which would return the inner product of a and b whether or not they are both row_vectors, both vectors, or one of each (in either order). Actually, I would prefer to just have "vector" and forget about "row_vector" but I could see (perhaps) how row_vector will help when we're doing matrix multiplication.

I'm not saying to go change everything right now; it just seems like it might be worth discussing. My worry about "*" is the potential confusion, I see the code as less readable because when I see "*" I will have to go back to declarations to see what is what. "inner(,)" seems more clear to me.

>
> Do you think we'll have cause to use matrix operations on
> integer-based matrices?

Possibly. This comes up in some models in psychometrics (and, I assume, in CS) with latent classes. Then again, we don't know how we will be fitting such models, as HMC doesn't work there!

> There aren't any angle brackets yet -- we're just making
> sure that the current language is extendable to something
> reasonable.

I thought you were proposing angle brackets for constrained parameters.

> We could start naming arguments. The only issue is verbosity.
> So we could also have an order-only usage like in R.
>
> A bigger issue is moving the matrix sizes into
> the array declaration. What we have now is:
>
> matrix[4,5] x[3];
>
> for a 3-array of 4x5 matrices to
>
> matrix x[3,4,5]

No, I am going back to your original plan. I think "matrix[4,5] x[3]" is fine! That was the point of my past email ("wow, cool!"): now I understand the maotivation for "matrix[4,5] x[3]" and I love it! Actually, I could imagine uses for "array[4,5,6] x[7,8]" but maybe we can wait on that…

What I'm getting at is, your notation "vector" and "matrix" has two distinct advantages. First is your motivation of computational efficiency: vector and matrix types allow you to implement fast vector computations. Second is for statistical notation. For instance, imagine a hierarchical model in which each person has a vector of parameters theta of length 4, and there are 40 people. Then you could do "vector[4] theta[40]".

It's the Pinnochio principle all over again: your comptuational idea, when implemented, suggested a statistical idea.

Now that I'm thinking about it, there are other ways of implementing the statistical idea. For example, suppose in my example we simply defined "real theta[40,4]" but then we wanted to give names to the rows or columns of theta. We could use a notation such as "theta_1" for the rows (thus theta_1 is a set of 40 vectors, each of length 4) and "theta_2" for the columns. This would generalize, for example if we have "real x[3,4,5]", then x_13 is a set of 4 matrices, each of dimension 3x5.

I could see how this general notation could be useful at some point but it doesn't seem worth worrying about now.

>
> This is possible, but I worry that it won't be clear where
> the matrix is. An alternative would be to put the named arguments
> on matrix:
>
> matrix[rows=4, cols=5] x[3];
>
> or
>
> matrix<rows=4, cols=5> x[3];

No, I _don't_ think that will be needed. It seems clear enough in "matrix[4,5] x[3]" that 4 and 5 are the dimensions of the matrix.

The only place where I wanted names was for bounds.

>
> If we do go with names, I'd prefer "lb" and "ub" to
> "floor" and "ceiling", because the latter have integer-rounding
> interpretations. So that'd give us something like:
>
> real[lb = -1, ub = 1] x[12,7];
>
> for a 12x7 array of reasl bounded in [-1,1]. Then we'd be
> able to write
>
> real[-1,1] x[12,7];

That could work, or instead of "lb" and "ub" (which to me seem too cute and unix-like, in the mode of "cat"), we could use "bounds", hence "real[bounds=[-1,1]] x" or even "real[bounds=-1,1] x" etc.

>
> as we do now. It'd also open the way for us to write:
>
> real[lb = -1] x[12,7];
>
> instead of
>
> real[-1,] x[12,7];
>
> And then maybe instead of writing:
>
> simplex[5] theta[3,4];
>
> for a 3x4 array of 5-simplexes, we could write:
>
> matrix[simplex,rows=5] theta[3,4];
>
> or
>
> matrix[constraint=simplex, rows=5] theta[3,4];
>
> and infer the number of columns matches the number of rows?
>

Shortcuts are fine but unambiguity is better, I think. The real issue, though, is that once we're fitting big models (e.g., the 12 predictors with all the 2-way interactions), I think we'll want some way for the user (i.e., me) to set up shortcuts, or subroutines, or macros, or whatever. The important gains will not be from removing keystrokes or characters within a line, the important gains will come from chunking in some way.

[Andrew Gelman]

I'd prefer to reserve "*" for component-wise multiplication and use "inner" for dot products.  But perhaps I'm being shortsighted if there are models with lots of matrix operations?  In my experience, I do a lot of componentwise multiplication and I don't mind typing "inner" in the relatively rare cases that I'm doing a dot product.  But maybe things are differnent for otehr users?

A

[Bob Carpenter]

On 8/9/12 3:36 AM, Andrew Gelman wrote:>

I prefer named bounds but not named rows and columns. Thus, for example:
>
> int[bounds=0,] c
> int[bounds=0,1] d
> real[bounds=-1,1] f[4,5]
> matrix[bounds=-1,1][4,5] m

I see the advantage of this syntax. As things stand now,
there are no names and bounds can't be declared on vectors,
row vectors, or matrices, only real and int. (Yes, it's on
the to-do list.)

Should I go through and change this before the first release?
I can add names later, but won't be able to get rid of them w/o
breaking backward compatibility of user models.

The only thing to remember here is which order to
put the bounds in relative to the sizes.

> Another reason I don't like the "rows" and "columns"
> thing is what happened to R: different functions use
> "row", "rows", "nrow", "nrows", "col.means", "colnames", etc

That's also a problem with naming functions.

> As a non-C-user, I'd go with square brackets everywhere so that users don't
> have to try to remember when they need angle brackets. But if you think
> that's silly, we can go with the angle brackets!

I like different syntax for different things, but there
are only so many symbols on an ASCII keyboard.

On 8/9/12 3:29 AM, Andrew Gelman wrote:
> I'd prefer to reserve "*" for component-wise multiplication and use "inner" for dot products. But perhaps I'm being
> shortsighted if there are models with lots of matrix operations? In my experience, I do a lot of componentwise
> multiplication and I don't mind typing "inner" in the relatively rare cases that I'm doing a dot product. But maybe
> things are differnent for otehr users?

I somehow missed a bunch of messages, so I'm only replying now.

We use ".*" for component-wise multiplication and "*"
for regular matrix multiplication.

If x is a row vector and y is a vector, then (x * y)
is a real and (y * x) is a matrix. In general, the
result is reduced to a smaller type.

- Bob

[Bob Carpenter]

On 8/9/12 3:25 AM, Andrew Gelman wrote:

>> The main reason for having these types is that we can
>> infer that for
>>
>> row_vector * vector
>>
>> the result type is real.
>
> If you think that's a good idea, we can go with it.

Great! I do think it's a good idea.

> My inclination (based on my background from R) is to _not_
> allow just plain "*" to do this. I would prefer something more explicitly such as
> inner (a, b);
> which would return the inner product of a and b whether or not they are both row_vectors, both vectors, or one of each (in either order). Actually, I would prefer to just have "vector" and forget about "row_vector" but I could see (perhaps) how row_vector will help when we're doing matrix multiplication.
>
> I'm not saying to go change everything right now;

I'd say it's now or never. If we change things later,
we'll break every model that a user has written.

> it just seems like it might be worth discussing. My worry about "*" is the
> potential confusion, I see the code as less readable because when I see "*" I
> will have to go back to declarations to see what is what. "inner(,)" seems more clear to me.

The way we did typing of vector, row vector and matrix
follows Eigen, which I really like. I also like that
you get more precise type inference for expressions.

Other people prefer a more R-like way of doing things
where all the types are very squishy and interconvertible
based on what's calling them.

We could've done that in Stan. We could still do it now,
but it's going to take A LOT of work, because it touches
every aspect of Stan from the model parser to the documentation
to the models.

I haven't tried to formulate a coherent alternative
typing. A lot of what R has to do with inner and what not
is because of the implicit type conversions they allow.

Our matrix ops are more like MATLAB.

".*" is componentwise multiplication, and * is always always
always matrix multiplication reducing to simplest type.

(matrix * matrix):matrix
(real * matrix):matrix
(real * vector):vector
(row_vector * vector):real
(vector * row_vector):matrix

etc.

The simpler types always recombine properly under
multiplication in compound expressions. But NOT in
assignment.

So I can't do

matrix(1,1) x;
row_vector(K) rv;
vector(K) v;
x <- rv * v; // (rv * v) is type real

This we could add easily enough.

We could get rid of row_vector and vector and just
have matrix. Then we'd treat row vectors and vectors
as 1 x N and M x 1 matrices. And when you multiply
a 1 x M matrix by an M x 1 matrix, you get a 1 x 1 matrix.

>
>>
>> Do you think we'll have cause to use matrix operations on
>> integer-based matrices?
>
> Possibly. This comes up in some models in psychometrics (and, I assume, in CS) with latent classes. Then again, we don't know how we will be fitting such models, as HMC doesn't work there!

That we can extend to later easily, too.
We'd just have to treat real as the default.

>> There aren't any angle brackets yet -- we're just making
>> sure that the current language is extendable to something
>> reasonable.
>
> I thought you were proposing angle brackets for constrained parameters.

Just for matrices, so to declare a as an M x N matrix
of non-negative real values, it'd be:

matrix<real[0,]>[M,N] a;

But I like

matrix[bounds=0,][M,N] a;

better, so let's just go with that (pending the decision
to completely redesign the language!).

- Bob

[Andrew Gelman]

Bob:
I defer to you on all of the items below. What will happen is, once I start fitting models, I'll come up with features I will want. Some of those will already exist in Stan, others I will want. I am hoping/expecting that we will be able to add these new things by adding new model classes or macros or whatever they are called, without needing to overhaul the language.
To put it another way, I don't think the language can handle all the things we'll want, so we'll bave to define new models or distributions as needed. One of my frustrations with Bugs was having to repeat blocks of code because there were no subroutines or macros.
A

> --
>

[Andrew Gelman]

See responses below.

On Aug 15, 2012, at 4:05 AM, Bob Carpenter wrote:

> On 8/9/12 3:36 AM, Andrew Gelman wrote:>
>
> I prefer named bounds but not named rows and columns. Thus, for example:
> >
> > int[bounds=0,] c
> > int[bounds=0,1] d
> > real[bounds=-1,1] f[4,5]
> > matrix[bounds=-1,1][4,5] m
>
> I see the advantage of this syntax. As things stand now,
> there are no names and bounds can't be declared on vectors,
> row vectors, or matrices, only real and int. (Yes, it's on
> the to-do list.)
>
> Should I go through and change this before the first release?
> I can add names later, but won't be able to get rid of them w/o
> breaking backward compatibility of user models.

What do the others think? My inclination is that names for bounds are a good thing. But I'm willing to be outvoted on this one.

> The only thing to remember here is which order to
> put the bounds in relative to the sizes.
>
> > Another reason I don't like the "rows" and "columns"
> > thing is what happened to R: different functions use
> > "row", "rows", "nrow", "nrows", "col.means", "colnames", etc
>
> That's also a problem with naming functions.
>
> > As a non-C-user, I'd go with square brackets everywhere so that users don't
> > have to try to remember when they need angle brackets. But if you think
> > that's silly, we can go with the angle brackets!
>
> I like different syntax for different things, but there
> are only so many symbols on an ASCII keyboard.

My worry with angle brackets is that it's just one more thing to remember: when to use angle brackets and when to use square brackets. I find it convenient to convey info through names rather than bracket types.

>
> On 8/9/12 3:29 AM, Andrew Gelman wrote:
>> I'd prefer to reserve "*" for component-wise multiplication and use "inner" for dot products. But perhaps I'm being
>> shortsighted if there are models with lots of matrix operations? In my experience, I do a lot of componentwise
>> multiplication and I don't mind typing "inner" in the relatively rare cases that I'm doing a dot product. But maybe
>> things are differnent for otehr users?
>
> I somehow missed a bunch of messages, so I'm only replying now.
>
> We use ".*" for component-wise multiplication and "*"
> for regular matrix multiplication.
>
> If x is a row vector and y is a vector, then (x * y)
> is a real and (y * x) is a matrix. In general, the
> result is reduced to a smaller type.

Persnally, I hate the idea of "row vectors" and "column vectors". I always find myself fumbling to keep these straight. But I get the impression that this matrix notation is a huge mess, so 'm guessing I should defer to you and Ben on this stuff.

A

[Bob Carpenter]

Changing the approach to matrix typing at this point
would be a HUGE hassle. I also have a strong preference
for the way it is. One thing we could do is add
other functions. For instance, there's already a dot
function that will do dot products of any pair of vectors
(not arrays). This again follows Eigen's lead.

I'm ambivalent about naming bounds or rows/columns.
Since Andrew feels strongly about naming just
bounds, why don't we do that? The change to the parser
is trivial and I can pretty automatically change the
models. Just to remind everyone, that's:

real[bounds=0,1] x;

instead of

real[0,1] x;

Basically, the syntax just inserts a "bounds" and then "="
token before the bounds. So that's trivial.

It also clears the way for doing matrices in a reasonable
way, with

matrix[bounds=-1,1][3,4] x;

I slightly prefer Marcus's suggested syntax:

matrix<real[0,1]>[3,4] x;

I don't have any trouble keeping track of the fact
that angle brackets are for type information and square
brackets for size info. But we do have to do something
so we don't wind up with the too-confusing notation:

matrix[0,1][3,4] x;

To allow matrices of int values, we'd probably use something
like:

matrix[scalar=int, bounds=1,5][3,4] n;

I don't like the overload of the comma here, but I can already
hear Andrew's voice saying he'll never be able to keep
the use of semicolons straight if we did something like:

matrix[scalar=int; bounds=1,5][3,4] n;

and that it's too much notation with more cleanly
delineated subexpressions as in:

matrix[scalar=int; bounds=[1,5]][3,4] n;

- Bob

[Andrew Gelman]

We could do real[lower=0,upper=1] x;
I just thought that real[bounds=0,1] x;
would be clearer to users. But I agree that there's somethign weird about that comma there.
A

> --
>

[Bob Carpenter]

On 8/15/12 5:10 PM, Andrew Gelman wrote:
> We could do real[lower=0,upper=1] x;
> I just thought that real[bounds=0,1] x;
> would be clearer to users. But I agree that there's somethign weird about that comma there.

Either way's easy to do from this end. I don't
have much preference. If we did [lower=0, upper=1],
people would probably expect us to allow the
order [upper=1, lower=0], too. This is a bit of
a pain, but not that much.

I sort of like the bounds=0,1. But it might get
confusing with the first two of:

real[bounds=,1] real[upper=1]
real[bounds=0,] real[lower=0]
real[bounds=0,1] real[lower=0,upper=1]

but I find

real[bounds=[,1]]
real[bounds=[0,]]
real[bounds=[0,1]]

- Bob

[Andrew Gelman]

B
Your last version is bracket-heavy but maybe ok. Stan will just be known as "the language with all the brackets."
A

> --
>