Matrix dimension limit?

I just got an error when compiling a model (with CmdStan) that makes me think there is limit on matrix dimension, namely that matrices must have 8 or fewer dimensions. Is that correct? The actual error was about the “to_vec” function and finding no matches that took 9 arguments, but matches for 8 and 7, etc. And the model compiles and runs once I reduce to an 8-dimensional matrix. I’m interpreting this as I said above but perhaps it originates in something I am trying to do with the matrix?

Is this difficult to change? I can work around it. The matrix in question is a post-stratification weight matrix for a model with several terms. But rather than use a multi-dimensional matrix, I could just use a long vector and indexes for each factor. That would make the json for the matrix simpler, at the cost of needing a lot of those new indexes. Is that the way people usually handle this?

Thanks!

Can you post an example of the code that’s failing?

Sure! It’s…not minimal. Or correctly formatted. I’m generating via a Haskell library so it’s harder to isolate just the bit that fails. I’m compiling using cmdStan 2.26.1, on OS X Big Sur with the built-in clang++.

Thanks!

data {
  int<lower=2> J_State;
  int<lower=1> N;
  int<lower=1> State[N];
  int<lower=2> J_Education;
  int<lower=1> Education[N];
  int<lower=2> J_Race;
  int<lower=1> Race[N];
  int<lower=2> J_Age;
  int<lower=1> Age[N];
  int<lower=2> J_WNH;
  int<lower=1> WNH[N];
  int<lower=2> J_Ethnicity;
  int<lower=1> Ethnicity[N];
  int<lower=2> J_CD;
  int<lower=1> CD[N];
  int<lower=2> J_WhiteNonGrad;
  int<lower=1> WhiteNonGrad[N];
  int<lower=2> J_Sex;
  int<lower=1> Sex[N];
  int<lower=0> T[N];
  int<lower=0> S[N];
  int K_CD;
  matrix[J_CD,K_CD] X_CD;
  int<lower=0> N_WI_ACS_WNH_State;
  real WI_ACS_WNH_State_wgts[N_WI_ACS_WNH_State,J_Education,J_Race,J_Age,J_WNH,J_Ethnicity,J_CD,J_WhiteNonGrad,J_Sex];
  int<lower=0> N_WI_ACS_NWNH_State;
  real WI_ACS_NWNH_State_wgts[N_WI_ACS_NWNH_State,J_Education,J_Race,J_Age,J_WNH,J_Ethnicity,J_CD,J_WhiteNonGrad,J_Sex];
  int<lower=0> N_NI_ACS_WNH_State;
  real NI_ACS_WNH_State_wgts[N_NI_ACS_WNH_State,J_Education,J_Race,J_Age,J_WNH,J_Ethnicity,J_CD,J_WhiteNonGrad,J_Sex];
  int<lower=0> N_NI_ACS_NWNH_State;
  real NI_ACS_NWNH_State_wgts[N_NI_ACS_NWNH_State,J_Education,J_Race,J_Age,J_WNH,J_Ethnicity,J_CD,J_WhiteNonGrad,J_Sex];
}
transformed data {
  vector[K_CD] mean_X_CD;
  matrix[J_CD,K_CD] centered_X_CD;
  for (k in 1:K_CD) {
    mean_X_CD[k] = mean(X_CD[,k]);
    centered_X_CD[,k] = X_CD[,k] - mean_X_CD[k];
  }
  matrix[J_CD,K_CD] Q_CD_ast = qr_thin_Q(centered_X_CD) * sqrt(J_CD - 1);
  matrix[K_CD,K_CD] R_CD_ast = qr_thin_R(centered_X_CD) / sqrt(J_CD - 1);
  matrix[K_CD,K_CD] R_CD_ast_inverse = inverse(R_CD_ast);
}
parameters {
  real alpha;
  vector[K_CD] thetaX_CD;
  real eps_Sex;
  real eps_WNH;
  real<lower=0> sigma_Race;
  vector[J_Race] beta_Race_raw;
  real eps_Ethnicity;
  real eps_Age;
  real eps_Education;
  real eps_WhiteNonGrad;
  real<lower=0> sigma_State;
  vector[J_State] beta_State_raw;
  real<lower=0> sigma_WNH_State;
  vector[J_State] eps_WNH_State_raw;
}
transformed parameters {
  vector[K_CD] betaX_CD;
  betaX_CD = R_CD_ast_inverse * thetaX_CD;
  vector[J_Race] beta_Race = sigma_Race * beta_Race_raw;
  vector[J_State] beta_State = sigma_State * beta_State_raw;
  vector[J_State] eps_WNH_State = sigma_WNH_State * eps_WNH_State_raw;
  vector[N] y_WNH_State;
  for (n in 1:N) {
    y_WNH_State[n] = {eps_WNH_State[State[n]],-eps_WNH_State[State[n]]}[WNH[n]];
  }
}
model {
  alpha ~ normal(0,2);
  thetaX_CD ~ normal(0,2);
  eps_Sex ~ normal(0,2);
  eps_WNH ~ normal(0,2);
  sigma_Race ~ normal(0,2);
  beta_Race_raw ~ normal(0,1);
  eps_Ethnicity ~ normal(0,2);
  eps_Age ~ normal(0,2);
  eps_Education ~ normal(0,2);
  eps_WhiteNonGrad ~ normal(0,2);
  sigma_State ~ normal(0,2);
  beta_State_raw ~ normal(0,1);
  sigma_WNH_State ~ normal(0,2);
  eps_WNH_State_raw ~ normal(0,1);
  S ~ binomial_logit(T,alpha + Q_CD_ast[CD] * thetaX_CD + to_vector({eps_Sex,-eps_Sex}[Sex]) + beta_Race[Race] + to_vector({eps_Ethnicity,-eps_Ethnicity}[Ethnicity]) + to_vector({eps_Age,-eps_Age}[Age]) + to_vector({eps_Education,-eps_Education}[Education]) + to_vector({eps_WhiteNonGrad,-eps_WhiteNonGrad}[WhiteNonGrad]) + beta_State[State] + y_WNH_State);
}
generated quantities {
  vector[N] log_lik;
  for (n in 1:N) {
    log_lik[n] = binomial_logit_lpmf(S[n] | T[n],alpha + Q_CD_ast[CD[n]] * thetaX_CD + {eps_Sex,-eps_Sex}[Sex[n]] + beta_Race[Race[n]] + {eps_Ethnicity,-eps_Ethnicity}[Ethnicity[n]] + {eps_Age,-eps_Age}[Age[n]] + {eps_Education,-eps_Education}[Education[n]] + {eps_WhiteNonGrad,-eps_WhiteNonGrad}[WhiteNonGrad[n]] + beta_State[State[n]] + {eps_WNH_State[State[n]],-eps_WNH_State[State[n]]}[WNH[n]]);
  }
  vector[N_WI_ACS_WNH_State] WI_ACS_WNH_State = rep_vector(0,N_WI_ACS_WNH_State);
  for (n in 1:N_WI_ACS_WNH_State) {
    real WI_ACS_WNH_State_WgtSum = 0;
    for (n_Education in 1:J_Education) {
      for (n_Race in 1:J_Race) {
        for (n_Age in 1:J_Age) {
          for (n_WNH in 1:J_WNH) {
            for (n_Ethnicity in 1:J_Ethnicity) {
              for (n_CD in 1:J_CD) {
                for (n_WhiteNonGrad in 1:J_WhiteNonGrad) {
                  for (n_Sex in 1:J_Sex) {
                    real pWI_ACS_WNH_State = inv_logit(alpha + Q_CD_ast[n_CD] * thetaX_CD + {eps_Sex,-eps_Sex}[n_Sex] + beta_Race[n_Race] + {eps_Ethnicity,-eps_Ethnicity}[n_Ethnicity] + {eps_Age,-eps_Age}[n_Age] + {eps_Education,-eps_Education}[n_Education] + {eps_WhiteNonGrad,-eps_WhiteNonGrad}[n_WhiteNonGrad] + beta_State[n] + {eps_WNH_State[n],-eps_WNH_State[n]}[n_WNH]);
                    WI_ACS_WNH_State_WgtSum += WI_ACS_WNH_State_wgts[n][n_Education, n_Race, n_Age, n_WNH, n_Ethnicity, n_CD, n_WhiteNonGrad, n_Sex];
                    WI_ACS_WNH_State[n] += pWI_ACS_WNH_State * WI_ACS_WNH_State_wgts[n][n_Education, n_Race, n_Age, n_WNH, n_Ethnicity, n_CD, n_WhiteNonGrad, n_Sex];
                  }
                }
              }
            }
          }
        }
      }
    }
    WI_ACS_WNH_State[n] /= WI_ACS_WNH_State_WgtSum;
  }
  vector[N_WI_ACS_NWNH_State] WI_ACS_NWNH_State = rep_vector(0,N_WI_ACS_NWNH_State);
  for (n in 1:N_WI_ACS_NWNH_State) {
    real WI_ACS_NWNH_State_WgtSum = 0;
    for (n_Education in 1:J_Education) {
      for (n_Race in 1:J_Race) {
        for (n_Age in 1:J_Age) {
          for (n_WNH in 1:J_WNH) {
            for (n_Ethnicity in 1:J_Ethnicity) {
              for (n_CD in 1:J_CD) {
                for (n_WhiteNonGrad in 1:J_WhiteNonGrad) {
                  for (n_Sex in 1:J_Sex) {
                    real pWI_ACS_NWNH_State = inv_logit(alpha + Q_CD_ast[n_CD] * thetaX_CD + {eps_Sex,-eps_Sex}[n_Sex] + beta_Race[n_Race] + {eps_Ethnicity,-eps_Ethnicity}[n_Ethnicity] + {eps_Age,-eps_Age}[n_Age] + {eps_Education,-eps_Education}[n_Education] + {eps_WhiteNonGrad,-eps_WhiteNonGrad}[n_WhiteNonGrad] + beta_State[n] + {eps_WNH_State[n],-eps_WNH_State[n]}[n_WNH]);
                    WI_ACS_NWNH_State_WgtSum += WI_ACS_NWNH_State_wgts[n][n_Education, n_Race, n_Age, n_WNH, n_Ethnicity, n_CD, n_WhiteNonGrad, n_Sex];
                    WI_ACS_NWNH_State[n] += pWI_ACS_NWNH_State * WI_ACS_NWNH_State_wgts[n][n_Education, n_Race, n_Age, n_WNH, n_Ethnicity, n_CD, n_WhiteNonGrad, n_Sex];
                  }
                }
              }
            }
          }
        }
      }
    }
    WI_ACS_NWNH_State[n] /= WI_ACS_NWNH_State_WgtSum;
  }
  vector[N_NI_ACS_WNH_State] NI_ACS_WNH_State = rep_vector(0,N_NI_ACS_WNH_State);
  for (n in 1:N_NI_ACS_WNH_State) {
    real NI_ACS_WNH_State_WgtSum = 0;
    for (n_Education in 1:J_Education) {
      for (n_Race in 1:J_Race) {
        for (n_Age in 1:J_Age) {
          for (n_WNH in 1:J_WNH) {
            for (n_Ethnicity in 1:J_Ethnicity) {
              for (n_CD in 1:J_CD) {
                for (n_WhiteNonGrad in 1:J_WhiteNonGrad) {
                  for (n_Sex in 1:J_Sex) {
                    real pNI_ACS_WNH_State = inv_logit(alpha + Q_CD_ast[n_CD] * thetaX_CD + {eps_Sex,-eps_Sex}[n_Sex] + beta_Race[n_Race] + {eps_Ethnicity,-eps_Ethnicity}[n_Ethnicity] + {eps_Age,-eps_Age}[n_Age] + {eps_Education,-eps_Education}[n_Education] + {eps_WhiteNonGrad,-eps_WhiteNonGrad}[n_WhiteNonGrad] + beta_State[n]);
                    NI_ACS_WNH_State_WgtSum += NI_ACS_WNH_State_wgts[n][n_Education, n_Race, n_Age, n_WNH, n_Ethnicity, n_CD, n_WhiteNonGrad, n_Sex];
                    NI_ACS_WNH_State[n] += pNI_ACS_WNH_State * NI_ACS_WNH_State_wgts[n][n_Education, n_Race, n_Age, n_WNH, n_Ethnicity, n_CD, n_WhiteNonGrad, n_Sex];
                  }
                }
              }
            }
          }
        }
      }
    }
    NI_ACS_WNH_State[n] /= NI_ACS_WNH_State_WgtSum;
  }
  vector[N_NI_ACS_NWNH_State] NI_ACS_NWNH_State = rep_vector(0,N_NI_ACS_NWNH_State);
  for (n in 1:N_NI_ACS_NWNH_State) {
    real NI_ACS_NWNH_State_WgtSum = 0;
    for (n_Education in 1:J_Education) {
      for (n_Race in 1:J_Race) {
        for (n_Age in 1:J_Age) {
          for (n_WNH in 1:J_WNH) {
            for (n_Ethnicity in 1:J_Ethnicity) {
              for (n_CD in 1:J_CD) {
                for (n_WhiteNonGrad in 1:J_WhiteNonGrad) {
                  for (n_Sex in 1:J_Sex) {
                    real pNI_ACS_NWNH_State = inv_logit(alpha + Q_CD_ast[n_CD] * thetaX_CD + {eps_Sex,-eps_Sex}[n_Sex] + beta_Race[n_Race] + {eps_Ethnicity,-eps_Ethnicity}[n_Ethnicity] + {eps_Age,-eps_Age}[n_Age] + {eps_Education,-eps_Education}[n_Education] + {eps_WhiteNonGrad,-eps_WhiteNonGrad}[n_WhiteNonGrad] + beta_State[n]);
                    NI_ACS_NWNH_State_WgtSum += NI_ACS_NWNH_State_wgts[n][n_Education, n_Race, n_Age, n_WNH, n_Ethnicity, n_CD, n_WhiteNonGrad, n_Sex];
                    NI_ACS_NWNH_State[n] += pNI_ACS_NWNH_State * NI_ACS_NWNH_State_wgts[n][n_Education, n_Race, n_Age, n_WNH, n_Ethnicity, n_CD, n_WhiteNonGrad, n_Sex];
                  }
                }
              }
            }
          }
        }
      }
    }
    NI_ACS_NWNH_State[n] /= NI_ACS_NWNH_State_WgtSum;
  }
  vector[N_WI_ACS_WNH_State] rtDiffWI = WI_ACS_WNH_State - WI_ACS_NWNH_State;
  vector[N_WI_ACS_WNH_State] rtDiffNI = NI_ACS_WNH_State - NI_ACS_NWNH_State;
  vector[N_WI_ACS_WNH_State] rtDiffI = rtDiffWI - rtDiffNI;
}

icic so yeah in the Stan compiler right now we limit it to 8 I’m pretty sure (fyi what you have above is an array of 8 dimensions). Are you getting an error when the stanc compiler compiles this or when the C++ compiles?

I’m pretty sure 8 was just an arbitrary number in the compiler. To lookup whether a function the user calls is accepted by the C++ we generate a big table of all the supported function signatures and I’m pretty sure the limit there is 8 just because the table gets to be pretty big after that. But if there’s a need for it I think we could increase that to 9 or 10.

I understand another compiler is generting this Stan code, though it may be good to talk to them as I think a more efficient form of this would just be a matrix with N_* rows and J_* columns

Thanks! I wrote the Haskell → Stan compiler–not nearly a complete one at this point, but helpful for managing the data and model bits all together. Anyway, I can re-write it to do as you suggest. Do you mean an index matrix, one with a row for every post-stratification cell and columns for the value of each group for that cell? So, in my case N x 8? Then a model term, e.g., for the State intercept, would look something like beta_State[PS_Indexes[n, Col_State]] ?

Is there some reason that would preferable to just generating a 1-D index array for each group the way I do for the modeled data? Thus leading to something more like beta_State[PS_State[n]] ?
I’d have to add each index to the data block but I wouldn’t need variables like Col_State to specify which column of the index matrix corresponds to each group.