Predicting the scoring time in hockey

TL;DR

This paper introduces a Bayesian density estimator based on a weighted beta prime distribution to predict the time until the r-th goal in hockey, leveraging past performance data and expert opinions.

Contribution

It proposes a novel Bayesian predictive density estimator for hockey scoring times using a gamma model and weighted beta prime distribution, outperforming existing methods.

Findings

Estimator outperforms existing models in prediction accuracy.

Efficient in both historical and current NHL seasons.

Uses ancillary information like past performance and expert opinions.

Abstract

In this paper, we propose a Bayesian predictive density estimator to predict the time until the r-th goal is scored in a hockey game, using ancillary information such as their performances in the past, points and specialists' opinions. To be more specific, we consider a gamma distribution as a waiting scoring model. The proposed density estimator belongs to an interesting new version of weighted beta prime distribution and outperforms the other estimator in the literature. The efficiency of our estimator is evaluated using frequentist risk along with measuring the prediction error from the old dataset, 2016-17, to the current season (2018-19) of the National Hockey League.

Figures

Tables

Table 1. Table 2: Bayesian predictive density estimators q ^ 1 subscript ^ 𝑞 1 \hat{q}_{1} , q ^ 0 subscript ^ 𝑞 0 \hat{q}_{0} , with and without ancillary information respectively, along with the truncated probability density function on ( 0 , 60 ) 0 60 (0,60) minutes as well as the mode, mean, and P i subscript 𝑃 𝑖 P_{i} , i × 100 𝑖 100 i\times 100 percentiles.

Estimator	Truncated Density on $(0,60)$	Mode,  Mean,  $P_{0.20}$,   $P_{0.50}$,   $P_{0.90}$
${\hat{q}}_0$	$1901470y_1^2/{\left(35.85+y_1\right)}^6$	17.92  28.35, 14.38, 26.62, 50.3
${\hat{q}}_1$	$y_1^2\left(0.055+\left(0.0004+{10}^{-6}{y}_1\right)y_1\right)/{\left(1.92+0.025y_1\right)}^8$	28.13,  33.12 19.06, 32.82, 53.48

Table 2. Table 1: Time elapsed (in minutes) until scoring the third goal in games corresponding to Toronto Maple Leafs (left) and Montreal Canadiens (right).

Toronto Maple Leafs v.s. Opponents	Time Elapsed
Winnipeg Jets	18.38
New York Rangers	11.12
Chicago Blackhawks	55.68
New Jersey Devils	53.57
Montreal Canadiens	32.67
Detroit Red Wings	15.75
Ottawa Senators	52.85
Los Angeles Kings	42.88
Carolina Hurricanes	27.17
Anaheim Ducks	58.48
St. Louis Blues	50.13
Vegas Golden Knights	15.03
Minnesota Wild	43.65
Boston Bruins	46.83
Montreal Canadiens	40.4
Carolina Hurricanes	31.6
Calgary Flames	41.88
Edmonton Oilers	13.08
Pittsburgh Penguins	12.9
Carolina Hurricanes	10.53
Buffalo Sabres	27.63
New York Rangers	31.35
Arizona Coyotes	11.4
Montreal Canadiens	57.17
Buffalo Sabres	57.93
Pittsburgh Penguins	31.57
Colorado Avalanche	58.07
Vegas Golden Knights	40.43
Dallas Stars	45.18
Buffalo Sabres	26.35
Montreal Canadiens	39.5
Ottawa Senators	52.43
Tampa Bay Lightning	35.67
Ottawa Senators	44.33
Dallas Stars	29.45
New York Islanders	30.52
New York Rangers	20.83
Anaheim Ducks	35.43
Ottawa Senators	11.48
Tampa Bay Lightning	31.57
Columbus Blue Jackets	28.05
Pittsburgh Penguins	35.73
Detroit Red Wings	59.48
New York Islanders	56.48
Boston Bruins	39.05
Tampa Bay Lightning	45.42
Detroit Red Wings	47.43
Florida Panthers	13.9
New York Islanders	37.28
Nashville Predators	36.72

Equations

$$p_{\lambda }(x)$$

$$L_{KL}(q_{\lambda },{\hat{q}}_1(\cdot ))={\int }_{\mathbb{R}}{q}_{\lambda }(y)\log \frac{q_{\lambda }(y)}{{\hat{q}}_1(y;x)}dy,$$

$$R_{KL}(\lambda ,{\hat{q}}_1)={\int }_{\mathbb{R}}{L}_{KL}\left(q_{\lambda }(y),{\hat{q}}_1(\cdot )\right)p_{\lambda }(x)dx.$$

$${\hat{q}}_{{\pi }_A}(y;x)={\int }_{\Lambda }{q}_{\lambda }(y)\pi (\lambda \mid x)d\lambda .$$

$$f_{a,b}(t)$$

$$F_{a,b}(t)$$

$$m(s_1,s_2)={\int }_0^{\infty }\pi (\xi ){\mathrm{IG}}_{s_1,s_2}(\xi )d\xi .$$

$$m(s_1,s_2)=\frac{s_1}{s_2}.$$

$$m(s_1,s_2;\bar{\xi })={\int }_0^{\bar{\xi }}\pi (\xi ){\mathrm{IG}}_{s_1,s_2}(\xi )d\xi ,$$

$$m(s_1,s_2;\bar{\xi })=\frac{\Gamma (s_1+1,\frac{s_2}{\bar{\xi }})}{s_2\Gamma (s_1)}.$$

$$C(k_1,k_2,s_1,s_2)={\int }_0^{\infty }\pi ({\xi }_1)m(s_1,s_2;{\xi }_1){\mathrm{IG}}_{k_1,k_2}({\xi }_1)d{\xi }_1.$$

$$\frac{k_1{k}_2^{k_1}{s}_2^{-k_1-2}\Gamma (k_1+s_1+2){}_2{\tilde{F}}_1\left(k_1+1,k_1+s_1+2;k_1+2;-\frac{k_2}{s_2}\right)}{\Gamma (s_1)},$$

$$\frac{1}{B(a,b)}\frac{1}{\sigma }\frac{(t/\sigma {)}^{a\gamma -1}}{{\left(1+(t/\sigma {)}^{\gamma }\right)}^{a+b}},$$

$${\hat{q}}_0(y_1;x_1)=\frac{1}{B(r_1-1,r_1^{\prime })}\frac{m(r_1+r_1^{\prime }-1,x_1+y_1)}{m(r_1-1,x_1)}\frac{x_1^{r_1-1}{y}_1^{r_1^{\prime }-1}}{(x_1+y_1{)}^{r_1+r_1^{\prime }-1}},y_1>0.$$

$$\pi ({\lambda }_1\mid x_1)=\frac{\pi ({\lambda }_1){\lambda }_1^{-r_1}{x}_1^{-r_1-1}{e}^{-x_1/{\lambda }_1}}{\Gamma (r_1-1)m(r_1-1,x_1)}.$$

$${\hat{q}}_0(y_1;x_1)$$

$$=\frac{1}{\Gamma (r_1^{\prime })\Gamma (r_1-1)m(r_1-1,x_1)}{\int }_0^{\infty }\pi ({\lambda }_1){\lambda }_1^{-(r_1+r_1^{\prime }-1)-1}{x}_1^{r_1-1}{y}_1^{r_1^{\prime }-1}{e}^{-(x_1+y_1)/{\lambda }_1}d{\lambda }_1$$

$$=\frac{1}{B(r_1-1,r_1^{\prime })}\frac{m(r_1+r_1^{\prime }-1,x_1+y_1)}{m(r_1-1,x_1)}\frac{x_1^{r_1-1}{y}_1^{r_1^{\prime }-1}}{(x_1+y_1{)}^{r_1+r_1^{\prime }-1}},$$

$$\frac{1}{B(r_1,r_1^{\prime })}\frac{x_1^{r_1}{y}_1^{r_1^{\prime }-1}}{(x_1+y_1{)}^{r_1+r_1^{\prime }}}{y}_1>0,$$

$${\hat{q}}_1(y_1;x_1,x_2)=\frac{C(r_1+r_1^{\prime }-1,x_1+y_1,r_2-1,x_2)}{C(r_1-1,x_1,r_2-1,x_2)}\frac{m(r_1-1,x_1)}{m(r_1+r_1^{\prime }-1,x_1+y_1)}{\hat{q}}_0(y;x_1),$$

$$\pi ({\lambda }_1,{\lambda }_2\mid x_1,x_2)=\frac{\pi ({\lambda }_1)\pi ({\lambda }_2)p_{r_1,{\lambda }_1}(x_1)p_{r_2,{\lambda }_2}(x_2)}{{\int }_0^{\infty }{\int }_0^{{\lambda }_1}\pi ({\lambda }_1)\pi ({\lambda }_2)p_{r_1,{\lambda }_1}(x_1)p_{r_2,{\lambda }_2}(x_2)d{\lambda }_{2}d{\lambda }_1},$$

$${\int }_0^{\infty }\pi ({\lambda }_1)p_{r_1,{\lambda }_1}(x_1){\int }_0^{{\lambda }_1}\pi ({\lambda }_2)p_{r_2,{\lambda }_2}(x_2)d{\lambda }_{2}d{\lambda }_1.$$

$$\frac{1}{r_2-1}{\int }_0^{{\lambda }_1}\pi ({\lambda }_2){\mathrm{IG}}_{r_2-1,x_2}({\lambda }_2)d{\lambda }_2=\frac{m(r_2-1,x_2;{\lambda }_1)}{r_2-1},$$

$$\frac{1}{(r_1-1)(r_2-1)}$$

$$=\frac{C(r_1-1,x_1,r_2-1,x_2)}{(r_1-1)(r_2-1)}.$$

$$\frac{\pi ({\lambda }_1)\pi ({\lambda }_2)\frac{{\lambda }_1^{-r_1}}{\Gamma (r_1-1)}\frac{{\lambda }_2^{-r_2}}{\Gamma (r_2-1)}{x}_1^{r_1-1}{x}_2^{r_2-1}{e}^{-x_1/{\lambda }_1}{e}^{-x_2/{\lambda }_2}}{C(r_1-1,x_1,r_2-1,x_2)},$$

$$\pi ({\lambda }_1\mid x_1,x_2)$$

$$=\frac{\pi ({\lambda }_1)\frac{{\lambda }_1^{-r_1}}{\Gamma (r_1-1)}{x}_1^{r_1-1}{e}^{-x_1/{\lambda }_1}}{C(r_1+-1,x_1,r_2-1,x_2)}{\int }_0^{{\lambda }_1}\pi ({\lambda }_2)\frac{{\lambda }_2^{-r_2}}{\Gamma (r_2-1)}{x}_2^{r_2-1}{e}^{-x_2/{\lambda }_2}$$

$$=\pi ({\lambda }_1)\frac{m(r_2-1,x_2;{\lambda }_1)}{C(r_1-1,x_1,r_2-1,x_2)}\frac{{\lambda }_1^{-r_1}}{\Gamma (r_1-1)}{x}_1^{r_1-1}{e}^{-x_1/{\lambda }_1},{\lambda }_1>0.$$

$${\hat{q}}_1(y_1;x_1,x_2)$$

$$=\frac{1}{\Gamma (r_1^{\prime })\Gamma (r_1-1)}\frac{x_1^{r_1-1}{y}_1^{r_1^{\prime }}}{C(r_1-1,x_1,r_2-1,x_2)}{\int }_0^{\infty }\pi ({\lambda }_1)m(r_1-1,x_1)e^{-(x_1+y_1)/{\lambda }_1}{\lambda }_1^{-(r_1+r_1^{\prime }-1)-1}d{\lambda }_1$$

$$=\frac{1}{B(r_1-1,r_1^{\prime })}\frac{C(r_1+r_1^{\prime }-1,x_1+y_1,r_2-1,x_2)}{C(r_1-1,x_1,r_2-1,x_2)}\frac{x_1^{r_1-1}{y}_1^{r_1^{\prime }}}{(x_1+y_1{)}^{r_1+r_1^{\prime }-1}},$$

$${\hat{q}}_1(y_1;x_1,x_2)=\frac{r_1{x}_2^{-r}{y}_1^r,\Gamma (r^{\prime }+r_1+r_2){}_2{F}_1\left(r^{\prime }+r_1,r^{\prime }+r_1+r_2;r^{\prime }+r_1+1;-\frac{x_1+y_1}{x_2}\right)}{(r^{\prime }+r_1)\Gamma (r^{\prime })\Gamma (r_1+r_2){}_2{F}_1\left(r_1,r_1+r_2;r_1+1;-\frac{x_1}{x_2}\right)}.$$

$$N_1\mid {\theta }_1\sim \mathrm{Po}\left({\theta }_1\frac{R_1}{R_2}\right),N_2\mid {\theta }_2\sim \mathrm{Po}\left({\theta }_2\frac{R_2}{R_1}\right),$$

Full text

Keywords and phrases: Ancillary information, Hockey, Predictive density estimation, Weighted beta prime distribution

1 Introduction

Predicting an unobserved random variable by finding the future density is often more meaningful for applications than conventional point and interval estimation of parameters, because it gives a richer prediction and all point or interval estimation, percentile, and even hypothesis testing can be obtained from the estimated density. Although in practice there is a plug–in approach in density estimation problems, Bayesian predictive distributions are a better alternative in the literature in many cases (see e.g. Marchand and Sadeghkhani 2018) and most of the time are easily calculated by using modern Monte Carlo techniques.

Often, there exists some ancillary information at our disposal which has been unduly ignored. Perhaps the most straightforward way to visualize this additional information is to place restrictions on the parameters of our model and translate them into parametric restrictions. In this paper, we focus on the gamma distribution as a waiting scoring model. Previous studies indicate estimating of future density based on current observed gamma population. For example, L’moudden et. al. (2017) studied Bayesian estimation of a future random variable of gamma distribution when the scale parameter is bound within an interval. This paper, uses ancillary information obtained from two independent gamma densities, to make a more accurate density estimation of one of them.

The remainder of the paper is organized as follows: In Section 2, we provide definitions and preliminary remarks which will be needed throughout this article. Section 3 discusses how to find a Bayesian predictive density estimator in gamma distribution with and without contemplating the ancillary information. In Section 4 we study an interesting application of proposed methods in estimating the future density of waiting time until scoring $r$-th goal in a hockey game. Finally, we make some concluding remarks in Section 5.

2 Main results

Let $X\mid\lambda\sim\operatorname{Gam}(r,\lambda)$, be a random variables (rv) from a gamma distribution with probability distribution function (pdf)

[TABLE]

where $r>0$ is known and $\lambda>0$ unknown. Suppose that we are interested in estimating the future density of unobserved rv $Y$ with pdf $q_{\lambda}(y)=\frac{y^{r^{\prime}-1}e^{-y/\lambda}}{\Gamma(r^{\prime})\,\lambda^{r^{\prime}}}\,,\,\,\,y>0$. Imagine that such an estimator can be obtained by $\hat{q}_{1}(\cdot\,;x)$. We use the Kullback–Leibler (KL) loss (equation 2.1) in order to compare the efficiency of the proposed estimator with the actual $q(\cdot)$.

[TABLE]

and the corresponding frequentist risk function is given as

[TABLE]

Previous studies (see, Corcuera and Giummolè 1999), indicates that under KL, the Bayes predictive density estimator for $Y$ based on observed $x$, prior and posterior $\pi(\cdot)$ and $\pi(\cdot\mid\mathbf{x})$ respectively, is given as

[TABLE]

The following contains some definitions and remarks that will be needed later on.

Definition 2.1**.**

The pdf and the cumulative density function (cdf) of a rv $T$ follows the inverse–gamma distribution, if we have

[TABLE]

where $\Gamma(m,n)$ is known as an upper incomplete gamma function and is defined as $\int_{n}^{\infty}t^{m_{1}}e^{-t}\,dt$.

Remark 2.1**.**

One can verify the following statements:

The marginal distribution of a random variable $\operatorname{IG}(s_{1},s_{2})$ associated with prior $\pi(\xi)$ is

[TABLE]

Note that if $\pi(\xi)=1/\xi$ for $\xi>0$, then

[TABLE] 2. 2.

[TABLE]

when $\pi(\xi)=1/\xi$ for $0\leq\xi\leq\bar{\xi}$, we have

[TABLE]

Noting $\bar{\xi}\to\infty$, $\Gamma(1+s_{1},0)=\Gamma(1+s_{1})$ and hence (2.8) and (2.9) are equal. 3. 3.

[TABLE]

In the case of $\pi(\xi)=\frac{1}{\xi_{1}\xi_{2}}$, for $\xi_{1}\geq\xi_{2}>0$, after some calculation we get

[TABLE]

*where ${}_{2}\tilde{F}_{1}$ is known as a regularized hypergeometric function. *

Definition 2.2**.**

A random variable $T$ is said to have a generalized beta prime density, whenever $T$ has the following distribution:

[TABLE]

where $\operatorname{B}(\cdot,\cdot)$ is the beta function, $\sigma>0$ is the scale parameter, and $a,b,\gamma>0$ are shape parameters. We denote it by $\operatorname{GB}^{\prime}(a,b,\gamma,\sigma)$. $\operatorname{GB}^{\prime}(a,b,\sigma,\gamma=1)$ is known as three parameter beta prime $\operatorname{B}^{\prime}(a,b,\sigma)$. Beta prime (also known as a beta distribution of the second kind) is obtained by $\operatorname{B}^{\prime}(a,b,\sigma=1)$ and denoted by $\operatorname{B}^{\prime}(a,b)$

3 Bayes predictive density estimation

The following provides the Bayes predictive density estimator under KL and prior density $\pi(\lambda)$.

Lemma 3.1**.**

The posterior predictive density for $Y_{1}\mid\lambda_{1}\sim\operatorname{Gam}(r^{\prime}_{1},\lambda_{1})$, based on $X_{1}\mid\lambda_{1}\sim\operatorname{Gam}(r_{1},\lambda_{1})$ and for for a prior $\pi(\lambda_{1})$, is given by

[TABLE]

Proof.

Using the notation of $\operatorname{Gam}_{r,\lambda}(t)$ for the pdf of rv $X$, which follows the gamma distribution with parameters of $r$ and $\lambda$, and the fact that $\operatorname{Gam}_{r,\lambda}(x)=\frac{\lambda}{x}\operatorname{IG}_{r,x}(\lambda)=\frac{1}{r-1}\operatorname{IG}_{r-1,x}(\lambda)$, for $x>0,r>0,\lambda>0$, along with equation (2.6), one may verify that $p(x_{1})=\int_{0}^{\infty}\frac{\pi(\lambda_{1})}{\Gamma(r_{1})\,\lambda_{1}^{r_{1}}}x_{1}^{r_{1}-1}e^{-\frac{x_{1}}{\lambda_{1}}}\,d\lambda_{1}=\frac{m(r_{1}-1,x_{1})}{r_{1}-1}$. The posterior density is given by

[TABLE]

Thus the posterior predictive density is given by

[TABLE]

and hence the proof. ∎

Remark 3.2**.**

Assuming the non–informative prior $\pi(\lambda)=\frac{1}{\lambda_{1}}$ for $\lambda_{1}>0$ (and therefore $m(s_{1},s_{2})=s_{1}/s_{2}$) to Lemma 3.1, the posterior predictive distribution can be simplified as

[TABLE]

which is known as a three parameter beta prime distribution, namely, $Beta^{\prime}(r^{\prime}_{1},r_{1},x_{1})$, as first observed by Aitchison (1975). L’moudden and Marchand (2017) also use Bayes predictive density estimator for a gamma model under Kullback–Leibler loss, when the scale parameter is bound in an interval $(a,b)$ for $b>a>0$, and show that restricted Bayes predictive density estimator dominates the unrestricted density estimator in (3.14).

3.1 Improving the Bayes predictive density estimation upon ancillary information

Often, we have at our disposal additional information from which we may gain a greater understanding, if correctly employed. In order to incorporate the ancillary information, let us assume the following scenario:

Let $X_{i}$ for $i=1,2$ be independent $\operatorname{Gam}(r_{i},\lambda_{i})$, $y_{1}$ has $\operatorname{Gam}(r^{\prime}_{1},\lambda_{1})$, and $y_{1}$ and $X_{i}$, are independent where $r_{i}>0$, $r^{\prime}_{1}$ are known and $\lambda_{i}>0$ is unknown but we have a prior information that the ratio of scale parameters satisfy $\frac{\lambda_{1}}{\lambda_{2}}\geq 1$ (i.e. $\lambda_{1}\geq\lambda_{2}>0$). The following theorem provides a posterior predictive density estimator subject to constraints placed on the ratio of scale parameters by assuming that $\lambda_{1}$ and $\lambda_{2}$ are independent; that is, $\pi(\lambda_{1},\lambda_{2})=\pi(\lambda_{1})\pi(\lambda_{2})$.

Theorem 3.1**.**

The posterior predictive density of $Y_{1}\sim\operatorname{Gam}(r^{\prime}_{1},\lambda_{1})$, based on $X_{1}\sim\operatorname{Gam}(r_{1},\lambda_{1})$ and $X_{2}\sim\operatorname{Gam}(r_{2},\lambda_{2})$, where $r_{i}>0$, $r^{\prime}_{1}>0$ and $\lambda_{i}>0$ are unknown associated with the prior information $\pi(\lambda_{1},\lambda_{2})=\pi(\lambda_{1})\pi(\lambda_{2})I_{(0,\infty)}(\lambda_{1})\,I_{(0,\lambda_{1})}(\lambda_{2})$ satisfy $\frac{\lambda_{1}}{\lambda_{2}}\geq 1$, is given by

[TABLE]

which is a weighted function of $\hat{q}_{0}(y;x_{1})$, the unrestricted posterior predictive density of $y_{1}$ based on $x_{1}$ obtained in equation (3.13), and $m(\cdot,\cdot)$ and $C(\cdot,\cdot)$ are given in equations (2.8) and (2.10), respectively.

Proof.

The joint posterior density is given by

[TABLE]

The denominator in equation (3.16) can be written as

[TABLE]

Applying notations in (2.8) and (2.10), the inner integral in the above quotation is equal to

[TABLE]

therefore the denominator is obtained by

[TABLE]

Therefore the posterior density (3.16) has the following form:

[TABLE]

so the marginal posterior is given by

[TABLE]

Substituting the above equation in the formula of posterior predictive density estimator, gives

[TABLE]

and replacing $\hat{q}_{1}(y_{1};x_{1})$ from (3.13) in the above equation thus completes the proof. ∎

Corollary 3.1**.**

In the case of non informative prior $\pi(\lambda_{1},\lambda_{2})=\frac{1}{\lambda_{1}\lambda_{2}}$, for $\lambda_{1}\geq\lambda_{2}\geq 0$, the Bayesian predictive density estimator given in equation (3.15) is in fact the weighted beta prime distribution and is given by

[TABLE]

4 Real examples and illustrations

In this section, we apply the proposed methods in Section 3 in order to find the Bayes predictive density estimators based on the National Hockey League’s dataset for the $2017-2018$ season to predict the upcoming $2018-2019$ season. We are concerned with the waiting time prediction until scoring $r$-th goal (in a real scenario $r$ rarely exceeds $6$) by team A as compared to team B. Similar explanations have been done by Suzuki et al. (2010) in predicting match outcomes football World Cup 2006. We shall assume that $N_{1}$ and $N_{2}$, the number of goals scored by team A and B, are two independent, random Poisson variables such that

[TABLE]

where $\theta_{1}$ and $\theta_{2}$ can be interpreted as the mean number of goals teams A and B score against their opponents, and $R_{1}$ and $R_{2}$ are the teams A and B points in the previous $2017-2018$ season. This assumption seems quite reasonable, since it is expected a that team with a higher point scores more goals rather than lower point team.

Equivalently, we can assume the waiting time to see team A scores $r_{1}$ goals to team B and team B scores $r_{2}$ goals to team A independently follow

[TABLE]

The goal is to estimate the future density for the expected time to seeing $r_{1}$ goal by team A to team B, $Y_{1}\mid\theta_{1},r^{\prime}_{1}\sim\operatorname{Gam}\left(r^{\prime}_{1},\theta_{1}\frac{R^{\prime}_{1}}{R^{\prime}_{2}}\right)$, using ancillary information $\theta_{1}\geq\theta_{2}$ (i.e., by contemplating that team A is more capable of scoring against team B based on previous records and$/$or the predictions of sports analysts.

As an example, we consider two Canadian hockey team, the Toronto Maple Leafs (A) and Montreal Canadiens (B). Based on the 2017–18 season teams’ points, available at www.nhl.com, $R_{1}=105$ and $R_{2}=71$, respectively. Table 1 (on page 10) shows the time elapsed (in minutes) until scoring the third goal in games ($r_{1}=r_{2}=3$) played by the Toronto Maple Leafs versus the Montreal Canadiens in the 2017--2018 National Hockey League season.

In Figure 1, the Bayesian predictive density estimators $\hat{q}_{0}$ and $\hat{q}_{1}$ (based on equations 3.14 and 3.18 with the exact densities given in Table 2) for the elapsed time until the Toronto Maple Leafs scores the $3^{rd}$ goal to the Montreal Canadiens in an upcoming game are depicted. Note that $\hat{q}_{1}$ is obtained by employing ancillary information that the performance of the Toronto Maple was better the Montreal Canadiens (specialists’ opinions or their points) i.e. $\theta_{1}/\theta_{2}\geq 1$ as well as contemplating their points. It is important to mention here that the densities in here have been truncated to $(0,60)$ minutes since NHL hockey games consist of three periods, each lasting $20$ minutes.

Table 2 contains the density of obtained estimator $\hat{q}_{0}$ and the Bayes estimators $q_{1}$, along with their modes, means, and $10^{th}$, $50^{th}$ and $90^{th}$ percentiles for the future density $Y_{1}$ of waiting time until scoring $3^{rd}$ score by the Toronto Maple Leafs based on the data from Table 1.

For instance, if we contemplate the ancillary information as we described above, we are anticipating to see the third goal of the Toronto Maple Leafs around minute $28$, rather than minute $18$ of the game, or we are expecting on average the Toronto Maple Leafs score the third goal versus their opponents around minute $28$ of each game. However, we apply the information at our disposal to improve our estimation, we hope to see the third goal of the Toronto Maple Leafs versus the Caniadiens around minute $33$ of the game, or in less than $20$ percent of matches Toronto scores the third goal at $14.38$. By applying the available information, this time moves to end of first period of the game.

4.1 Prediction errors

We calculate the the prediction error (pe) of $\hat{q}_{0}$ and $\hat{q}_{1}$ in order to investigate their performance under KL loss given in equation (2.1). We verify their pe’s based on the $2016-17$ season and see their performance by comparing the KL loss (distance) Bayesian density estimator and exact density for waiting time until the Toronto Maple Leafs scores the $3^{rd}$ goal. This can be calculated by

[TABLE]

where $Y_{1}$, waiting time until the Toronto Maple Leafs scores the $3^{rd}$ goals based on the $2016-17$ dataset, follows a truncated gamma $\operatorname{Gam}(3,18.3)$ on $(0,60)$ (with an expected value of $35.8$). It can be verified that the prediction error of our estimator when the ancillary information is taken into account $pe(\hat{q}_{1})=0.04$ while $pe(\hat{q}_{0})=0.45$ and $\hat{q}_{1}(y_{1})$ outperforms other estimator, in estimating the density of waiting time until the $3^{rd}$ goal occurs in a match Toronto Maple Leafs versus the Montreal Canadiens in a National Hockey League’s game.

However, the frequentist risks of $\hat{q}_{1}(y_{1})$ and $\hat{q}_{0}(y_{1})$, from equation (2.2), can be used to verify the dominance result of our proposed estimator $\hat{q}_{1}(y_{1})$ in order to make a prediction the result of upcoming match the Toronto Maple Leafs versus the Montreal Canadiens in the $2018-19$ season based on data from the $2017-18$ season. Figure 2 depicts the performance of $\hat{q}_{1}(y_{1})$ and $\hat{q}_{0}(y_{1})$ via the frequentist risk functions. Since there is no ancillary information (i.e. $\theta_{1}=\theta_{2}$ and $R_{1}/R_{2}=1$), both risks are equal and as long as the ratio of $\theta_{1}/\theta_{2}$ increases (gain from ancillary information), the risk of $\hat{q}_{1}(y_{1})$ decreases (the risk of $\hat{q}_{0}[y_{1}]$ is constant as we anticipated, since it is an MRE estimator) and both risk functions converge to each other due to lack of information.

5 Concluding remarks

In summation, our results demonstrate that applying additional ancillary information can be used to come up with a more accurate prediction via a better Bayesian predictive density estimator for the future density of a gamma. Our proposed density estimator for the waiting time until scoring the $r$-th goal outperforms the usual MRE estimator. This was verified by comparing the prediction errors of two density estimators as well as their frequentist risk functions.

Acknowledgement

The authors thank the Stathletes company specially Meghan Chayka and Jeff Goeree for providing the data and helpful comments to this manuscript. The Natural Sciences and the Engineering Research Council of Canada, and Ontario Centre of Excellence supported the research of Professor S. Ejaz Ahmed.