Creating a linear predictor

Modeling
1

I have a linear predictor for a ragged array: subjects have two kinds of response, but varying numbers of each. For each response there are a set of predictors.
data {
int<lower=0> N; // number of observations
int<lower=0> M; //number of persons (subjects)
int<lower=0> P; // number of variables
real y[N]; // observed values
int outcome[N]; // which outcome for each y
real x[N,P]; // covariates
}
parameters {
matrix[P,2] beta;
}

Below is the code I’d like, which will allow for general P. What I currently have is a separate .stan file for P=2, 3, and 4 with the inner loop written out. I assume there is a way to do this, which will be obvious once someone points it out.

transformed parameters {
vector [N] lp; // linear predictor

for (i in 1:N) {
lp[i] =0;
for (j in 1:P) { lp[i] = lp[i] + beta[j, outcome[i]] * X[i,j];}
}
}

I’ve left off the remainder of the code, which includes a subject/response specific random intercept, prior for the parameters, etc.

Terry T.

2

That looks fine to me, though I don’t see where it is ragged.

Said another way, it looks like you have two groups, and so two sets of betas, one for each group. Is that correct?

If you declared x like:

matrix[N, P] x;

Then you could write:

transformed parameters {
  vector [N] lp; // linear predictor

  for(i in 1:N) {
    lp[i] = x[i,] * beta[,outcome[i]];
  }
}

for short.

If you surround text with three tick marks on each side you get a code block which formats things nicely (that’s the preformatted text button in the little post toolbar).

3

That is obvious once written down (but wasn’t obvious to me before).
Thank you very much.

Terry T.

4

A follow-up to my follow-up. Stan now complains that “the expression is ill formed” with respect to the line lp[i] = x[i,] * beta[,outcome[i]];

x is an [N,P] matrix and beta is a [P,2] matrix. x[i,] is a row vector and beta[P,2] is a column vector, mathematically, so it should be unambiguous, but STAN is unhappy. How do I convince it of my good intentions?

5

And a last: I was wrong: it was real x[N,P] in the declaration. A matrix of values is not a matrix I guess. The code now works.

1 Like
6

Whoops, yeah, Stan is picky about what is what. There are functions to convert around between types.