Hi,
I’m trying to use on a model which uses logistic regression to model the time varying ratings of multiple teams competing in a tournament. For a given match n, we assume that the score between teams follows a binomial logit model i.e. the number of rounds won by a team is dependent on the rating difference between the two teams on a logit scale. For simplicity, I assume that teams ratings can vary with time according to a Brownian process. It would probably be better to use a Gaussian process, but a Brownian process is easier to evaluate.
Unfortunately, Stan has real issues sampling from my model . I assume this is because we are sampling from the ratings directly, but logistic regression only cares about the difference between the ratings (centered parameterization problem). I tried incorporating a “mu” parameter to allow the mean of the ratings to vary, so I could sample from just a Gaussian distribution, but this doesn’t seem to work too well either.
I think the most efficient way to obtain team scores would be to sample from an M by M matrix of the differences between team ratings, with zeros along the diagonal. However, I’m not sure how to constrain a matrix to be antisymmetric in this way in stan, and I’m not sure what kind of prior I would use for this kind of matrix. I do think this might be the most efficient sampling method though.
Here’s my current code. If you have any suggestions on how to avoid this issue of centered parameterization, please let me know.
data {
int<lower=0> N;
int<lower=0> M;
int roundscores[N];
int totalrounds[N];
int bluteamnumber[N];
int redteamnumber[N];
real times[N];
}
parameters {
matrix[M,N] team_scores;
real<lower=0> sd_change;
real <lower=0> sd_teams;
real<lower=-6,upper=6> mu[N];
}
transformed parameters {
vector[N] score_diff;
matrix[M,N] team_scores_star;
for (n in 1:N){
team_scores_star[,n]=team_scores[,n]+mu[n];
score_diff[n] = team_scores_star[bluteamnumber[n]][n]-team_scores_star[redteamnumber[n]][n];
}
}
model{
sd_change~ normal(0,0.1);
sd_teams~ normal(0,2);
for (m in 1:M){
team_scores[m]~normal(0,sd_teams);
}
for (n in 2:N){
team_scores_star[,n]~normal(team_scores_star[,n-1],(times[n]-times[n-1])*sd_change);
}
roundscores~binomial_logit(totalrounds,score_diff);
}