The use of curly braces in stan functions

I’d like to ask about the use of curly braces in stan.

My understanding is that for example, when we want to write a for loop with more than one arguments in the loop, we’d want to use {}. However, I’ve seen multiple times that these braces get used in the middle of a function, or in the model block for no apparent reason. Below is an example.

vector theta_pred_rng(array[] vector x, array[] vector x_pred, matrix LK,
                        vector theta, real alpha, real l, real epsilon) {
		int n = size(x);
		int n_pred = size(x_pred);
    vector[n_pred] theta_pred;
		{
      matrix[n,n_pred] K12 = gp_matern32_cov(x,x_pred,alpha,l);
      matrix[n,n_pred] LKinv_K12 = mdivide_left_tri_low(LK,K12);
      vector[n_pred] LKinv_theta0 = mdivide_left_tri_low(LK,theta);
      vector[n_pred] theta_pred_mu = LKinv_K12' * LKinv_theta0;
      matrix[n_pred,n_pred] K22 = add_diag(gp_matern32_cov(x_pred,alpha,l), epsilon);
      matrix[n_pred,n_pred] theta_pred_cov = K22 - LKinv_K12'*LKinv_K12;
      theta_pred = multi_normal_rng(theta_pred_mu, theta_pred_cov);
		}
		return theta_pred;
	}

At first, I thought this is only for clarify/ ease to read of the codes. However, later I realized that the braces have actual impacts on the execution of the codes … E.g., I noticed that if instead of declaring theta_pred before the braces, I put vector[n_pred] theta_pred=multi_normal_rng() in the braces, then I wouldn’t be able to return this variable. I’d get an error message saying “Identifier ‘theta_pred’ not in scope.”.

So I wonder, are the braces explicitly creating a local environment? What are the pros and cons of using these braces and when should I use them? I couldn’t find any documentation on this. Would appreciate if someone points me to some references of could explain a bit more.

Thanks!

Yes, these local blocks create their own scope (In more technical terms, Stan is a “lexically block scoped language”).

There are a bunch of reasons someone might want to use these, but in Stan I think the most common uses are in the transformed parameters and generated quantities blocks, for if you want to create a helper variable but don’t want it in the output. A silly example:

generated quantities {
 vector[N] foo; // will be in output
  {
    vector[N] bar; // won't be in output
    // some calculation to create bar
    foo = abs(bar); // or whatever
  }    
}

Note that using { } after for is technically optional - the for loop will just run the next statement, including if the next statement is a local block:


//equivalent

for (i in 1:10)
  print(i);

for (i in 1:10) {
  print(i);
}

We recommend the second form for style reasons

3 Likes

I’ll add, as I learned from personal experience, for GPs if you don’t use the local scoping the Stan output will include the full covariance matrix and eventually crash R after consuming all the available RAM.

1 Like

I noticed the same thing … even if I explicitly ask not to output those variables by using include=T, pars=c().

Now I’m using the generated quantity block to get posterior predictives. With n=100, and predicting for 121 points, I noticed that sampling finishes pretty fast, but after that it’s been taking more than an hour, but still couldn’t finish the run. I don’t understand why it’s taking such a long time …

I have a follow-up question. Say if I don’t use these local blocks in the codes, but when I do sampling, I use include=T, pars=c() to specify which variables I want in the output. Does that have the same effects? Or I’d still be better off to use these local blocks to save memory during the run?

I believe RStan is implemented in such a way that those items are dropped from the output early to save memory (I’m far from an expert in RStan, though!). They still get passed around a bunch, but I don’t think space for all of them is allocated

Using local blocks when possible is probably still a good idea for both clarity and to be more portable to other Stan interfaces, though!

1 Like

This is from munging and analyzing the output. @bgoodri could provide a definitive answer about what’s going on. cmdstanr will finish faster, because it will stream output to the file system (disk or SSD), assuming you’re not bottlenecked on file system writes like in some clusters. But then it will have to read the data into R to do the posterior analyses.

In contrast, CmdStanPy does all this analysis in C++ from disk using standalone executables. That is, unless you use ArviZ to get some more functionality, at which point things start working like cmdstanr.

Thanks for the reply! Here’s what I’ve noticed so far:

  1. same codes, on the same environment, same setup, could sometimes take > 2 hours and still not finished, or can finish within 10 minutes. In both cases, the sampling stage were quickly done
  2. if I delete the generated quantity block, then it seems the run won’t get stuck
  3. if I run everything on the server instead of on my laptop, it seems the run can finish without a problem so far. This made me wonder whether it has anything to do with the memory limit setup of my Rstudio
  4. I thought about using cmdstanr, however, I decided to go with rstan because I want to use include=T, pars=c() to only save outputs for selected parameters. Otherwise my output file will be too large …

Thanks!!