Help with rewriting multinomial distribution to dirichlet-multinomial distribution

Hi,

In mixed martial arts a fighter can aim for head, body or legs. I would like to model the fighting style of different fighters.

My first model looks as follows:

data {
    int<lower=1> n_fights;
    int<lower=2> n_fighters;
    int<lower=0> n_strikes_head[n_fights];
    int<lower=0> n_strikes_body[n_fights];
    int<lower=0> n_strikes_legs[n_fights];
    int<lower=0> n_strikes[n_fights];
}


parameters {
    vector<lower=0, upper=1>[n_fights] phi_head;
    vector<lower=0, upper=1>[n_fights] phi_body;
    vector<lower=0, upper=1>[n_fights] phi_legs;
}

model {
    n_strikes_head ~ binomial(n_strikes, phi_head);
    n_strikes_body ~ binomial(n_strikes, phi_body);
    n_strikes_legs ~ binomial(n_strikes, phi_legs);
}

The first problem I came across is that phi_head + phi_body + phi_legs is not equal to 1.
I tried to solve this by incorporating an unit simplex and bundling the strikes together in a single matrix

data {
    int<lower=1> n_fights;
    int<lower=2> n_fighters;
    int<lower=1, upper=n_fighters> fighter_id[n_fights];
    array[n_fights, 3] int<lower=0> n_strikes;                 // n_strikess_head, n_strikes_body, n_strikes_legs
}

parameters {
    vector<lower=0>[3] alpha;
    simplex[3] phi[n_fights];
}

model {
    for (fight in 1:n_fights) {
        phi[fight] ~ dirichlet(alpha);
        n_strikes[fight] ~ multinomial(phi[fight]);
    }

}

Now, I want to give each fighter his own phi / fighting skill.
Is this done in the following code?

data {
    int<lower=1> n_fights;
    int<lower=2> n_fighters;
    int<lower=1, upper=n_fighters> fighter_id[n_fights];
    array[n_fights, 3] int<lower=0> n_strikes;                 // n_strikess_head, n_strikes_body, n_strikes_legs
}


parameters {
    vector<lower=0>[3] alpha;
    simplex[3] phi[n_fighters];
}

model {
    for (fight in 1:n_fights) {
        phi[fighter_id[fight]] ~ dirichlet(alpha);
        n_strikes[fight] ~ multinomial(phi[fighter_id[fight]]);
    }
}

Thank you

You should set the dirichlet distribution once for each fighter rather than for each fight. You also need a prior for alpha; the default (improper uniform) is too wide.

model {
    alpha ~ cauchy(0, 1);
    for (fighter in 1:n_fighters) {
        phi[fighter] ~ dirichlet(alpha);
    }
    for (fight in 1:n_fights) {
        n_strikes[fight] ~ multinomial(phi[fighter_id[fight]]);
    }
}

You might want to try the built-in dirichlet_multinomial distribution: Multivariate Discrete Distributions