\section{Introduction}
Quantitative analysis in sports analytics has received a great deal of attention in recent years. Traditional statistical research for sports analytics mainly focused on game result prediction, such as predicting the number of goals scored in soccer matches \citep{dixon1997modelling, karlis2003analysis,baio2010bayesian}, and the basketball game outcomes \citep{carlin1996improved,caudill2003predicting,cattelan2013dynamic}. More recently, fast development in game tracking technologies has greatly improved the quality and variety of collected data sources \citep{albert2017handbook}, and in turn substantially expanded the role of statistics in sports analytics, including performance evaluation of players and teams  \citep{cervone2014pointwise,franks2015characterizing,wu2018modeling,hu2020bayesian,hu2021cjs}, commentator's in-game analysis and coach's decision making \citep{fernandez2018wide,sandholtz2019measuring}.

The home field advantage is a key concept that has been studied across different sports for over a century. The first public attempt to statistically quantify the home field advantage was about 40 years ago \citep{physicaleducation}, where \citet{schwartz1977home} studied the existence of home advantage in several sports and found that its effect was most pronounced in the indoor sports such as ice hockey and basketball.
Since then, statisticians and psychologists have dug their heels into trying to map out this anomaly, and there is some consensus on it. The percentages found through different studies in general float around 55-60\% of an advantage to the home team (wherein 50\% would mean there is no advantage either way). \citet{courneya1992home} found that different sports in general held different advantage percentages: 57.3\% for American football versus 69\% for soccer, for example. In a study of basketball home court advantage, \citet{calleja2018brief} found that the visiting team scored 2.8 points less, attempted 2 fewer free throws, and made 1 less free throw than their season averages. Season length is also shown to have a negative correlation with the strength of home field advantage, where sports with 100 games or more per season were shown to have significantly less advantage than a sport with fewer than 50 games per season \citep{Metaanalysis}. This finding can be explained by the fact that each game becomes less important on average as the number of games increases. Some studies have shown that the home field advantage becomes less significant in the modern era. For example, \citet{Metaanalysis} found that the home field advantage was the strongest in the pre-1950s, as opposed to any other twenty-year block leading up to 2007, the year he conducted the study. In another work, \citet{sanchez2009analysis} showed that between 2003 and 2013, the home field advantage for the UEFA Champions League decreased by 1.8\%, supporting the claim that as time goes on, the home field advantage decreases to a certain level and remains still. In general, the home field advantage tends to deteriorate as the athletes play longer on an away field.

In principle, there are four major factors associated with the home field advantage: crowd involvement, travel fatigue, familiarity of facilities, and referee bias that benefits the home teams. Several studies \citep{schwartz1977home,greer1983spectator} have shown a positive correlation between the size of the crowd and the effect of the home field advantage, an advantage that can get as high as 12\% over a team's opponents. Travel fatigue is another important factor for home field advantage \citep{Metaanalysis,calleja2018brief,NFLtimezone}, since the athletes have an overall reduction in mean wellness, including a reduction of sleep, self-reported feelings of jet lag and energy reduction; the importance of sleep and proper recovery cannot be over-exaggerated, and it can be difficult to do either of those properly in unfamiliar settings (e.g., a hotel room or a bus). Furthermore, familiarity of facilities gives the home team a strong advantage over a visiting team \citep{MovingNewStadium}. A home team not only does not have to travel to an unfamiliar facility, be away from home, and spend hours transporting themselves; they also get to use their own facility, locker room, and field. It cannot be underestimated the value of knowing and knowing well the nooks and crannies of a facility. Ultimately, giving rule advantages to the home team is very common in certain sports \citep{courneya1992home}. For example, in baseball the home team gets the advantage of batting last, giving them the final opportunity to score a run in the game. In hockey, the last line change goes to the home team.

Soccer has proven to have one of the highest home field advantages among major sports in a majority of current research. \citet{MovingNewStadium} states that the home field advantage in soccer is equivalent to 0.6 goals per game and the visual cues that come from being intimately familiar with a facility can be exponentially helpful in fast paced sports such as soccer, which could perhaps explain why the home field advantage is less pronounced in slower, stopping sports such as baseball. \citet{sanchez2009analysis} observes that the better a soccer team is, the more often the home field advantage appears, after studying home field advantage across groups of variably ranked soccer teams.

Despite the vast amount of existing work on home field advantage quantification in soccer and other sports \citep{leitner2020no,benz2021estimating,fischer2021does}, very little is known about its {\it causal effect} on team performance; and it is our goal in this paper to fill this gap. In particular, we study soccer games by analyzing a data set collected from the English Premier League (EPL) 2020-2021 season, where 380 games were played, that is, two games between each pair of 20 teams. We choose eleven team-level summary statistics as the main outcomes that represent team's performance on defensive and offensive sides, and the referee bias. We then develop a new causal inference approach for assessing the causal effect of home field advantage on these outcomes. More details about our data application is provided in Section \ref{sec:data}.

In causal inference literature (for observational studies), causal effect estimation of a binary treatment is a classic topic that has been intensively studied \citep{licausal,imbens2004nonparametric}.
Over the past few decades, a number of methods in both statistics and econometrics have been proposed to identify and estimate the average treatment effects of a particular event/treatment \citep[see a recent overview in][]{athey2017estimating} under the assumption of ignorability (also unknown as unconfounded treatment assignment) \citep{rosenbaum1983central}, including regression imputation, (augmented) inverse probability weighting \citep[see e.g., ][]{horvitz1952generalization, rosenbaum1983central, robins1994estimation, bang2005doubly, cao2009improving}, and matching \citep[see e.g., ][]{rubin1973matching, rosenbaum1989optimal, heckman1997matching, hansen2004full, rubin2006matched, abadie2006large, abadie2016matching}. However, the data structure in soccer games is distinguished from those in the existing literature, leading to {\it two  unique challenges} in identifying the causal effects of home field advantage. First, the EPL data set is obtained neither from a randomized trial or observational studies as standard in casual inference literature, but instead is based on a collection of pair-wise matches between every two teams. For one season, each team will have one home game and one away game with each of the other 19 opponents. All matches are pre-scheduled. Thus, propensity-score-based methods can hardly gain efficiency with such design. Secondly, for each match, there is one team with home field advantage and the other without, i.e., there are no matches in neutral field. In fact, this is the common practice for other major professional sports leagues such as NBA, NFL, and MLB. In other words, there is technically no control group where both teams have no home field advantage. Hence, we cannot rely on matching-based method to estimate the casual effects. To overcome these difficulties, in this paper, we establish a  hierarchical causal model to characterize the underlying true causal effects at the league and team levels, and propose a novel causal estimation approach for home field advantages with inference procedures.

Our proposed method is unique in the following aspects. First, the idea of {\it pairing home and away games} for solving causal inference problem is novel. In fact, this idea and our proposed approach are widely applicable to general sports applications such as football, baseball and basketball studies, and provide a valuable alternative to the existing literature that mainly relies on propensity score. Secondly, under the proposed  hierarchical  model framework, both league level and team level causal effects are identifiable and can be conveniently estimated. Moreover, our inference procedure is developed based on linear model theory that is accessible to a wide audience including first-year graduate students in statistics. Thirdly, our real data analysis results reveal several interesting findings from England Primer League, which may provide new insights to practitioners in sports industry. It is our hope that the data application presented in this paper as well as the developed statistical methodology can reach to a wide range of audience, be useful for educational purpose in statistics and data science classes, and stimulate new ideas for more causal analysis in sports analytics.

The rest of the paper is organized as follows. In Section \ref{sec:data}, we
give an overview of our motivating data application and introduce several representative summary outcomes. In Section~\ref{sec:model}, we introduce the causal inference framework and propose a  hierarchical causal model for both team-level and league-level home field advantage effect estimation. Extensive simulation studies are presented in Section~\ref{sec:simu} to investigate the empirical performance of our approach.
We apply our method to analyze the 2020-2021 England Primer League in Section~\ref{sec:app} and conclude with a discussion of future directions in Section~\ref{sec:disc}.

\section{Motivating Data}\label{sec:data}
We are interested in studying the causal effect of home field advantage for in-game performance during English Premier League 2020-2021 season. Despite as obvious as it sounds, the home field advantage is in fact not evident at first glance through the usual descriptive statistics such as the game outcomes and number of scored goals. For example, among the 380 games in that season, the number of home wins is 144 (37.89\%), which is even less than the number of away wins, 153 (40.26\%). Among the total of 1024 goals scored in that season, 514 (50.19\%) were scored by the home team, and 510 (49.81\%) by the away team.

To better understand this phenomenon, we study a data set provided by Hudl \& Wyscout, a company that excels at soccer game scouting and match analysis. The data is collected from 380 games played by 20 teams in EPL, and includes
a number of statistics collected from each game that range across the defensive and offensive capabilities of the home and away teams and players.
We choose to focus on eleven in-game statistics as follows:
\begin{itemize}
    \item \textbf{Attacks w/ Shot} - The number of times that an offensive team makes a forward move towards their goal (a dribble, a pass, etc) followed by a shot;
    \item \textbf{Defence Interceptions} - The number of times that a defending team intercepts a pass;
    \item \textbf{Reaching Opponent Box} - The number of times that a team moves the ball into their opponent's goal box;
    \item \textbf{Reaching Opponent Half} - The number of times that a team moves the ball into their opponent's half;
    \item \textbf{Shots Blocked} -  The number of times that a defending team deflects a shot on goal to prevent scoring;
    \item \textbf{Shots from Box} - The number of shots taken from the goal box;
    \item \textbf{Shots from Danger Zone} -  The number of shots taken from the "Danger Zone", which is a relative area in the center of the field approximately 18 yards or less from the goal;
    \item \textbf{Successful Key Passes} - The number of passes that would have resulted in assists if the resulting shot had been made;
    \item \textbf{Touches in Box} - Number of  passes or touches that occur within the penalty area;
    \item \textbf{Successful Key Passes} - The number of passes that would have resulted in assists if the resulting shot had been made;
    \item \textbf{Expected Goals (XG)} - The average likelihood a goal will be scored given the position of the player over the course of a game;
    \item \textbf{Yellow Cards} - The number of yellow cards given to a team in a game.
\end{itemize}

Table \ref{tab:summary_statistics} contains the eleven summary statistics from the raw data, where each row shows the statistic, its primary role in the game (defense, offense, or referee), the means for the home team and for the away team, respectively, and the overall standard deviation. These in-game statistics are chosen because they are most relevant to studying the home field advantage and also sufficient to cover different aspects of the soccer games (e.g., team offensive and defensive performance).
\begin{table}[th]
    \centering
        \caption{Selective summary statistics for 2020-2021 England Primer League.}
    \label{tab:summary_statistics}
    \begin{tabular}{  c  c c  c  c  }
  \toprule
  \textbf{Statistic} & \textbf{Role} & \textbf{Mean Home} & \textbf{Mean Away} & \textbf{Overall SD}\\

  \midrule
   Attacks w/ Shot & Offense & 11.547 & 9.979 & 4.910 \\

  Defence Interceptions & Defense & 42.639 & 44.637 & 11.287  \\

  Reaching Opponent Box & Offense &  14.779 & 13.097 & 6.135\\

  Reaching Opponent Half & Offense & 58.561 & 54.968 &  12.566\\

  Shots Blocked & Defense & 3.382 & 2.95 & 2.273 \\

  Shots from Box & Offense & 7.5 & 6.434 & 3.594 \\

  Shots from Danger Zone & Offense & 5.058 & 4.271 &  2.622\\

  Successful Key Passes & Offense & 3.584 & 3.037 & 2.328\\

  Touches in Box & Offense & 19.468 & 17.216 & 8.883 \\

  Expected Goals (XG) & Offense & 1.577 & 1.395 & 0.867\\

  Yellow Cards & Referee & 1.447 & 1.474 & 1.151\\
  \bottomrule
\end{tabular}

\end{table}

To better understand the role of these selected statistics. We choose three of them (representing offense, defense, and referee), calculate the difference in these statistics between home and away games for each team and present their distributions in Figure~\ref{fig:reaching} to Figure~\ref{fig:yellow}.
The offensive statistic is Reaching Opponent Half, and its distribution for different teams is illustrated in Figure \ref{fig:reaching}. The difference between home and away teams, in theory, should be positive overall if the home field advantage is present, negative if there is actually an advantage towards the away team, and 0 if there is no advantage either way. From the picture, we can see a clear trend for home advantage, e.g., 17 of the 20 teams have a positive mean, which implies that they reach their opponents' half more often on average when they play at home as opposed to away. This finding is in fact fairly consistent across all offensive variables. The three teams that do not show evidence of home field advantage for Reaching the Opponent's Half are Aston Villa, Chelsea, and Southampton.  The top three teams from that season (Manchester City, Manchester United, and Liverpool) exhibit similar patterns in the picture.

\begin{figure}
    \centering
   \includegraphics[width=0.95\textwidth]{fig/Reaching_opp_half.png}
    \caption{Net difference of Reaching Opponent Half for 20 teams.}
    \label{fig:reaching}
\end{figure}

 The Defense Interceptions is chosen as a representation for the defensive statistic; and we present its distribution of the difference between home and away teams in Figure \ref{fig:Interception}. Since defence interceptions negatively impact a team, the home field advantage will be seen here if, oppositely to the offensive statistic, the distribution of the difference is skewed negatively. Five teams do not show evidence of home field advantage: Chelsea, Crystal Palace, Everton, Fulham, and Southampton. Interestingly, the top three teams from that season (Manchester City, Manchester United, and Liverpool) all have negative means but show a different level of variability in this statistic.
\begin{figure}
    \centering
   \includegraphics[width=0.95\textwidth]{fig/Defence_Interceptions_1.png}
    \caption{Net difference of Defence Interceptions for 20 teams.}
    \label{fig:Interception}
\end{figure}

Finally, the referee bias is represented by the team-level number of yellow cards in each match. The difference between a team's yellow card calls at home versus away is shown in Figure~\ref{fig:yellow}. If there was a home field advantage effect in yellow cards given, the boxplot would trend negatively. Different from the offensive and defensive statistics, we see a majority of teams centralizing their distributions around the zero line, which is an indication that the referee calls are similar whether the team is home or away. Out of 20 teams, 14 have means about the zero, and the other 6 have a mix of positive and negative means.

\begin{figure}
    \centering
   \includegraphics[width=0.95\textwidth]{fig/Yellow_Cards_1.png}
    \caption{Net difference of Yellow Cards for 20 teams.}
    \label{fig:yellow}
\end{figure}

\section{Method}\label{sec:model}
\subsection{Framework and Assumptions}
Suppose there are $n$ different teams in the league, denoted by $\{T_1,\cdots, T_n\}$, with a total number of $N\equiv n(n-1)$ matches between any pair of two teams, $T_i$ and $T_j$, where $i\not = j$ and $i,j \in \{1,\cdots,n\}$. Without loss of generality, we assume $n\geq 3$.
Therefore, each team $T_i$ has $(n-1)$ matches with home field advantage and another $(n-1)$ without.
Define a treatment indicator $\delta_{i}=1$ if team $T_i$ has the home field advantage in a match against team $T_j$ for $j\not = i$, and $\delta_{i}=0$ otherwise.
By definition, we have $\delta_{i} + \delta_{j} =1$ for any match between $T_i$ and $T_j$.
The main outcome is defined as the net difference in the outcome of interest, i.e., $Y_{i,j} =(c_i - c_j)\in\mathbb{R}$, where $c_i$ and $c_j$ are one of the eleven statistics that we described in Section \ref{sec:data} for team $T_i$ (with home field advantage) and team $T_j$ (without home field advantage), respectively. Following the potential outcome framework \citep[see e.g., ][]{rubin1974estimating}, we define the potential outcome $Y^*(\delta_{i}=a,\delta_{j}=1-a)$ as the outcome of interest that would be observed after the match between team $T_i$ and team $T_j$, where $a=1$ or $a=0$ corresponds to that team $T_i$ or team $T_j$ has the home field advantage, respectively. As standard in the causal inference literature \citep[see e.g., ][]{rosenbaum1983central}, we make the following assumptions for any pairs of $i \neq j$.

\noindent \textbf{(A1)}. Stable Unit Treatment Value Assumption:
\begin{equation*}
\begin{split}
Y_{i,j}= \delta_{i}(1- \delta_{j})Y^*(\delta_{i}=1,\delta_{j}=0)  +(1- \delta_{i}) \delta_{j} Y^*(\delta_{i}=0,\delta_{j}=1).
\end{split}
\end{equation*}

\noindent \textbf{(A2)}. Ignorability:
\begin{equation*}
\{Y^*(\delta_{i}=1,\delta_{j}=0),Y^*(\delta_{i}=0,\delta_{j}=1)\}\independent \{\delta_{i},\delta_{j}\}.
\end{equation*}

\noindent \textbf{(A3)}. Non-monotonicity: If team A dominates team B and team B dominates team C, team C can still dominate team A with a certain probability.

Assumptions (A1) and (A2) are standard in causal inference literature \citep[see e.g.,][]{athey2017estimating}, to ensure that the causal effects are estimable from observed data. By game design, each team will compete the rest teams with and without home field advantage once respectively, hence the ignorability assumption holds automatically in our study. Assumption (A3) is to rule out the correlation between different matches that involves the same team. In reality, there are style rivalries between teams in soccer games. Hence the transitive relation does not always hold.

\subsection{Hierarchical Causal Modeling and Estimation}
 In this section, we detail the proposed hierarchical causal model and its estimation and inference procedures. Specifically, we are interested in team level and league level estimators of
causal effects. To this end, we define the causal effect of home field advantage associated with team $T_i$ as
\begin{eqnarray}\label{def_betai}
\beta_i = {\mbox{E}}_j \{Y^*(\delta_{i}=1,\delta_{j}=0)-Y^*(\delta_{i}=0,\delta_{j}=0)\},~~\text{for}~i=1,\ldots,n,
\end{eqnarray}
and the causal effect of home field advantage for the entire league as
\begin{eqnarray}\label{def_delta}
\Delta = {\mbox{E}}_i \{\beta_i\}.
\end{eqnarray}
Given finite number of teams in the league, we are interested in two estimands, the average home field advantage of the league $\Delta\equiv \sum_{i=1}^n \beta_i /n $ and the team-specific home field advantage $\beta_i$.
Yet, since there is no match in neutral field in major professional sports, we always have $\delta_{i} + \delta_{j} =1$, i.e.,  $Y^*(\delta_{i}=0,\delta_{j}=0)$ can never be observed. To address this difficulty and estimate $\beta_i$ from the observational studies, we propose to decompose the outcome function into two parts, one corresponding to the home field advantage and the other representing the potential outcome in a hypothetical neutral field, via a mixed two-way ANOVA design. To be specific, the outcome of a match between team $T_i$ (with home field advantage) and team $T_j$ (without home field advantage) is modeled by
\begin{eqnarray}\label{main_model_beta}
 Y_{i,j}=\alpha_{i,j}+ \beta_{i}+\epsilon_{i,j},
\end{eqnarray}
where $\alpha_{i,j}= {\mbox{E}} \{Y^*(\delta_{i}=0,\delta_{j}=0)\}$ is the expected net outcome between team $T_i$ and $T_j$ in a hypothetical neutral field (i.e., if there is no team taking home field advantage) based on Assumption (A3),
and $\epsilon_{i,j}$ is random noise with $\mathcal{N}(0,\sigma^2_0)$. The factorization in \eqref{main_model_beta} enables us to unravel the  home field advantage at team level by utilizing the pair-wise match design in soccer games that we will discuss shortly. Assumption (A3) is required for \eqref{main_model_beta} so that we can allow independent baseline effect $\alpha_{i,j}$ for any pair of two teams in a hypothetical neutral field. In addition, we adopt a hierarchical model \citep{berry2013bayesian,chu2018bayesian,geng2020mixture} to characterize the relationship between the home field advantage of individual teams and that of the whole league as
\begin{eqnarray}\label{main_model_alpha}
 \beta_{i}\sim \mathcal{N}(\Delta,\sigma^2),
\end{eqnarray}
where $\sigma^2$ describes the variation of home field advantages across different teams.  Based on \eqref{main_model_beta}, oppositely, we have the model for the match between team $T_j$ (with home field advantage) and team $T_i$ (without home field advantage) as
\begin{equation}\label{main_model_op}
    Y_{j,i}=-\alpha_{i,j}+  \beta_{j}+\epsilon_{j,i},
\end{equation}
where the net score without home field advantage satisfies $\alpha_{i,j} = - \alpha_{j,i}$. Thus, combining \eqref{main_model_beta} and \eqref{main_model_op}, we have \begin{equation}\label{hte}
 Y_{i,j} + Y_{j,i}= \beta_{i}+\beta_{j}+\epsilon_{i,j}+\epsilon_{j,i}.
\end{equation}
In other words, $\alpha_{i,j}$'s are treated as nuisance parameters and do not need to be estimated in our model. By repeating \eqref{hte} over all paired matches between team $T_i$ and team $T_j$,
for $i\not = j$ and $i,j \in \{1,\cdots,n\}$, we have
\begin{equation}\label{est3}
\underbrace{\begin{bmatrix}
1 & 1 & 0&\cdots &0&0&0\\
1 & 0 & 1&\cdots &0&0&0\\
\vdots & \vdots  & \vdots &\cdots &\vdots &\vdots &\vdots \\
0 & 0 & 0&\cdots &1&0&1\\
0 & 0 & 0&\cdots &0&1&1\\
\end{bmatrix}}_{\bm{H}_{N/2\times n}}
\underbrace{\begin{pmatrix}
  {\beta}_1\\
  {\beta}_2\\
 {\beta}_3\\
\vdots\\
 {\beta}_{n-2}\\
 {\beta}_{n-1}\\
 {\beta}_n
\end{pmatrix}}_{\widehat{\bm{\beta}}_{n\times 1}}
+ \underbrace{\begin{pmatrix}
\epsilon_{1,2} + \epsilon_{2,1}\\
\epsilon_{1,3} + \epsilon_{3,1}\\
\vdots\\
\epsilon_{(n-2),n} + \epsilon_{n,(n-2)}\\
\epsilon_{(n-1),n} + \epsilon_{n,(n-1)}
\end{pmatrix}}_{\bm{\epsilon}_{N/2\times 1}}
=
\underbrace{\begin{pmatrix}
Y_{1,2} + Y_{2,1}\\
Y_{1,3} + Y_{3,1}\\
\vdots\\
Y_{(n-2),n} + Y_{n,(n-2)}\\
Y_{(n-1),n} + Y_{n,(n-1)}
\end{pmatrix}}_{\bm{Y}_{N/2\times 1}}.
\end{equation}
This motivates us to estimate the team-specific home advantage effects through the above linear equation system. Specifically, since $\bm{H}$ is a full rank matrix under $n\geq 3$, we use the following estimate of $\bm{\beta}=[\beta_1, \cdots,\beta_n]^\top$,
\begin{equation*}
\widehat{\bm{\beta}} = (\bm{H}^\top \bm{H})^{-1} \bm{H} \bm{Y}.
\end{equation*}
Based on the definition that $\Delta=\sum_i \beta_i /n$, we have an estimator of  $\Delta$ as
\begin{equation}\label{est1}
    \widehat{\Delta} =  \sum_i \widehat{\beta}_i /n,
\end{equation}
where $\widehat{\beta}_i$ is the $i$-th element in $\widehat{\bm{\beta}} $. Following the standard theory for linear regression, we establish the normality of $\widehat{\bm{\beta}} $ in the following proposition.

\begin{prop}
\label{beta_variance}
Assuming the noise terms $\epsilon_{i,j}$ are independent and identically distributed (i.i.d.) Gaussian random variables, i.e., $\epsilon_{i,j}\sim \mathcal{N} (0,\sigma_0^2)$. Under (A1)-(A3) with $n\geq3$, we have
\begin{equation*}
      \widehat{\bm{\beta}} - \bm{\beta}   \sim \mathcal{N}_{n}(\bm{0},\Sigma_{\bm{\beta}}),
\end{equation*}
where $\bm{\beta}=[\beta_1,\cdots,\beta_n]^\top$ denotes the true causal effect for $n$ teams. The covariance matrix $\Sigma_{\bm{\beta}}$ can be estimated by $\widehat{\sigma}^2_{\beta}(\bm{H}^{\top}\bm{H})^{-1}$, where $\widehat{\sigma}^2_{\beta}=2||\bm{H}\widehat{\bm{\beta}}-\bm{Y}||^2_{2}/N$.
\end{prop}
Proposition \ref{beta_variance} holds for every $n \geq 3$ since the normality of $\widehat{\beta}$ is exact. A two-sided $(1-\alpha)$ marginal confidence band  of $\bm{\beta}$ can be obtained as
\begin{equation}
    \widehat{\bm{\beta}}\pm z_{\alpha/2}\sqrt{\text{diag}\{\widehat{\sigma}^2_{\beta}(\bm{H}^{\top}\bm{H})^{-1})\}},
\end{equation}
where $z_{\alpha/2}$ denotes the upper $\alpha/2-$th quantile of a standard normal distribution,
using the estimation of $\Sigma_{\bm{\beta}}$ from Proposition \ref{beta_variance}. Since $\widehat{\Delta}$ is the sample mean of $\widehat{\beta}_i$ as indicated in \eqref{est1}, the normality of the estimated average treatment effect can be obtained immediately from Proposition \ref{beta_variance}
\begin{prop}
\label{prop:sigma_delta}
Suppose that the same set of assumptpions in Proposition \ref{beta_variance} hold. Let $\Delta = \sum_i \beta_i/n$. Then
\begin{equation*}
    \sqrt {n}\{\widehat{\Delta} -\Delta\}  \sim  \mathcal{N}(0,\sigma^2),
\end{equation*}
where an unbiased estimator for the variance $\sigma^2$ is \begin{equation}\label{eq:est}
\widehat{\sigma}^2 ={1\over n -1} \sum_{i=1}^n (\widehat{\beta}_i - \widehat{\Delta})^2 -  \sum_{i=1}^n \text{diag}\{\widehat{\sigma}^2_{\beta}(\bm{H}^{\top}\bm{H})^{-1})\}_i /n
\end{equation}
and $\text{diag}(W)_i$ is the $i$-th diagonal element of matrix $W$.
\end{prop}
The first part of Proposition \ref{prop:sigma_delta} is obvious since $\widehat{\Delta}$ is a linear combination of $\widehat{\bm{\beta}}$, which follows a multivariate normal distribution. However, the estimation of $\sigma^2$ is non-trivial because the off-diagonal entries in the covariance estimator for $\Sigma_{\beta}$ do not perform well due to the high-dimensionality, i.e., number of the free parameters in $\Sigma_{\beta}$
is $O(n^2)$ and our data have a sample size of the same order. Therefore, the usual quadratic form estimator for $\sigma^2$, $n^{-1} \widehat{\sigma}^2_{\beta} \bm{1}^{\top} (\bm{H}^{\top}\bm{H})^{-1})\bm{1}$, is biased, where $\bm{1}$ is the vector of $n$ ones. To solve this problem, the key observation is to use the law of total variance as suggested in the hierarchical modeling literature \citep{berry2013bayesian,chu2018bayesian,geng2020mixture} by considering
\begin{equation*}
 \sigma^2=\text{Var}\{{\mbox{E}}_i(\widehat{\beta_i}|\bm{O})\}= \text{Var}_i(\widehat{\beta_i})-{\mbox{E}}_i\{\text{Var}(\widehat{\beta_i}|\bm{O})\},
\end{equation*}
where $\bm{O}$ is the observed data. The first term $\text{Var}_i(\widehat{\beta_i})$ can be estimated by the sample variance $  \sum_i(\widehat{\beta}_i - \widehat{\Delta})^2/(n -1)$, and the second term can be estimated by $\sum_i\text{diag}\{\widehat{\sigma}^2_{\beta}(\bm{H}^{\top}\bm{H})^{-1})\}_i /n$. Both estimators are unbiased. Hence \eqref{eq:est} holds.
The two-sided $(1-\alpha)$ confidence interval of $\Delta$ thus can be constructed as
\begin{equation}
    \widehat{\Delta}\pm z_{\alpha/2} \sqrt{\widehat{\sigma}^2/n},
\end{equation}
based on Proposition \ref{prop:sigma_delta}.
\begin{remark}
The normality of $\beta_i$ in \eqref{main_model_alpha} can be relaxed to other distributions with mean $\Delta$ and variance $\sigma^2$, which leads to an asymptomatic normality of $\widehat{\Delta}$ as $\sqrt {n}\{\widehat{\Delta} -\Delta\}  \stackrel{d}{\rightarrow} \mathcal{N}(0,\sigma^2)$ when $n\to \infty$, by the central limit theorem. Yet, in reality, we have finite and usually a fairly small number of teams in one league, such as $n=20$ for EPL. Therefore we choose to keep the normality assumption.
\end{remark}
\section{Simulation}\label{sec:simu}
\subsection{Simulation Setup}
In this section, simulation studies are conducted to evaluate the empirical performance of estimators for both $\bm{\beta}$'s and $\Delta$. We consider two data generation scenarios for $\alpha_{ij}$ in our simulation. In the first scenario, we generate ability difference between team $T_i$ and team $T_j$ as $\alpha_{ij}\overset{\text{i.i.d}}{\sim} \mathcal{N}(0,2^2), i=1,\ldots,n,j=1,\ldots,n$ and $\alpha_{ij}=-\alpha_{ji}$. In second scenario, we firstly generate ability of each team  $\text{Ab}_{i}\sim \mathcal{N}(0,2^2),i=1,\ldots,n$, and then set $\alpha_{ij}=\text{Ab}_{i}-\text{Ab}_{j}$ for $i=1,\ldots,n,\, j=1,\ldots,n$. We fix $\Delta=1$ for all the simulation studies, $\beta_i\overset{\text{i.i.d}}{\sim} \mathcal{N}(1,0.3^2)$, and $\epsilon_{ij}\overset{\text{i.i.d}}{\sim} \mathcal{N}(0,\sigma^2_0), i=1,\ldots,n,\, j=1,\ldots,n$. We choose three different values for the variance of $\epsilon_{ij}$, $\sigma^2_0=0.5,1,2$, and four different values for the number of teams, $n=10,20,40,80$. In total, we generate $1,000$ independent replicates for each scenario. The performance of the estimates are evaluated by the bias, the coverage probability (CP), the sample variance of the 1,000 estimates (SV), and the mean of the 1,000 variance estimates (MV) in the following way; take $\Delta$ as an example:
\begin{equation*}
\begin{split}
&\text{Bias} = \frac{1}{1000} \sum_{i=1}^{1000}
\widehat{\Delta}^{(i)} - \Delta^{(0)} , \\
&\text{CP} =
\frac{1}{1000}1\left\{\Delta^{(0)}\in (\Delta^{(i)}_{\alpha/2},\Delta^{(i)}_{1-\alpha/2})\right\},\\
&\text{SV} =  \frac{1}{1000-1}
\sum_{i=1}^{1000}
 \left(\widehat{\Delta}^{(i)} -\bar{\widehat{\Delta}} \right)^2,\\
 &\text{MV} =  \frac{1}{1000}
\sum_{i=1}^{1000}
\widehat{\text{Var}}(\Delta^{(i)}) ,
   \end{split}
\end{equation*}
where $\Delta^{(0)}$ is the true value for $\Delta$, $\widehat{\Delta}^{(i)}$, $ (\Delta^{(i)}_{\alpha/2},\Delta^{(i)}_{1-\alpha/2})$, and $\widehat{\text{Var}}(\Delta^{(i)})$ are the point estimate, $(1-\alpha)$-level confidence bands, and variance estimate from the $i$th replicated simulation data, respectively. We set $\alpha = 0.05$ throughout the simulation.
\subsection{Simulation Results}
We first examine the simulation results for $\Delta$, as well as the
variance estimator of $\Delta$. Table \ref{tab:delta_results} summarizes the performance of the $\Delta$ estimator in different data generation scenarios, and reports the bias, coverage probability, sample variance
of the 1,000 estimates and mean of the 1,000 variance estimates. From this table, it is
clear that our proposed estimator for $\Delta$ has a very small bias that decreases quickly as $n$ increases. The coverage probabilities are also very close to the desired $95\%$ nominal coverage level, especially when $n \geq 20$. Note that the standard error of the estimated coverage probabilities is $\sqrt{.05\times.95/1,000}=0.0069$. The last two columns also confirm the accurate variance estimation for $\widehat{\Delta}$ under all the simulation scenarios; and the estimation accuracy becomes better as $n$ increases.

\begin{table}[!thp]

    \centering
     \caption{Summary of simulation results for $\widehat{\Delta}$. (Column (1): mean bias, (2): coverage probabilities for estimators, (3): sample variance of the 1,000 estimates, and (4): mean of the 1,000 variance estimates).}\label{tab:delta_results}
    \begin{tabular}{cclcccc}
    \toprule
      &  & & (1)& (2) & (3)& (4)\\
     \midrule
  Scenario 1 &	$n=10$& $\sigma_0^2=0.5$&-0.0016&0.899&0.0143&0.0143\\
    &     & $\sigma_0^2=1$&-0.0042&0.913&0.0195&0.0198\\
    &            & $\sigma_0^2=2$&0.0073&0.893&0.0323&0.0345\\
                   	  \midrule
   &	$n=20$& $\sigma_0^2=0.5$&-0.0006&0.923&0.0061&0.0055\\
  &       & $\sigma_0^2=1$&-0.0059&0.914&0.0075&0.0073\\
    &            & $\sigma_0^2=2$&0.0055&0.950&0.0092&0.0106\\
                    \midrule
   &	$n=40$& $\sigma_0^2=0.5$&0.0002&0.928&0.0027&0.0025\\
   &      & $\sigma_0^2=1$&0.0023&0.941&0.0030&0.0030\\
   &             & $\sigma_0^2=2$&-0.0002&0.936&0.0039&0.0038\\
  \midrule
      &	$n=80$& $\sigma_0^2=0.5$&-0.0013&0.947&0.0012&0.0012\\
   &      & $\sigma_0^2=1$&0.0004&0.949&0.0013&0.0013\\
   &             & $\sigma_0^2=2$&0.0016&0.947&0.0015&0.0015\\
         \midrule
         Scenario 2 &	$n=10$& $\sigma_0^2=0.5$&0.0036&0.882&0.0149&0.0139\\
    &     & $\sigma_0^2=1$&0.0032&0.910&0.0192&0.0208\\
    &            & $\sigma_0^2=2$&-0.0004&0.913&0.0321&0.0344\\
                    \midrule
   &	$n=20$& $\sigma_0^2=0.5$&0.0011&0.927&0.0057&0.0055\\
  &       & $\sigma_0^2=1$&0.0003&0.928&0.0074&0.0072\\
    &            & $\sigma_0^2=2$&0.0016&0.933&0.0100&0.0107\\
                    \midrule
   &	$n=40$& $\sigma_0^2=0.5$&0.0019&0.949&0.0025&0.0025\\
   &      & $\sigma_0^2=1$&0.0008&0.940&0.0028&0.0029\\
   &             & $\sigma_0^2=2$&0.0027&0.945&0.0036&0.0038\\
    \midrule
      &	$n=80$& $\sigma_0^2=0.5$&-0.0009&0.954&0.0018&0.0018\\
   &      & $\sigma_0^2=1$&0.0004&0.949&0.0013&0.0013\\
   &             & $\sigma_0^2=2$&-0.0010&0.955&0.0014&0.0015\\
   \bottomrule
    \end{tabular}

\end{table}
Next we examine the simulation results for $\bm{\beta}$'s, as well as the
coverage probability of $\bm{\beta}$'s. Figures \ref{fig:bias_beta1} to \ref{fig:bias_beta4} summarize the boxplots of all $\bm{\beta}$s' estimators in different data generation process to show estimation performance. It is
clear from these figures that the estimator is essentially unbiased for all scenarios considered in the simulation. Furthermore, we present the empirical coverage probabilities for the confidence intervals in Figure~\ref{fig:cover_beta}. All the coverage probabilities are fairly close to the nominal level of 0.95, especially when $n \geq 20$. The coverage probability in general becomes more accurate
as the sample size $n$ increases.

In summary, the simulation results confirm the excellent performance (e.g., bias, variance estimate, and coverage probabilities) of our method in terms of estimating both $\Delta$ and $\bm{\beta}$'s. Note that our simulation scenarios also include the situation when $\alpha_{ij}$ are not $i.i.d.$ generated. Our estimator still performs well in this case because our estimation procedure does not require modeling assumptions or estimation of  $\alpha_{ij}$. The simulation results validate the theoretical results, especially the variance formula, in Propositions \ref{beta_variance} and \ref{prop:sigma_delta}.
\begin{figure}[!h]
    \centering
    \includegraphics[width=0.95\textwidth]{fig/box_N10.pdf}

    \caption{Simulation: Bias of $\bm{\beta}$'s for $n=10$.}
    \label{fig:bias_beta1}
\end{figure}
\begin{figure}[!h]
    \centering
    \includegraphics[width=0.95\textwidth]{fig/box_N20.pdf}

    \caption{Simulation: Bias of $\bm{\beta}$'s for $n=20$.}
    \label{fig:bias_beta2}
\end{figure}
\begin{figure}[!h]
    \centering
    \includegraphics[width=0.95\textwidth]{fig/box_N40.pdf}

    \caption{Simulation: Bias of $\bm{\beta}$'s for $n=40$.}
    \label{fig:bias_beta3}
\end{figure}
\begin{figure}[!h]
    \centering
    \includegraphics[width=0.95\textwidth]{fig/box_N80.pdf}

    \caption{Simulation: Bias of $\bm{\beta}$'s for $n=80$.}
    \label{fig:bias_beta4}
\end{figure}

\begin{figure}[!h]
    \centering
    \includegraphics[width=0.95\textwidth]{fig/point.pdf}
    \caption{Simulation: Coverage Probabilities for $\bm{\beta}$'s in different scenarios.}
    \label{fig:cover_beta}
\end{figure}

\section{Analysis of EPL Data}\label{sec:app}
We focus on the eleven selected in-game statistics collected from 20 teams in 2020-2021 English Premier League as described in Section \ref{sec:data}. For each statistic, we treat it as the main outcome and apply the proposed causal inference method to calculate the team-specific home field advantage $\widehat{\beta}_i$ and the league wide home field advantage $\widehat{\Delta}$  to further analyze the implications of the home field advantage on a team-by-team basis and a statistic-specific basis.
Figure~\ref{fig:cover_beta1} shows the 95\% confidence intervals for the $\widehat{\beta}_i$. A $\beta=0$ would indicate no home field advantage. Brighton is clearly a standout team here, with a wide confidence interval across all statistics.

\begin{figure}[!t]
    \centering
    \includegraphics[width=0.95\textwidth]{fig/All_plots_beta.pdf}
    \caption{EPL data analysis: estimated $\bm{\beta}$ for all in-game statistics.}
    \label{fig:cover_beta1}
\end{figure}

In Figure~\ref{fig:cover_beta1}, each of the eleven statistics had a range of 3-6 teams out of 20 that were significant for that particular statistic. While this is clearly not a majority, it does give us some data on the type of statistic that retains the home field advantage, which we will discuss more in the next section. Conversely, we can look at the significance of the $\widehat{\beta}_i$ by team, which will continue to build on to the story. The Figure~\ref{fig:pvalue_beta} shows the p-values of the $\widehat{\beta}_i$'s for each team (x-axis and color) and for each in-game statistic (denoted by shape).

\begin{figure}[!t]
    \centering
    \includegraphics[width=0.95\textwidth]{fig/Pvalues_all_teams.pdf}
    \caption{EPL data analysis: p-values for $\bm{\beta}$ estimators of 20 teams.}
    \label{fig:pvalue_beta}
\end{figure}

While a majority of the statistics are above the 0.05 standard threshold for significance, there are some teams that have a majority of the in-game statistics being significant. We now look at the makeup of those particular teams and begin to make inference on what makes those specific teams susceptible to the phenomenon that is the home field advantage. The teams that had the highest number of significant statistics were Fulham, Brighton, Newcastle United and Wolverhampton Wanderers. This prompted us to ask what these teams had in common that the other teams did not that caused them to reap the benefits of the home field advantage. The most obvious and telling common factor amongst the teams was their records- they all were ranked in the bottom $50\%$ at the end of the season, with Fulham being one of the three relegation teams in the English Premier League that year. This could imply that less performant teams reap higher benefits of home-field advantage than teams that excel. Teams that perform in the top rankings of the English Premier League can dominate their opponents regardless of their match location. Perhaps their talent and skill prevails such that the difference in their statistics when they are home versus when they are away is nearly indiscernible. In contrast, teams that are not as talented or skilled retain every advantage from being at home - from increased confidence due to crowd involvement, to familiarity, and to a lack of travel fatigue.

The $\widehat{\Delta}$ values tell us the estimated home field advantage that a specific statistic awards across the league. Significant $\widehat{\Delta}$ values give us insight to the statistics that have the highest impact on the home field advantage and how much they contribute. Table~\ref{tab:delta_estimate} shows the estimated $\Delta$ values for each of the in-game statistics, their estimated standard deviations, and their p-values for significance.

\begin{table}[!t]
    \centering
        \caption{EPL data analysis: estimated $\Delta$ for eleven summary statistics.}\label{tab:delta_estimate}
    \begin{tabular}{  c  c c  c }
  \toprule
  \textbf{Statistic} & $\widehat{\Delta}$ & $\widehat{\sigma}$ & \textbf{P-value} \\

   \midrule
   Attacks w/ Shot & 1.568 & 0.556 & 0.005 \\
  \hline
  Defence Interceptions & -1.997 & 1.200 & 0.096  \\
  \hline
  Reaching Opponent Box & 1.732 &  0.608 & 0.004\\
  \hline
  Reaching Opponent Half & 3.592 & 0.949 & 0.000\\
  \hline
  Shots Blocked & 0.350 & 0.329 & 0.287\\
  \hline
  Shots from Box & 1.066 & 0.379 & 0.005 \\
  \hline
  Shots from Danger Zone & 0.786 & 0.279 & 0.005\\
  \hline
  Successful Key Passes &  0.489 & 0.237 & 0.039\\
  \hline
  Touches in Box & 2.253 & 1.043 & 0.031 \\
  \hline
  XG & 0.232 & 0.092 & 0.011 \\
  \hline
  Yellow Cards & -0.026 & 0.126 & 0.834 \\
  \bottomrule
\end{tabular}

\end{table}

Seven of the eleven statistics we chose to analyze proved to significantly have an effect on the net increase (or decrease) of a statistic in favor of the home team.
We notice that offense based statistics, such as attacks with shot and reaching opponent box, are significant at $\alpha=0.05$; while defense based statistics, such as shots blocked and defence interceptions, are not significant, as well as the referee based statistic, yellow cards. This could be an indication that teams retain an advantage offensively when they compete at their home field. Going back to the causes of home-field advantage and examining the statistics, it appears that familiarity of the field would have an impact on the outcomes of the offensive statistics. For example, attacks with shot, reaching the opponent's half or box, shots from box or danger zone, and successful key passes are statistics that are based on the players' ability to get into scoring position, which would increase their likelihood of making a goal. The statistic goal itself did not show to hold a significant home field advantage, which allows us to infer that perhaps the number of goals scored by the home team may not be more than that of the away team, but the opportunities that are presented due to the causal factors of the home field advantage phenomenon are significantly greater. That is, the quality of play is better for the home team than the away team.

\section{Discussion}\label{sec:disc}

Our research found significant measurement of the home field advantage in various soccer statistics. The home field advantage resides more heavily in offensive-based statistics than it does in defensive or referee based statistics. This does not illuminate the importance of defense over offense, but rather that the home field advantage phenomenon is not as prominent in the defensive side of soccer. We found no home field advantage exhibited by the officials of the game, a positive indication of unbiasedness in the sport.  We elected to deeply analyze the eleven statistics that we chose based on their $\widehat{\Delta}$ values and $\widehat{\bm{\beta}}$ values. The statistics that showed significance had less to do with the stand-out statistics of soccer goals, free kicks, fouls, etc. These statistics were more based in quality of play. For example, successful key passes and attacks with shots are statistics based on the quality of the action. What we can derive from this is that the home team is not necessarily benefiting in terms of the obvious, but in the details. The home team takes advantage in being familiar enough with their field such that they get more shots off of attacks and passes, and reach the opponent's box more times. While this may not show up in the box score, it nevertheless can give home team an edge.
Further, we discovered that teams that performed poorly in the season and had lower rankings in the English Premier League retained higher home field advantage. Less successful teams would likely be more confident and comfortable playing in their home environment, while the highly successful teams would be confident playing anywhere.

There are several possible directions for further investigations. First, considering the specific factors such as crowd involvement, familiarity of facilities, and travel fatigue that make up the home field advantage denotes an interesting future direction. Furthermore, building multiple comparison procedure and multivariate response model would merit future research from both methodological and applied perspectives. In addition, our current model depends on the normality assumption for responses. Proposing distribution-free estimators as well as inference procedure will broaden the applications in different areas. From application point of view, studying the in-game statistics over the Covid-19 season in comparison to previous years would help to identify the strength of crowd involvement in the home field advantage.

   \bibliographystyle{agsm}