Finding maximum likelihood

Question

This is a model that is used to model soccer scores, so $i$ and $j$ are, respectively, home and away teams. Random variables $(x,y)$ are the goals scored by the home and away teams, respectively. Parameter $\lambda$ is a known mean goals scored by the home team and $\mu$ is the mean goals scored by the away team. I have managed to fix all the other parameters except for $\rho$, which I have to estimate via MLE.

$$Pr(X_{i,j}=x, Y_{i,j}=y)=\tau_{\lambda, \mu}(x,y)\frac{\lambda^x \text{exp}(-\lambda)}{x!}\frac{\mu^y\text{exp}(-\mu)}{y!}$$ where $$\lambda=\alpha_{i}\beta_{j}\gamma$$ $$\mu=\alpha_{j}\beta_{i}$$ and $$\tau_{\lambda,\mu}(x,y)=\left\{\begin{array}{cc} 1-\lambda\mu\rho &\text{if $x=y=0$,} \\ 1+\lambda\rho &\text{if $x=0,y=1$,}\\ 1+\mu\rho &\text{if $x=1,y=0$,}\\ 1-\rho &\text{if $x=y=1$,}\\ 1 &\text{otherwise}\end{array} \right.$$

Based on the above equations, all the parameters $(\lambda, \mu, \alpha, \beta, \gamma)$ are known constants.

So, now, the problem that I am having is that I have no clue on how to estimate $\rho$ using the maximum likelihood function since a piece-wise equation is involved.

Also, it will be great if anyone can do this using R.

Answer 1

If I understand it, your data is $\{(i_m,j_m,x_m,y_m) : m=1,\dots,n\}$, in which, for the $m$-th match, $i_m$ is the index of the home team, $j_m$ is the index of the away team ($i_m\neq j_m$), $x_m$ is the number of goals scored by the home team, and $y_m$ is the number of goals scored by the away team.

If we have $t$ teams disputing the $n$ matches, the likelihood $L(\rho,\gamma,\alpha_i,\beta_i, i=1,\dots,t)$ is proportional to $$ \prod_{m=1}^n \tau_{\lambda_m,\mu_m}(x_m,y_m)\,\lambda_m^{x_m} \,\mu_m^{y_m}\,e^{-(\lambda_m+\mu_m)} \, , $$ in which $\lambda_m = \alpha_{i_m} \beta_{j_m}\gamma$ and $\mu_m = \alpha_{i_m} \beta_{j_m}$. This likelihood is a function of $2(t+1)$ parameters. You will probably need to implement some conjugate gradient method to solve this constrained optimization problem. I wonder what happens if we put some priors on the parameters and MCMC the posterior with a random walk Metropolis. You said that you know the values of the parameters, except for $\rho$ (and I wonder how...). Hence, your first step is to write an R function that computes the product above for each value of $\rho$. Then, plot it and use the optimizer to see what happens.