<p>Handball has received growing interest during the last years, including academic research for many different aspects of the sport. However, modelling the game’s outcome has attracted less interest mainly because of the additional challenges. For example, it has been observed that the number of goals exhibit variance smaller than the mean, namely they are underdispersed relative to a Poisson distribution. In this paper, we focus on modelling the difference of scores in a handball game. We fit a model directly on the goal difference to avoid the observed underdispersion. In the present paper, we propose some models starting from a Skellam regression model and also considering zero-inflated versions as well as other discrete distributions. Furthermore, we develop some new bivariate copula models to estimate the conditional distribution of the score difference in the second half given the score difference in the first half, thus providing novel insights into the game. Data from the German Bundesliga are used to show the potential of the new models.</p>