Parametrised integrals

General
1

Is it possible to define a function in stan that returns a parametrized integral? For example, an integrand parametrized by t_i,

I_i = \int_{a}^{b} f(s,t_{i}) \,\mathrm{d}s \,.

As far as I understand this is not possible, since the integrand has a strict definition

real integrand(real x, real xc, real[] theta,real[] x_r, int[] x_i)

that does not accept parametrized functions. Or the parametrization of the integrand should be fixed and dependent only on x_r and x_i.

My questions are:

  • Is my understanding correct?
  • Is there a workaround? I guess one workaround is to solve the integral using an ode solver.
  • Is there a particular reason the integrand has a strict definition? Since in the definition of the ode systems the user is free to pass any parameters.

Thank you for your time and your assistance.

2

That is what the real[] theta argument to the integrand is for.

3

Isn’t this variable only meant for variables of the block parameters?

It is not clear from my post, but in my case, t_i is not a parameter to be inferred, it is one of the data block.

4

If t is a real number, then pass it to the integrand function via x_r. If t is an integer, then use x_i.

5

I don’t think this can work since the integral is inside a loop and depends on t_i. I can pass the whole vector t through x_r but not the specific value t_i that I need on the i-th iteration.

6

Can’t you pass i through x_i, pass the whole vector \mathbf{t} through x_r, and pick out real t_i = x_r[x_i[1]]?

7

I tried this. I tried to change the value of x_i[1] inside the loop, i.e., x_i[1]=i and I got an error.

Is it allowed to change the value of transformed data inside the model block? If so, maybe I misinterpreted the error.

8

Integer arrays should always behave like data (the relevant distinction is not the block it’s defined in but whether the sampler needs to calculate a derivative). Something like this

real integrand(real x, real xc, real[] theta, real[] x_r, int[] x_i) {
  real t_i = x_r[x_i[1]];
  return f(x, ti);
}
...
real x_r[N];
...
real theta[0];
for (i in 1:N) {
  ... = integrate_1d(integrand, a, b, theta, x_r, { i });
}

Note the curly braces that create a new array. Actually, this also works

integrate_1d(integrand, a, b, theta, {x_r[i]}, {1})

because anything built out of data-only expressions is still “data-like”.

1 Like
9

Hi Niko. Thanks for the reply. I tried it and I get the error:

   -------------------------------------------------
    51:        s[i] += 10.^gn * (x2[i]-1/x[i]);
    52:
    53:        real q = integrate_1d( J1_I, 0, pi(), {b,lw,th}, {x[i]}, {1} );
                        ^
    54:
    55:        s[i] += cm * q;
   -------------------------------------------------

Ill-typed arguments supplied to function 'integrate_1d'. Available signatures:
((real, real, real[], data real[], data int[]) => real, real, real, real[], data real[], data int[]) => real
((real, real, real[], data real[], data int[]) => real, real, real, real[], data real[], data int[], data real) => real
Instead supplied arguments of incompatible type: (real, real, real[], real[], int[]) => real, int, real, real[], real[], int[].

Same error with the integrate_1d(integrand, a, b, theta, x_r, { i }); approach.

From the error message, I understand that x_r a x_i must be data.

10

That’s weird because this compiles for me.

functions {
  real integrand(real x, real xc, real[] theta, real[] x_r, int[] x_i) {
    real t_i = x_r[x_i[1]];
    return x + t_i;
  }
}
data {
  int N;
  real x_r[N];
  int x_i[0];
}
parameters {
  real b;
  real lw;
  real th;
}
model {
  for (i in 1 : N) {
    target += integrate_1d(integrand, 0, pi(), {b, lw, th}, {x_r[i]}, {1});
  }
}

Is the integrate_1d call inside a function? If so you must mark x as data in the signature.

Anyway, theta argument does not have to be from the parameters block.

integrate_1d( J1_I, 0, pi(), {b,lw,th,x[i]}, {1.0}, {1} )
11

Thank you Niko. This worked.

Yes, the integrate_1d is called inside another function and I have to add data in front of all the variables I entered in place of x_r in the integrate_1d call.