Dear All,

I am trying to fit a random-effects stochastic frontier model (basically same as a regression but with inefficiency component u) with spatial element based on the specification in the link https://rstudio-pubs-static.s3.amazonaws.com/243277_01730c1f0a984132bce5d5d25bec62aa.html

I have 20 units with each of the units being observed 32 times, which gives me a grand total of 640 observations.

I am trying to model the spatial relationship for the standard deviation of the inefficiency component. My spatial weights matrix is 20 x 20, but of course, when using multi_normal distribution it is demanding that my spatial weights matrix be 640 x 640, which is of course not needed as it will be 31 repetitions for each unit.

I am sure that there is a way to code the 20 x 20 matrix into the multi_normal distribution so it knows that unit id ‘1’ in weights matrix relates to 32 ‘1s’ in the data.

Here is my code and the error I am getting says **"multi_normal_lpdf: Size of random variable (640) and rows of covariance parameter (20) must match in size "**

Any help will be greatly appreciated:

```
data {
int<lower=0> N;// Number of observations
int<lower=0> J; // Number of groups
int<lower=0> P; // Number of predictors in the regression
real Y[N]; // dependent variables
matrix[N,P] X;
int<lower=1,upper=J> dhb_id[N]; // specifying 20 groups for random effects
matrix<lower=0>[J,J] W; // 20 x 20 Spatial weight matrix
matrix<lower=0,upper=1>[J,J] I; // 20 x 20 identity matrix
}
parameters {
real alpha;
vector[P] beta;
real<lower=0> sigma_e;
real<lower=0> sigma_u; // SD of variance
real<lower=-1,upper=1> lambda; // measures the spatial correlation
vector<lower = 0>[N]mean_u;
vector<lower = 0>[N]u; // inefficiency component
vector[J]RE; // vector capturing random effects
}
transformed parameters {
vector[N]yhat;
yhat = alpha + X*beta +u + RE[dhb_id] ;
}
model {
sigma_e~exponential(1);
alpha ~ normal(0,1);
beta ~ normal (0,2);
lambda ~ normal (0,2);
sigma_e~exponential(1);
mean_u ~ normal(0,2);
u ~ multi_normal(mean_u, **tcrossprod(I - lambda * W)*(sigma_u*sigma_u)**);
RE ~ normal(0,2);
Y ~ normal(yhat, sigma_e);
}
```

With much gratitude

Antony