A simple matrix completion model for fixed K in Stan

matrix_completion.stan

data {
  int N; // number of individuals
  int T; // number of time periods
  matrix[N, T] Y; // outcome matrix; missing entries set to -9.0
  int K; // rank of matrix
}
parameters {
  matrix[N, K] M; // individual loadings
  matrix[T, K] U; // time factors
  real<lower = 0> sigma; // error scale
  positive_ordered[K] gamma; // vector of individual prior scales on time shocks
  real<lower = 0> gamma_scale; // hyperscale
}
transformed parameters {
  matrix[N, T] Theta = M * U';
}
model {
  sigma ~ inv_gamma(1, 1);
  gamma_scale ~ exponential(1);
  gamma ~ normal(0, gamma_scale);
  
  for(t in 1:T) {
    U[t] ~ normal(0, gamma);
  }
  
  // likelihood
  for(i in 1:N) {
    M[i] ~ normal(0, 1);
    for(t in 1:T){
      if(Y[i,t]!= -9.0) {
        Y[i,t] ~ normal(Theta[i,t], sigma);
      }
    }
  }
}