map_rect functions aren't able to make use of vectorized forms

[seantalts]

Summary:

map_rect needs a vectorized version, essentially.

Description:

Right now, map_rect specifies that it takes in a function that operates on a single row of a dataset, even if a node is executing that code on a much larger chunk of data at once. We can find speed improvements if we create a version that specifies a user defined function that operates on a chunk of data instead of a single row.

Expected Output:

A new map_rect type signature where the UDF takes in a vector and likely has controllable chunk sizes or chunking broadly.

Tagging @wds15 and @yizhang-cae - either of you have a semi-realistic but fake data set and model that is useful for something like this? I guess it'd be anything with a ton of observations on a vectorized lpdf.

Current Version:

v2.18.0

[wds15]

As discussed, I think this should give us more performance... intuitively. Before we embark on this we should probably first do a simple test run where we pack the data into a suitable form for the current map_rect (and unpack in the UDF). With this test we can quickly find out if chunking will help us (it really should).

(More utilities for packing and unpacking may give us a better ROI overall)

[seantalts]

Agreed. I was asking for fake data and model for that purpose but I can just make something simple for a POC.

[bgoodri]

Perhaps when we have tuples, then map_rect can take a std::vector of tuples and loop over them in parallel

…

On Thu, Jul 26, 2018 at 4:12 PM seantalts ***@***.***> wrote: Yep. I was asking for fake data and model for it but I can just make something simple for a POC. — You are receiving this because you are subscribed to this thread. Reply to this email directly, view it on GitHub <#955 (comment)>, or mute the thread <https://github.com/notifications/unsubscribe-auth/ADOrqkp4rZ4T1fQeTLe9WOQsdQhM8XXaks5uKiLDgaJpZM4ViMWG> .

[bob-carpenter]

This issue isn't concrete enough to act on. It needs

• signature for function being proposed
• semantics

It's not even clear to me what the point of the issue is other than a hint that it'll lead to speed improvements. The "Mapped Logistic Regression" section of the user manual uses vectorized bernoulli logit in each shard. Just so nobody has to look too far afield, here is its:

functions {
  vector lr(vector beta, vector theta, real[] x, int[] y) {
    real lp = bernoulli_logit_lpmf(y | beta[1] + to_vector(x) * beta[2]);
    return [lp]';
  }
} 
data {
  int y[12];
  real x[12]; 
}
transformed data {
  int ys[3, 4] = { y[1:4], y[5:8], y[9:12] };
  real xs[3, 4] = { x[1:4], x[5:8], x[9:12] };
  vector[0] theta[3];
}
parameters {
  vector[2] beta;
}
model {
  beta ~ std_normal();
  target += sum(map_rect(lr, beta, theta, xs, ys));
}

There are three shards, each of size four, and each using a vectorized probability distribution.

[seantalts]

@wds15 I think had a specific example where this was harder - I can imagine working this out with a vectorized multinormal to be even more annoying, for example.

There is another interesting thing here that I just realized - I've always been trained to write an issue that solely describes a problem and avoids prescribing a solution. I guess that's exactly the opposite of what we ask for on Stan. Woops.

[bob-carpenter]

I agree that the packing and unpacking will be painful, but the alternative is anticipating what shapes of things people need. Arrays of tuples might be good for that if we could have something like a generic interface. Writing an issue up that's a bug that needs to be solved is fine. The functional spec is that it should work the way it is documented. Writing up a feature request is trickier, as it's often not clear what we want. I think they need to be at least at the level of functional specs in the sense of what would qualify as a solution. What's happening now is that we are running into problems where volunteers are picking up these issues and trying to solve them and then they submit PRs that don't match what we had in mind. That's why I want the issues to be more precise.

[bgoodri]

Can someone remind me again what is stopping us from being able to write generalized lambda functions in the model block that are not tied to a particular type? Something like

model {
  function lr = [&beta](theta, x, y) {
    return [bernoulli_logit_lpmf(y | beta[1] + x * beta[2])]';
  }
  target += sum(map_rect(lr, theta, xs, ys));
}

where theta, xs, and ys are anything you can iterate over rather than specifically vector [], real [,], and int[,].

[wds15]

This looks like a state-full lambda construct which we don't support for the moment being with map_rect. The reason is that we use only the type to encode what functor get's called on the remote.

However, what you have coded up can be adapted to what we have. You would need to cast beta to be a shared parameter vector. The theta, x & y need to evaluate for each element to an equally sized container.

[bgoodri]

I realize we can't do it right now, as is, but I thought there was no real obstacle to relaxing the constraints that the functors in Stan impose. On that note, why do theta, x, and y fundamentally need to anything stricter than iterable over their first dimension and having the same number of elements on their first dimension?

…

On Sat, Oct 27, 2018 at 2:11 PM wds15 ***@***.***> wrote: This looks like a state-full lambda construct which we don't support for the moment being with map_rect. The reason is that we use only the type to encode what functor get's called on the remote. However, what you have coded up can be adapted to what we have. You would need to cast beta to be a shared parameter vector. The theta, x & y need to evaluate for each element to an equally sized container. — You are receiving this because you commented. Reply to this email directly, view it on GitHub <#955 (comment)>, or mute the thread <https://github.com/notifications/unsubscribe-auth/ADOrqlXFBERLh7tD-tHDLEZW4QohsTx3ks5upKHmgaJpZM4ViMWG> .

[bob-carpenter]

Mainly because it's hard to code in the parser. It'll also be a bit tricky on the C++ math side, but I think that's doable. This is one of the main motivations to move to a new language for the parser and code generator. So we probably won't get to lambdas (closures or not) for six or twelve months. But they're very high priority when we do get a new parser. Polymorphic (templated) lambdas got introduced in C++14 with auto and I just tested that they got into C++1y. For MPI, ODE solvers, etc., I think we'll be able to relax some of the need to pack everything into vectors with parameter-pack type code on the C++ side. But we're still going to need to call out data vs. parameters in MPI (need to know data to send over and derivatives to send back) and ODEs (need to know which Jacobian is necessary). We can also implement other functions on top of map_rect that do their own index fiddling or packing. Suggestions welcome.

[seantalts]

+1 this is a high priority feature! I think in this case we'd need to have stanc actually dissect the lambda and generate a templated or data-specific actual non-lambda function for this example to work, precisely because as Sebastian noted this lambda is capturing state with &beta and I don't think C++ has a way of sending capturing lambdas over MPI.

[bgoodri]

Ah, OK. I had threads in mind and am less familiar about what is needed for MPI.

…

On Wed, Oct 31, 2018 at 8:22 AM seantalts ***@***.***> wrote: +1 this is a high priority feature! I think in this case we'd need to have stanc actually dissect the lambda and generate a templated or data-specific actual non-lambda function for this example to work, precisely because as Sebastian noted this lambda is capturing state with &beta and I don't think C++ has a way of sending capturing lambdas over MPI. — You are receiving this because you commented. Reply to this email directly, view it on GitHub <#955 (comment)>, or mute the thread <https://github.com/notifications/unsubscribe-auth/ADOrqorIITts7kcz2FLWc7lQMtt_58bVks5uqZYCgaJpZM4ViMWG> .