\section{Introduction}}
Label shift \cite{QuioneroCandela2009WhenTA} occurs predominantly 
in classification tasks, in domains like computer vision \citep{barandela2003restricted,xiao2010sun}, medical diagnosis \citep{grzymala2004approach,mac2002problem}, fraud detection \citep{philip1998toward} and others \citep{haixiang2017learning,radivojac2004classification}.
Label shift often stems from, for example, a non-stationary environment and a biased way that the training and test data sets are collected. 


The basic assumption in label shift \citep{lipton2018detecting} is that, while the class prior changes, the conditional distributions of data within a class are maintained in training and testing. 
A well known special case of label shift is learning with imbalanced data  \citep{cao2019learning,zhao2021active} where the training are remarkably imbalanced due to some sampling bias, while the test data have a more balanced prior on the labels, e.g., uniform prior. It is commonly believed that training with imbalanced data can significantly undermine the overall performance of the trained classifiers \citep{japkowicz2002class,mazurowski2008training}. 

Typically, the test label distribution is unknown and many methods have been proposed to estimate the test priors \citep{Storkey09whentraining,chan05estimation}. When the test label distribution has been estimated, learning under label shift \citep{bishop1995nn,elkan2001costsensitive} reduces to the problem of resampling the training data. Common techniques include oversampling the minority class, downsampling the majority class \citep{branco2015survey,he2009learning,weiss2004mining}, and reweighting \citep{huang2016learning,wang2017learning}. For example, \citet{seiffert2009rusboost} integrate downsampling with boosting for classification with imbalanced training data. The performance (measured by the area under the ROC curve) improves compared to base line methods (e.g., adaboost). However, these balancing method comes with several drawbacks. When the data imbalancing is extreme, downsampling incurs significant loss of information, and oversampling can lead to over-fitting \citep{chawla2002smote,cui2019class}.  Reweighting methods tend to make the optimization of deep models difficult \citep{huang2016learning,cui2019class}.

The aforementioned findings primarily focus on the underparametrized regime, where the sample size is much larger than the data dimension. Recently, significant progress have been made in training overparametrized models, such as deep neural networks. Such progress have stimulated empirical and theoretical studies on overparametrized models, whose statistical properties surprisingly challenge the conventional wisdom. For example, the typical U-shaped bias-variance trade-off curve is complemented by the double descent phenomenon observed in various models (see more in the related work section, \citet{belkin2019reconciling,hastie2019surprises,bartlett2020benign,mei2019generalization}).  

In this paper, we study binary classification with label shift. The classifier is taken as Fisher Linear Discriminant Analysis (LDA, \citet{fisher1936use,bishop:2006:PRML}). 
We consider a Gaussian mixture data model, where $x\in \RR^p$
is the feature, $y \in \{0, 1\}$ is the label, and $x \given y$ is Gaussian distributed. Suppose $n$ training data are collected under certain prior with $n_{\ell}$ samples in class $\ell$ such that $n_0 +n_1 =n$. Specifically, data imbalance refers to the case where the majority and minority class priors do not match, i.e., $n_1/n_0 \neq 1$. The test prior on each class is denoted as $\pi_{\ell} = \PP(y=\ell)$ for $\ell =0,1$. We assume the test priors are known and different from the training priors. When the test priors are not known, we can still estimate them from the empirical label distributions in the test data.

{\bf Our contributions}. We provide a theoretical analysis on the performance of LDA under label shift, in both the under- and over-parametrized regime. We explicitly quantify the misclassification error in the proportional limit of $n \rightarrow \infty$ and $p/n_{\ell}  \rightarrow \gamma_{\ell}$ for $\ell=0,1$, where $\gamma_{\ell} > 0$ is a constant.
Our theory shows a peaking phenomenon when the sample size is close to the data dimension.

We demonstrate a {\it phase transition} phenomenon about data imbalance: The misclassification error exhibits different behaviors as the two-class ratio $n_1/n_0$ varies, depending on the value of $\gamma_0$.
In particular, when $\gamma_0$ is fixed and the ratio $n_1/n_0$ increases from $1$, we observe the following three phases:

\begin{itemize}
\setlength\itemsep{0em}

\item In the underparametrized regime (e.g., $\gamma_0 = 0.5$), the misclassification error first decreases then increases as $\gamma_1$ decays, yet the error decrease is marginal.
  
\item In the lightly overparametrized regime (e.g., $\gamma_0 = 2.5$), the misclassification error first increases then decreases as $\gamma_1$ decays. 
  
\item In the overparametrized regime (e.g., $\gamma_0 = 5$), the misclassification error first decreases then increases, and finally decreases again as $\gamma_1$ decays.
\end{itemize}
  
Such a phase transition suggests that LDA trained with imbalanced data can outperform the counterpart trained with reduced balanced data, in certain overparametrized regime.
  
Moreover, we investigate the impact of the $\ell_2$ regularization on the performance of LDA under label shift: The aforementioned phase transition vanishes when the regularization is sufficiently strong. While the phase transition persists when the regularization is weak.

{\bf Related work}. 
 In literature, many methods have been developed to handle classification under label shift. The sampling techniques include the informed upsampling \citep{liu2008exploratory,kubat1997addressing}, synthetic oversampling \citep{chawla2002smote}, cluster-based oversampling \citep{jo2004class}. Cost-sensitive methods \citep{sun2007cost, huang2016learning} use a cost matrix to represent the penalty of
classifying examples from one class to another. Examples are cost-sensitive decision trees \citep{maloof2003learning} and cost-sensitive neural networks \citep{kukar1998cost}. Kernel-based methods are developed in \citep{liu2006boosting,wang2010boosting,hong2007kernel}.
Despite the empirical success, there are limited theories about how the classification results are affected by data imbalance. 

Statistical properties of LDA has been well established in existing works (\citet{anderson1962introduction}, \citet[Section 10.2]{fukunaga2013introduction}, \citet{velilla2005consistency, zollanvari2011analytic, zollanvari2015generalized, sifaou2020high}). LDA with balanced training and test data is studied in \citet{raudys1998expected}  in the overparametrized case, while the assumption is more restrictive than ours, e.g., the feature vector is Gaussian with the identity covariance matrix in \citet{raudys1998expected}. In the asymptotic regime (i.e., $p, n \rightarrow \infty$), \citet{bickel2004some} show that, when $p/n\rightarrow \infty$, LDA tends to random guessing. Later, \citet{wang2018} consider the proportional scenario where $p / n \rightarrow \gamma \in (0, 1)$. Our theory is more general, and covers both $0<\gamma<1$ (underparametrization) and $\gamma>1$ (overparametrization). 
The misclassification error of Regularized LDA is analyzed in \citet{elkhalil2020large}. Our error analysis on LDA  can not be implied from \citet{elkhalil2020large} by taking the limit of the regularization parameter to $0$ since the covariance matrix is not invertible in the overparametrized case.
We note a parallel line of work studying LDA in sparsity constrained high-dimensional binary classification problems \citep{cai2011direct,shao2011sparse,mai2012direct}.


Our theory demonstrates a peaking phenomenon of LDA, which has been recognized in history \citep{hughes1968mean,duin1995small,hua2005optimal,sima2008peaking} and recently for neural networks \citep{belkin2019reconciling}. This phenomenon has been justified for linear regression \citep{hastie2019surprises,bartlett2020benign,belkin2020two,muthukumar2020harmless,bibas2019new}, random feature regression \citep{mei2019generalization}, logistic regression \citep{AMDDHBLC},  max-margin linear classifier \citep{GEMMLC}, and others \citep{xu2019number,derezinski2019exact,nakkiran2019more}. To our knowledge, we are the first to provide a theoretical justification of the peaking phenomenon for LDA under label shift.

The rest of the paper is organized as follows: Section \ref{sec:BinaryClassificationUsingLDA} introduces LDA; Section \ref{sec:PeakingPhenomenonOfLDA} presents an asymptotic analysis of the misclassification error for LDA, and the phase transition phenomenon under data imbalance; Section \ref{sec:RegularizationImpactOnMisclassification} presents the impact of regularization; Section \ref{sec:Experiments} presents real-data experiments; Section \ref{sec:proof} presents a proof of our main results; Section \ref{sec:ConclusionAndDiscussion} discusses binary classification with extremely imbalanced data and in the highly overparametrized regime. We also discuss future directions.

{\bf Notation}: Given a vector $v \in \RR^p$, we denote $\norm{v}_\Sigma=\sqrt{v^\top\Sigma v}$ for a positive definite matrix $\Sigma$. Given a matrix $M$, we denote $M^\dagger$ as its pseudo-inverse. Let $\Phi(\cdot)$ be the CDF of the standard normal distribution. For two random variables $X$ and $Y$, we denote $X\mathop{=}\limits^d Y$ as $X$ and $Y$ having the same distribution. For a sequence of random variables $\{X_n\}$, we denote $X_n\overset{\rm a.s.}{\longrightarrow}X$ as the almost sure convergence.


\save{\section{Binary Classification using LDA}} \label{sec:BinaryClassificationUsingLDA}

Binary classification aims at classifying an input feature $x\in \RR^p$ into two classes labeled by $y \in \{0,1\}$. A linear classifier achieves this goal by predicting the label based on a linear decision boundary in the form of $\beta^\top (x-\alpha) = b$ with $\alpha, \beta \in \RR^p$ and $b\in \RR$.

To be specific, a linear classifier gives the label $y$ of the feature $x$ by
\begin{align}
f^b_{\alpha,\beta}(x) = 
\begin{cases}
    0, \quad \text{if }\beta^\top (x-\alpha)>b,\\
    1, \quad \text{if }\beta^\top (x-\alpha)\leq b.
\end{cases}
\label{eq:linearclassifier}
\end{align}
Given the class priors $\pi_0 + \pi_1 = 1$, we determine $\alpha, \beta$ and $b$ by minimizing the misclassification error defined as
\begin{equation}\label{eq:miscls-err}
\begin{split}
    \mathcal{R} \left(f^b_{\alpha,\beta}\right)  & = \pi_0 \PP\left(f^b_{\alpha,\beta}(x)=1\given y=0\right) 
    \\  &+ \pi_1 \PP\left(f^b_{\alpha,\beta}(x)=0\given y=1\right).
\end{split}
\end{equation}
LDA approaches the binary classification problem by assuming that the conditional distribution of $x$ given label $y=0$ (resp. $y=1$) is Gaussian $\cN(\mu_0, \Sigma)$ (resp. $\cN(\mu_1, \Sigma)$). Accordingly, the optimal classifier in LDA (also known as the Bayes rule) takes 
\begin{equation}
\label{eq:alphabetastar}
\alpha^* = \frac{\mu_0+\mu_1}{2}, ~~\beta^* = \Sigma^{-1} \left(\mu_0-\mu_1\right), ~~ b^* = \ln \frac{\pi_1}{\pi_0}.
\end{equation}
The decision boundary $\beta^*$ coincides with the Fisher linear discriminant rule \citep{fisher1936use}, which maximizes the ratio of between-class variance and within-class variance:
\begin{equation}
  \argmax_{\beta \in \RR^p}~   \frac{ (\beta^\top \mu_1- \beta^\top \mu_0)^2}{\beta^\top  \Sigma \beta}.
  \label{eqmaxfisher}
\end{equation}
The optimal ratio in \eqref{eqmaxfisher} at $\beta^*$ is defined as the Signal-to-Noise Ratio (SNR), i.e., $\textrm{SNR} = \norm{\beta^*}_\Sigma^2 = (\beta^*)^{\top} \Sigma \beta^*$.

In practice, we receive $n_0$ and $n_1$ i.i.d training data points from class $0$ and $1$, respectively. We denote the per-class data points as $\{x^\ell_i\}_{i=1}^{n_\ell}$ for $\ell = 0, 1$. The total number of samples is $n = n_0 + n_1.$ We obtain the empirical Fisher linear discriminant classifier $f_{\hat\alpha,\hat \beta}^{\hat b}$ with
\begin{align}\label{eq:estimator}
\hat{\alpha} = \frac{\hat{\mu}_0 + \hat{\mu}_1}{2}, \quad \hat{\beta} = \hat{\Sigma}^{\dagger}(\hat{\mu}_0 - \hat{\mu}_1), \quad \hat{b} = \ln \frac{n_1}{n_0},
\end{align}
where $\hat{\mu}_0, \hat{\mu}_1$ and $\hat{\Sigma}$ are empirical estimators of $\mu_0, \mu_1, \Sigma$:
\begin{align*}
\hat{\mu}_\ell = \frac{1}{n_\ell} \sum_{i=1}^{n_\ell} x_i^\ell, ~~ \hat{\Sigma} = \frac{1}{n-2} \sum_{\ell = 0, 1}\sum_{i=1}^{n_\ell} (x_i^\ell - \hat{\mu}_\ell)(x_i^\ell - \hat{\mu}_\ell)^\top.
\end{align*}


\save{\section{Phase Transition of LDA under Label Shift}}\label{sec:PeakingPhenomenonOfLDA}


In this section, we present our main results on the misclassification error analysis of LDA, which covers both the under- and over-parametrized regime. 
\save{\subsection{Error Analysis of LDA}}\label{sec:peaklda}
We first introduce a data model for our theoretical analysis.
\begin{assumption}\label{asm:training data and test data} For both the training and the test data, the conditional distribution of $x$ given $y=\ell$ is Gaussian, i.e., $$x \given(y=\ell) \sim \cN(\mu_\ell,\Sigma), \quad\textrm{for} \quad \ell=0,1.$$ The training data $\{x_i^\ell\}_{i=1}^{n_\ell}$ are i.i.d. sampled for class $\ell = 0, 1$, respectively. The test data have priors $\pi_0$, $\pi_1$ such that $\pi_0+\pi_1=1$. \end{assumption}


Assumption \ref{asm:training data and test data} allows arbitrarily imbalanced training data with $n_0 \neq n_1$ where the test label distribution can be different from the training data.
Under Assumption \ref{asm:training data and test data}, we prove an asymptotic behavior of the misclassification error of LDA in the limit of $n,p\to\infty$.

\begin{theorem}
\label{theorem1}
Let $\gamma_0,\gamma_1$ and $\Delta^2$ be positive constants and set $\gamma=\frac{\gamma_0\gamma_1}{\gamma_0+\gamma_1}$.
Under Assumption \ref{asm:training data and test data}, we let $n,p\to\infty$ with
\begin{equation*}
p/n_0 \to \gamma_0,\ p/n_1 \to \gamma_1,\ p/n \to \gamma, \text{ and } \norm{\beta^*}_\Sigma^2\to\Delta^2.
\end{equation*} 
Then the misclassification error $\cR(f_{\widehat\alpha,\widehat\beta}^{\hat b})$ converges to a limit almost surely when $\gamma \in (0,1)\cup (1,\infty)$, i.e.,
\begin{equation}\label{eq:theorem-underparametrized}
    \mathcal{R} (f_{\widehat\alpha, \widehat\beta}^{\hat b}) \overset{\rm a.s.}{\rightarrow} \sum_{\ell=0,1}\pi_\ell \Phi\left(\frac{\frac{g(\gamma_0, \gamma_1, \ell)}{1-\gamma}+(-1)^\ell\ln{\frac{\gamma_0}{\gamma_1}}}{k(\gamma_0, \gamma_1)\frac{1}{(1-\gamma)^{3/2}}}\right)
\end{equation}
for $0 < \gamma < 1$, and 
\begin{equation}\label{eq:theorem-overparametrized}
    \mathcal{R}(f_{\widehat\alpha, \widehat\beta}^{\hat b}) \overset{\rm a.s.}{\rightarrow}\sum_{\ell=0,1}\pi_\ell \Phi\left(\frac{\frac{g(\gamma_0, \gamma_1, \ell)}{\gamma(\gamma-1)}+(-1)^\ell\ln{\frac{\gamma_0}{\gamma_1}}}{k(\gamma_0, \gamma_1)\frac{1}{(\gamma-1)^{3/2}}}\right)
\end{equation}
for $\gamma > 1$, where
\begin{align*} 
g(\gamma_0, \gamma_1, \ell) &= -\frac{1}{2}(\Delta^2 + (-1)^\ell (\gamma_0 - \gamma_1)) \\
 k(\gamma_0, \gamma_1) &= \sqrt{\Delta^2 + \gamma_0 + \gamma_1}.
 \end{align*}
\end{theorem}

Theorem \ref{theorem1} is proved in Section \ref{appendix:pfthm1}.
Theorem \ref{theorem1} demonstrates a peaking phenomenon of LDA. For simplicity, we illustrate this peaking phenomenon for balanced training ($\gamma_0=\gamma_1$) and test data ($\pi_0 = \pi_1 = 0.5$). In this case, the limit of $\mathcal{R}(f_{\widehat\alpha, \widehat\beta}^{\hat b})$ can be simplified to
\begin{equation*}
    \mathcal{R}(f_{\widehat\alpha,\widehat\beta}^{\hat b})\overset{\rm a.s.}{\longrightarrow}    \begin{cases}
        \Phi\left(-\frac{\Delta^2 \sqrt{1-\gamma}}{2\sqrt{\Delta^2 + 4\gamma}}\right), &\ 0<\gamma<1,\\
        \Phi\left(-\frac{\Delta^2 \sqrt{\gamma-1}}{2\gamma \sqrt{\Delta^2 + 4\gamma}}\right), &\ \gamma>1.
    \end{cases}
     \label{eq.R.balance}
\end{equation*}
Figure \ref{fig:limiting_risk} (left panel) shows the asymptotic misclassification error as a function of $\gamma$ with various SNR. 

A peak occurs as $\gamma$ approaches $1$, when the training sample size is approximately equal to the data dimension. In the {\it underparametrized} regime ($0<\gamma<1$), the misclassification error increases with respect to $\gamma$. In the {\it overparametrized} regime ($\gamma>1$), the misclassification error has a local minimum such that the error first decreases and then increases. This peaking phenomenon exists for all levels of SNR. 




 We interpret the peaking phenomenon as an interplay between the conditioning of the covariance matrix and the variance in the statistical estimation of the means and the covariance matrix. With balanced training data, the training model and the test model are the same so there is no model mismatch.
  In the underparametrized regime ($0<\gamma<1$), the covariance matrix is full rank. In this case,  variance dominates the estimation error, and variance decreases as the sample size increases, since a larger number of samples yield better estimations of $\mu_0, \mu_1$, and $\Sigma$.   Therefore, the misclassification error decreases as $\gamma$ decreases. In the overparametrized regime ($\gamma>1$), the covariance matrix is rank-deficient.  According to Bai-Yin theorem \citep{bai1993}, the condition number of $\Sigma^\dagger$ is proportional to $(1-\sqrt{1/\gamma})^{-2}$, which decreases as $\gamma$ increases from $1$. As a result, the misclassification error decreases as $\gamma$ increases from $1$. When $\gamma$ further increases, variance dominates the error due to limited number of samples, so the error increases again. This explains the local minimum when $\gamma>1$.
  
\begin{figure}[ht]
   
    \centering
    \includegraphics[width=0.9\textwidth]{figures_final/figure1.pdf}
    \caption{Misclassification error of LDA as a function of $\gamma$ when $n_1/n_0=1$ (left panel) and $n_1/n_0=8$ (right panel), with balanced test data ($\pi_0 = \pi_1 = 0.5$), for various levels of asymptotic SNR. The solid curve represents the theoretical error given by Theorem \ref{theorem1}, and the scatters denote synthetic data errors, with $n=200,p=\lceil \gamma n \rceil$.}
    \label{fig:limiting_risk}
   
\end{figure}
  More interestingly, under the label shift with imbalanced training data ($\gamma_0 \neq \gamma_1$) and balanced test data ($\pi_0=\pi_1=0.5$), there are trade-offs among three factors: (1) the conditioning of the pseudo-inverse of the covariance matrix; (2) the variance in the statistical estimation of the means and the covariance matrix; (3) an additional model mismatch. The misclassification error exhibits intricate behaviors with multiple local minima in the error curve. See Figure \ref{fig:limiting_risk} (right panel) for an example with $n_1/n_0=8$.
  
\save{\subsection{Phase Transition under Label Shift}}\label{sec:phaselda}
\vspace{-1mm}
A critical question for classification under label shift is: When the class priors vary between the training data and the test data, is it beneficial to correct the training data distribution?

In order to evaluate the overall performance of classifier trained with imbalanced data, we usually consider a special case of label shift with balanced test data ($\pi_0=\pi_1=0.5$). The question above is reduced to whether it is beneficial to downsample the majority class. Theorem \ref{theorem1} suggests an interesting dichotomy to the question above. Specifically, we fix $\gamma_0$ and investigate how the two-class ratio $n_1/n_0 \in [1, 10]$ affects the performance of LDA. We identify three distinct behaviors of LDA depending on the value of $\gamma_0$. We demonstrate the three behaviors in Figure \ref{fig:lda_asym_risk_imbalanced123}:

\noindent$\bullet$ {\it Behavior I:} When $\gamma_0=0.5$, i.e., the class of $\ell = 0$ is underparametrized, the misclassification error first decreases and then increases as a function of $n_1/n_0$; 

\noindent$\bullet$ {\it Behavior II:} When $\gamma_0 =2.5$, i.e., the class of $\ell = 0$ is slightly overparametrized, the misclassification error first increases and then decreases as a function of $n_1/n_0$;


\noindent$\bullet$ {\it Behavior III:} When $\gamma_0 = 5$, i.e., the class of $\ell = 0$ is overparametrized, the misclassification error first decreases and then increases, and finally decreases again, as a function of $n_1/n_0$. 

We obtain rich insights on training with imbalanced data from Figure \ref{fig:lda_asym_risk_imbalanced123}. When the data imbalance is moderate, for example $n_1/n_0 \in [1, 3]$, training with imbalanced data can outperform the counterpart of using reduced balanced data as in {\it Behavior I} and {\it Behavior III}. Nonetheless, the improvement in {\it Behavior I} is only marginal. On the contrary, {\it Behavior II} indicates that downsampling the majority class improves the performance of LDA.

\begin{figure}[h]
   
    \centering
    \includegraphics[width=0.9\textwidth]{figures_final/figure2.pdf}
    \caption{We demonstrate three behaviors of the misclassification error of LDA as a function of $n_1/n_0$ for various $\gamma_0$ and a fixed asymptotic SNR $\Delta^2=9$. The solid curve represents the error based on our theoretical analysis in Theorem \ref{theorem1}. The scatters denote synthetic data errors, with $n_0=40, p=\lceil \gamma_0n_0 \rceil$. Upper left panel: Plots of an underparametrized ($\gamma_0 =0.5$), a slightly overparametrized ($\gamma = 2.5$) and an overparametrized ($\gamma = 5$) case. Upper right panel: Behavior I in the underparametrized cases; Bottom Left panel: Behavior II in the slighly overparametrized cases. Bottom right panel: Behavior III in the overparametrized cases.}
    \label{fig:lda_asym_risk_imbalanced123}
   
\end{figure}

As the data imbalance becomes more significant, for example, $n_1/n_0 \in [8, 10]$, {\it Behavior II} and {\it III} both indicate that downsampling the majority class hurts the performance. Such behavior is expected, since the downsampling incurs severe information loss.

Theorem \ref{theorem1} also characterizes the misclassification error when $n_1/n_0$ is extremely large. We can check that the misclassification error converges to $0.5$ as $n_1/n_0 \rightarrow \infty$, regardless of the value of $\gamma_0$ (see Appendix \ref{appendix:pfphasetransition}). This indicates that extreme label shift renders the trained classifier suffering from the model mismatch. Nonetheless, in such an extreme imbalanced case, the minority group is prone to be outliers, and detection of outliers is also of great interest.

In the sequel, we formally characterize three phases corresponding to the aforementioned different behaviors. We explicitly identify two phase transition knots $\gamma_a$ and $\gamma_b$ (derived in Appendix \ref{appendix:pfphasetransition}):
\begin{equation*}
    \gamma_a = 2 \quad \textrm{and} \quad \gamma_b = \frac{1}{8}(12-\Delta^2+\sqrt{\Delta^4+40\Delta^2+144}).
\end{equation*}
We claim three phases depending on the value of $\gamma_0$.



\noindent$\bullet$ {\bf Phase I}: If $\gamma_0\in(0,\gamma_a)$, the misclassification error has {\it Behavior I};

\noindent$\bullet$ {\bf Phase II}: If $\gamma_0\in(\gamma_a,\gamma_b)$, the misclassification error has {\it Behavior II};

\noindent$\bullet$ {\bf Phase III}: When $\gamma_0\in(\gamma_b,\gamma_c)$ for some $\gamma_c>0$,  the misclassification error has {\it Behavior III}.


We observe that the first transition $\gamma_a$ appears at the exact parametrized case, i.e., $p = n$. The second transition $\gamma_b$ depends on the SNR and is always larger than $\gamma_a$. 

We remark that in {\bf Phase III}, we cut off $\gamma_0$ at some threshold $\gamma_c$. If $\gamma_0$ is extremely large, i.e., the problem is highly overparametried, we can observe a fourth behavior on the misclassification error. In fact, the misclassification error has multiple local maxima, reflecting a complex interaction between the limited information in the training data and the mismatch of the training and test model. We discuss the highly overparametried regime in Section \ref{sec:ConclusionAndDiscussion}.

\save{\section{Regularization Impact on LDA}}\label{sec:RegularizationImpactOnMisclassification}
In machine learning, regularization is commonly used to stabilize the computation and improve the generalization performance. 
In this section, we study regularized LDA \citep{friedman1989regularized,guo2007regularized} and analyze its asymptotic misclassification error. 

\save{\subsection{Error Analysis of Regularized LDA}}\label{sec:rldaerror}
When an $\ell_2$ regularization term is added on $\beta$, we consider the following optimization problem based on \eqref{eqmaxfisher}:
\begin{equation*}
    \argmin_{\beta \in \RR^p}~ -(\beta^\top \mu_1- \beta^\top \mu_0)^2 + \lambda \norm{\beta}_2^2, \ \  \textrm{s.t.} \ \ \beta^\top  \Sigma \beta = 1.
\end{equation*}
This gives rise to an optimal solution $\beta^*_\lambda$. For simplicity, we formalize it in an equivalent form as
\begin{equation}\begin{split}
    \beta^*_\lambda = (\Sigma + \lambda I )^{-1} \big({\mu}_0-{\mu}_1\big).
\end{split}\end{equation}
We denote the regularized Fisher linear discriminant classifier by $f_{\alpha^*,\beta_\lambda^*}^{b^*}$ where $\alpha^*$ and $b^*$ are given in \eqref{eq:alphabetastar}.


The empirical counterpart of the regularized LDA is given by $f_{\hat\alpha,\hat\beta_\lambda}^{\hat b}$, where the empirical parameters $\hat\alpha$, $\hat\beta_\lambda$ and $\hat b$ are computed through $\widehat\mu_0, \widehat\mu_1$ and $\widehat\Sigma$ according to \eqref{eq:estimator} and
$$
   \hat \beta_\lambda = (\widehat\Sigma + \lambda I )^{-1} \big({\widehat\mu}_0-{\widehat\mu}_1\big).
$$

Similar to Theorem \ref{theorem1}, we prove an asymptotic behavior of the misclassification error of the regularized LDA.

\begin{theorem}\label{theorem2regularized}
Let $\gamma_0,\gamma_1$ and $\Delta^2$ be positive constants and set $\gamma=\frac{\gamma_0\gamma_1}{\gamma_0+\gamma_1}$.
Under Assumption \ref{asm:training data and test data}, we let $n,p\to\infty$ with
\begin{equation*}
    p/n_0 \to \gamma_0,\  p/n_1 \to \gamma_1,\  p/n \to \gamma, \text{ and } \norm{\beta^*}_\Sigma^2 \to \Delta^2.
\end{equation*} The misclassification error $\cR(f^{\hat{b}}_{\widehat\alpha, \widehat{ \beta}_\lambda})$
converges to a limit almost surely when $\gamma \in (0,1) \cup (1,\infty)$, i.e.,
\begin{equation*}
\cR(f^{\hat{b}}_{\widehat\alpha, \widehat{ \beta}_\lambda}) \overset{\rm a.s.}{\rightarrow}\sum_{\ell=0,1} \pi_\ell  \Phi\left(\frac{g(\gamma_0, \gamma_1, \ell)m(-\lambda)+(-1)^\ell\ln{\frac{\gamma_0}{\gamma_1}}}{k(\gamma_0, \gamma_1)\sqrt{m'(-\lambda)}}\right),
\end{equation*}
where $g(\gamma_0, \gamma_1, \ell)$ and $k(\gamma_0,\gamma_1)$ are defined in Theorem \ref{theorem1}, and $m(\lambda) = \int 1/(s - \lambda) dF_\lambda(s)$ with $F_\lambda$ denoting the Marchenko-Pastur law.
\end{theorem}

Theorem \ref{theorem2regularized} is proved in Appendix \ref{appendixtheorem2regularized}. 
Theorem \ref{theorem2regularized} implies that regularization has a smoothing effect on the peaking phenomenon. For simplicity, we show this smoothing effect for balanced training ($\gamma_0=\gamma_1$) and test data ($\pi_0 = \pi_1 = 0.5$).

Figure \ref{fig:rlda_asym_risk_lam} shows the asymptotic misclassification error of regularized LDA as a function of $\gamma$. 
When the regularization is weak, e.g., $\lambda = 10^{-4}$, we observe a similar peaking phenomenon as in Figure \ref{fig:limiting_risk}, while the peak is lower than that in Figure \ref{fig:limiting_risk}. Compared to the unregularized classifier $f^{\hat b}_{\hat \alpha, \hat \beta}$, regularization improves the conditioning of the estimated covariance matrix, i.e., $\hat{\Sigma} + \lambda I$ is never singular, which in turn mitigates the performance degradation when $\gamma \approx 1$.

\begin{figure}[htb!]
    \centering
    \includegraphics[width=0.9\textwidth]{figures_final/figure3.pdf}
    \caption{Misclassification error of regularized LDA as a function of $\gamma$ for various levels of asymptotic SNR. The left panel shows the impact of weak regularization $\lambda = 10^{-4}$, and the right panel uses a strong regularization $\lambda = 1$. The solid line represents the theoretical error given by Theorem \ref{theorem2regularized}, and the scatters denote synthetic data errors, with $n=200,p=\lceil \gamma n \rceil$.}
    \label{fig:rlda_asym_risk_lam}
\end{figure}




\begin{figure}[htb!]
\begin{minipage}[c]{\textwidth}
    \centering
        \includegraphics[width=0.45\textwidth]{figures_final/figure4_upperleft.pdf}
        \includegraphics[width=0.45\textwidth]{figures_final/figure5.pdf}
\end{minipage}

    \caption{We demonstrate three behaviors of the misclassification error of regularized LDA as a function of $n_1/n_0$ for various $\gamma_0$ when we use a weak regularization with $\lambda = 10^{-4}$ (left panel). They disappears with a strong regularization $\lambda = 1$ (right panel). We fix the asymptotic SNR $\Delta^2=9$. The scatters denote synthetic data errors, with $n_0=40, p=\lceil \gamma_0n_0 \rceil$.}
    \label{fig:rlda_asym_risk_imbalanced123}
    
\end{figure}



When the regularization is strong, e.g., $\lambda = 1$, the peaking phenomenon disappears in Figure \ref{fig:rlda_asym_risk_lam} (right panel). In this case, the matrix $\hat{\Sigma}+\lambda I$ is always well-conditioned. The error is dominated by the variance in the statistical estimation of the means and the covariance matrix. As a result, the error increases as $\gamma$ increases. We remark that proper regularization greatly reduces the misclassification error when the problem is approximately exactly parametrized ($\gamma \approx 1$).

\save{\subsection{Phase Transition of Regularized LDA}}\label{sec:rldaphase}

In this section, we study the impact of the regularization on the phase transition phenomenon discussed in Section \ref{sec:phaselda}. We discuss the impact of weak and strong regularization separately, as they lead to very different behaviors.



When the regularization is weak, e.g., $\lambda = 10^{-4}$, we observe a similar phase transition phenomenon as in Section \ref{sec:phaselda}. We depict the misclassification error curves as a function of $n_1/n_0$ in Figure \ref{fig:rlda_asym_risk_imbalanced123}.


When the regularization is strong, e.g., $\lambda = 1$, the phase transition phenomenon disappears (See a formal justification in Appendix \ref{appendix:rldaphasetransition}). Figure \ref{fig:rlda_asym_risk_imbalanced123} shows that the asymptotic misclassification error of regularized LDA as a function of $n_1/n_0$ for various $\gamma_0$. We observe that the error curve consistently first decreases and then increases. In this case, the matrix $\hat{\Sigma}+\lambda I$ is always well-conditioned. Therefore the misclassification error is the consequence of the trade-off between two factors: 1) the model mismatch and 2) the variance in the statistical estimation of the means and the covariance matrix.





\save{\section{Real-Data Binary Classification}}\label{sec:Experiments}

We connect our theoretical findings to real-data binary classification tasks. We consider Neural Network (NN) classifiers for the MNIST (CC-BY 3.0) and CIFAR-10 (MIT) datasets \citep{lecun1998gradient, krizhevsky2009learning}. We focus on the overparametried regime, which is the working regime for neural networks.

\textbf{MNIST dataset and NN classifier}. The MNIST dataset consists of handwritten digits of resolution $28 \times 28$. We train a  neural network classifier with one hidden layer to distinguish digits $3$ and $8$. We vary the number of hidden units in $\{32, 64, 128\}$. The activation function is ReLU, i.e., $\sigma(\cdot) = \max\{\cdot, 0\}$. We use Adam \citep{kingma2014adam} for training, with default hyperparameters in Pytorch.

\textbf{CIFAR-10 dataset and NN classifier}. The CIFAR-10 dataset consists of RGB images of resolution $32 \times 32$ from $10$ categories. We pick two similar categories, e.g., horse v.s. deer, for binary classification. We downsample the data to a resolution of $16 \times 16$. We also train an NN classifier with one hidden layer. The number of neurons in the hidden layer is $10$, and the activation function is ReLU. We use momentum SGD for training, with momentum coefficient $0.9$ and learning rate $0.001$.

\begin{figure}[htb!]
    \centering
    \includegraphics[width=0.9\textwidth]{figures_final/figure6.pdf}
    \caption{Misclassification error of neural network classifiers on MNIST and CIFAR-10 for binary classification.}
    \label{fig:realdata}
\end{figure}
In both tasks, during training, we fix the number of samples from one category ($500$ in MNIST and $1000$ in CIFAR-10), and vary the samples in the other category. The total number of training epochs is $50$ for MNIST and $20$ for CIFAR-10. After training, the classifier is tested on a balanced test set, which consists of $900$ samples from each class in MNIST, and $1000$ samples per-class in CIFAR-10. The test error is averaged over $10$ independent runs for MNIST and $5$ independent runs for CIFAR-10 with random seeds.

\textbf{Result}. The misclassification error in MNIST and CIFAR-10 as a function of $n_1/n_0$ is plotted in Figure \ref{fig:realdata}. In both the MNIST and CIFAR-10 experiments, the neural network is overparametrized. Therefore, we expect the misclassification error exhibits {\it Behavior III} in Figure \ref{fig:lda_asym_risk_imbalanced123}. This is corroborated in Figure \ref{fig:realdata} as the misclassification error first decreases and then increases as $n_1/n_0$ grows.

More importantly, these experiments consistently indicate that in the overparametrized regime, downsampling the majority class may hurt the model performance (cf. $n_1/n_0 \approx 3$ in MNIST and $n_1/n_0 \approx 1.5$ in CIFAR-10). Meanwhile, when the data imbalance is relatively severe, downsampling the majority class can be beneficial, due to its mitigation on label shift between training and testing.

 
\save{\section{Proof of Theorem \ref{theorem1}}} \label{sec:proof}

In this Section, we prove our main Theorem \ref{theorem1}. The proof for Theorem \ref{theorem2regularized} about the regularized LDA is similar to the proof of Theorem \ref{theorem1}. The proof of Theorem \ref{theorem2regularized} is given in Appendix \ref{appendixtheorem2regularized}.

\subsection{Lemmas to be used for the proof of Theorem \ref{theorem1}}
\label{subsec:lemma}
Our analysis relies on change of variables to exploit the independence between the sample mean estimator and the covariance estimator.

\begin{lemma}[Independence between sample mean and sample covariance]\label{lemma:independent} 
Let $x_i \overset{\rm i.i.d.}{\sim} \cN(\mu, \Sigma)$ for $i = 1, \dots, n$, be samples from a Gaussian distribution. We denote $\hat{\mu} = \frac{1}{n} \sum_{i=1}^n x_i$ as the estimator of the sample mean, and $\hat{\Sigma} = \frac{1}{n-1} \sum_{i=1}^n (x_i - \hat{\mu})(x_i - \hat{\mu})^\top$ as the estimator of the sample covariance. Then $\hat{\mu}$ and $\hat{\Sigma}$ are independent.
\end{lemma}

We also utilize the asymptotic characterization of the spectrum of Wishart matrix.
\begin{lemma}[Isotropicity of wishart matrix]\label{lemma3:isotropicity}
Assume $\frac{Z^\top Z}{n-2} \sim \mathcal{W}(I_p , n-2)$. For any vector $z \in \RR^p$ independent of $Z$, we have
\begin{align}\label{apeq:wishart}
{z}^\top \left(\frac{Z^\top Z}{n-2} \right)^{\dagger}{z} \overset{d}{=} \frac{1}{p}\norm{z}^2_2 \tr\left(\left(\frac{Z^\top Z}{n-2}\right)^\dagger\right).
\end{align}
\end{lemma}

Lemma \ref{lemma:independent} and \ref{lemma3:isotropicity} are proved in Appendix \ref{appendix:pflma}.
We next present some helper lemmas.




\begin{lemma}\label{lemma:trace}
Given a matrix $Z$ with i.i.d. standard normal distributed entries, we have 
\begin{align}\label{apeq:trace}
\tr\left(\left(Z^\top Z\right)^\dagger\right) =\tr\left(\left(ZZ^\top \right)^\dagger\right).
\end{align}
\end{lemma}


\begin{lemma}[Strong law of large numbers]\label{lemma:SLLN}
Assume $z\sim\mathcal{N}(0,I_p)$ and $\mu_d$ is a non-random p-dimensional vector such that $\norm{\mu_d}_2^2 \overset{\rm a.s.}{\longrightarrow} \Delta^2$ and $n_\ell$ satisfies Binomial distribution $B(n_0+n_1,\pi_\ell)$ then we have
\begin{equation*}\begin{split}
    \frac{1}{\sqrt{p}}\mu_d^\top z \overset{\rm a.s.}{\longrightarrow}0,\quad
    \frac{1}{p}z^\top z \overset{\rm a.s.}{\longrightarrow} 1,\quad
    \frac{n_\ell}{n_0+n_1} \overset{\rm a.s.}{\longrightarrow} \pi_\ell.
\end{split}\end{equation*}
\end{lemma}

\begin{lemma}[Marchenko-Pastur law]\label{lem.MPLaw}
Let $F_{\gamma}$ be the Marchenko-Pastur (MP) law of $0<\gamma<1$. Then for any real number $\zeta\leq 0$, the Stieltjes transform of MP law and its derivative at $\zeta$ are given as
\begin{align}
    m(\zeta)&=\int \frac{1}{s-\zeta}dF_{\gamma}(s)=\frac{1-\gamma-\zeta-\sqrt{(\zeta-\gamma-1)^2-4\gamma}}{2\gamma\zeta}
    \label{eq.MP.m}
\end{align}
and 
\begin{align}
    \frac{d}{d\zeta}m(\zeta)&=\int \frac{1}{(s-\zeta)^2}dF_{\gamma}(s) =\frac{2\gamma(\gamma-1)-\frac{2\gamma\left(\zeta(\gamma+1) -(\gamma-1)^2\right)}{\sqrt{(\zeta-\gamma-1)^2 -4\gamma}} }{4\gamma^2\zeta^2}.
    \label{eq.MP.dm}
\end{align}

In particular, we have
\begin{align}
    &m(0)=\frac{1}{1-\gamma}  \quad \mbox{ and } \quad \frac{d}{d\zeta}m(0)=\frac{1}{(1-\gamma)^3}. \nonumber
\end{align}
\end{lemma}

Lemma \ref{lemma:trace}, \ref{lemma:SLLN} and \ref{lem.MPLaw} are proved in Appendix \ref{subsec:lemmaproof2}.


\save{\subsection{Proof of Theorem \ref{theorem1}}}\label{appendix:pfthm1}

\begin{proof}[Proof of Theorem \ref{theorem1}]

To begin with, we recall the misclassification error of binary classification in \eqref{eq:miscls-err}. Substituting the Fisher linear discriminant classifier $f^{\hat{b}}_{\hat{\alpha}, \hat{\beta}}$ (with $\hat{\alpha}, \hat{\beta}, \hat{b}$ given in \eqref{eq:estimator}) and prior $\pi_0 = \pi_1 = 1/2$ into \eqref{eq:miscls-err}, we derive
\begin{align}
\cR(f^{\hat{b}}_{\hat{\alpha}, \hat{\beta}})& = \frac{1}{2}\PP \left(\hat \beta ^\top (x-\hat\alpha) \leq \ln\frac{n_1}{n_0} \given x\sim \cN(\mu_0,\Sigma)\right) \nonumber \\
& \quad + \frac{1}{2} \PP \left(\hat \beta ^\top (x - \hat \alpha) > \ln \frac{n_1}{n_0} \given x \sim \cN(\mu_1,\Sigma)\right)\nonumber\\
&= \frac{1}{2} \PP \left(\frac{\hat \beta ^ \top (x - \mu_0)}{\| \hat \beta \|_\Sigma} \leq  \frac{\hat \beta ^ \top (\hat \alpha - \mu _0) + \frac{n_1}{n_0}}{\|\hat \beta \|_\Sigma} \given x\sim \cN(\mu_0,\Sigma)\right) \nonumber \\
& \quad + \frac{1}{2} \PP \left(\frac{\hat \beta ^\top(x-\mu_1)}{\| \hat \beta \|_\Sigma} > \frac{\hat \beta ^\top(\hat \alpha - \mu_1)+\ln \frac{n_1}{n_0}}{\|\hat \beta\|_\Sigma}\given x \sim \cN(\mu_1,\Sigma)\right) \nonumber\\
&= \frac{1}{2} \Phi\left(\left[\widehat\beta^\top  (\widehat\alpha - \mu_0)+\ln{\frac{n_1}{n_0}}\right]/\| \widehat\beta \|_\Sigma\right) + \frac{1}{2} \Phi\left(\left[\widehat\beta^\top  (\mu_1 - \widehat\alpha)+\ln{\frac{n_0}{n_1}}\right]/\|\widehat\beta\|_\Sigma\right).\label{apeq:risk}
\end{align}
Therefore, it suffices to find the limits of 
\begin{align*}
q_0 = \left(\widehat\beta^\top  (\widehat\alpha - \mu_0)+\ln\frac{n_1}{n_0}\right)/\| \widehat\beta \|_\Sigma \quad \textrm{and} \quad
q_1 = \left(\widehat\beta^\top  (\mu_1 - \widehat\alpha)+\ln\frac{n_0}{n_1}\right)/\|\widehat\beta\|_\Sigma,
\end{align*}
since the Gaussian CDF $\Phi$ is continuous.

To further aid our analysis, we characterize the distributions of $\hat{\alpha}$, $\hat{\beta}$, and $\hat{\Sigma}$. Specifically, we make the following change of variables:
\begin{equation}\label{apeq:changeofvar}
\begin{split}
\hat\mu_0 &\mathop{=}\limits^d \frac{1}{\sqrt{n_0}}\Sigma^{\frac{1}{2}}z_0 + \mu_0,\\
\hat\mu_1 &\mathop{=}\limits^d\frac{1}{\sqrt{n_1}}\Sigma^{\frac{1}{2}}z_1 + \mu_1,\\
\hat{\Sigma}~ &\mathop{=}\limits^d \frac{1}{n-2}\Sigma^{\frac{1}{2}}Z^\top Z\Sigma^{\frac{1}{2}},
\end{split}
\end{equation}
where $z_0,z_1 \sim \cN(0,I_p)$ and $Z\in\RR^{(n-2)\times p}$ with each element $Z_{i,j}\sim \cN(0,1)$. Note that $z_0, z_1$ and $Z$ are independent with each other and $Z^\top Z$ is a Wishart matrix by Lemma \ref{lemma:independent}.



\noindent $\bullet$ {\bf Case 1. $0 < \gamma < 1$}. We present in detail how to characterize the asymptotic limit of $q_0$, and $q_1$ follows a similar argument. We tackle the numberator and denominator of $q_0$ separately. Using Lemma \ref{lemma:SLLN}, we check
\begin{align*}
\ln \frac{n_1}{n_0} \overset{\rm a.s.}{\longrightarrow} \ln \frac{\pi_1}{\pi_0} = \ln \frac{\gamma_0}{\gamma_1}.
\end{align*}
Therefore, we temporarily omit the threshold term $\ln (n_1/n_0)$ in both $q_0, q_1$ to ease the presentation. By the change of variables formula in \eqref{apeq:changeofvar} and some manipulation, we deduce
\begin{align}
\hat{\beta}^\top (\hat{\alpha} - \mu_0) & = - \frac{1}{2} (\mu_d + \frac{1}{\sqrt{n_0}}z_0 - \frac{1}{\sqrt{n_1}}z_1)^\top \left(\frac{Z^\top Z}{n-2}\right)^\dagger (\mu_d - \frac{1}{\sqrt{n_0}} z_0 - \frac{1}{\sqrt{n_1}} z_1)\nonumber\\
& =\frac{1}{2n_0} z_0^\top \left(\frac{Z^\top Z}{n-2}\right)^\dagger z_0  - \frac{1}{2}(\mu_d - z_1/\sqrt{n_1})^\top \left(\frac{Z^\top Z}{n-2}\right)^\dagger (\mu_d - z_1/\sqrt{n_1}) \label{apeq:q0numeratorstep1}\\
& \overset{d}{=} \underbrace{\left[ -\frac{1}{2}\|\mu_d -\tfrac{1}{\sqrt{n_1}}z_1\|_2^2 + \frac{1}{2n_0} \|z_0\|_2^2 \right]}_{A} \times \underbrace{\frac{1}{p} \tr\left(\frac{Z^\top Z}{n-2}\right)^{\dagger}}_{B}. \label{apeq:q0numeratorstep2}
\end{align}
where $\mu_d = \Sigma^{-1/2}(\mu_0 - \mu_1)$ and the last equality follows from the isotropicity of the Wishart matrix $Z^\top Z$ in lemma \ref{lemma3:isotropicity},
In the sequel, we establish the limits of terms $A$ and $B$, since they are independent.

\noindent {\bf Asymptotic convergence of $A$}. We show
\begin{align}\label{apeq:A}
A \overset{\rm a.s.}{\longrightarrow} -\frac{1}{2}(\Delta^2 + \gamma_1 - \gamma_0).
\end{align}
To see the result above, we expand term $A$ as
\begin{align*}
A=-\frac{1}{2}\norm{\mu_d}^2 + \frac{1}{\sqrt{n_1}}\mu_d^\top z_1 -\frac{1}{2n_1}\norm{z_1}^2 + \frac{1}{2n_0}\norm{z_0}^2.
\end{align*}
Recall that in Lemma \ref{lemma:SLLN}, we verify that $\frac{1}{n_1}\norm{z_1}_2^2 = \frac{p}{n_1} \frac{1}{p} \norm{z_1}_2^2 \overset{\rm a.s.}{\longrightarrow} \gamma_1$ and $\frac{1}{n_0}\norm{z_0}_2^2 \overset{\rm a.s.}{\longrightarrow} \gamma_0$ by the strong law of large numbers again. Using the Borel-Cantelli lemma, we claim $\frac{1}{\sqrt{n_1}} \mu_d^\top z_1 \overset{\rm a.s.}{\longrightarrow} 0$. Invoking the assertion in Theorem \ref{theorem1} that $\norm{\beta^*}_{\Sigma}^2 = \norm{\mu_d}_2^2 = \Delta^2$, we deduce the desired convergence of $A$.

\noindent {\bf Asymptotic convergence of $B$}. We next show
\begin{align}\label{apeq:B}
B \overset{\rm a.s.}{\longrightarrow} \frac{1}{1-\gamma}.
\end{align}
The convergence result above utilizes the Stieljes transformation of the Marchenko-Pastur law. Specifically, by the Bai-Yin theorem \citep{bai1993}, for $\gamma < 1$, 
\begin{equation*}
    \zeta_{\min}\left(\frac{Z^\top Z}{n-2}\right)\geq \frac{1}{2}(1-\sqrt{\gamma})^2,
\end{equation*}
which implies that $Z^\top Z$ is almost surely invertible.
Conditioned on $Z^\top Z$ being invertible, we rewrite $B$ as 
\begin{align*}
B = \frac{1}{p}\tr\left(\frac{Z^\top Z}{n-2}\right)^{-1} = \frac{1}{p} \sum_{i=1}^p \frac{1}{s_i} = \int \frac{1}{s} d F_{\frac{Z^\top Z}{n-2}} (s),
\end{align*}
where $s_i$'s denote the eigenvalues of $Z^\top Z/(n-2)$ and $F_{M}(a) = \frac{1}{p} \sum_{i=1}^p \mathbbm{1}\{\lambda_i(M) \leq a\}$ is the empirical measure of eigenvalues of $M$.

Now apply the Marchenko-Pastur theorem \cite{marchenko-pastur}, which says that $F_{Z^\top Z/(n-2)}$ converges weakly, almost surely, to the Marchenko-Pastur law $F_\gamma$ (depending only on $\gamma$). Invoking the Portmanteau theorem \citep{marchenko-pastur}, weak convergence is equivalent to the convergence in expectation of all bounded functions $h$, that are continuous except on a set of zero probability under the limiting measure. Defining $h(s)=1/s \cdot \mathds{1}\{s\geq a/2\}$, where we abbreviate $a=(1-\sqrt{\gamma})^2$, it follows that as $n,p \to \infty$, almost surely,
    \begin{equation*}\begin{split}
        \int_{a/2}^\infty \frac{1}{s} dF_{\frac{Z^\top Z}{n-2}}(s) \mathop{\to} & \int_{a/2}^\infty \frac{1}{s} dF_\gamma(s)
    \end{split}\label{apeq.dF}
    \end{equation*}
We can remove the lower limit of integration on both sides above; for the right-hand side, this follows since support of the Marchenko-Pastur law $F_\gamma$ is $[a,b]$, where $b = (1+\sqrt{\gamma})$; for the left-hand side, this follows again by the Bai-Yin theorem \cite{bai1993} (which as already stated, implies the smallest eigenvalues of $Z^\top Z/(n-2)$) is almost surely greater than $a/2$ for large enough n). Thus the last display implies that as $n,p\to \infty$, almost surely,
\begin{equation*
B \to \int \frac{1}{s} F_\gamma(s)
\end{equation*}

Thanks to the Stieljes transformation of the Marchenko-Pastur law, we can explicitly compute $\int \frac{1}{s} d F_\gamma(s)$. In particular, the Stieljes transformation of $F_\gamma$ is defined as
\begin{align*} 
m(\zeta) = \int \frac{1}{s - \zeta} d F_\gamma(s) \quad \textrm{with} \quad \zeta \in \CC.
\end{align*}
Taking $\zeta \rightarrow 0$, we obtain $m(0) = \lim_{\zeta \rightarrow 0} m(\zeta) = \frac{1}{1-\gamma}$ by Lemma \ref{lem.MPLaw}. Consequently, we establish $B \overset{\rm a.s.}{\longrightarrow} \frac{1}{1-\gamma}$.

Substituting \eqref{apeq:A} and \eqref{apeq:B} into \eqref{apeq:q0numeratorstep2} yields
\begin{align}\label{apeq:numerator}
\hat{\beta}^\top (\hat{\alpha} - \mu_0) \overset{\rm a.s.}{\longrightarrow} -\frac{1}{2}(\Delta^2 + \gamma_1 - \gamma_0) \frac{1}{1-\gamma}.
\end{align}

Next we consider the denominator $\|\hat{\beta}\|_{\Sigma}$ in $q_0$. Applying the change of variables in \eqref{apeq:changeofvar} and using Lemma \ref{lemma3:isotropicity}, analogous to \eqref{apeq:q0numeratorstep2}, we derive
\begin{align}
\|\widehat\beta \|_\Sigma^2
= & (\widehat{\mu}_0-\widehat{\mu}_1)^\top \widehat\Sigma^\dagger \Sigma \widehat\Sigma^\dagger(\widehat{\mu}_0-\widehat{\mu}_1)\nonumber  \\
= & ({\mu}_d + \frac{1}{\sqrt{n_0}}{z}_0 - \frac{1}{\sqrt{n_1}}{z}_1)^\top  \left(\left(\frac{Z^\top Z}{n-2}\right)^{\dagger}\right)^2 ({\mu}_d + \frac{1}{\sqrt{n_0}}{z}_0 - \frac{1}{\sqrt{n_1}}{z}_1)\nonumber \\
\overset{d}{=}& \underbrace{\norm{\mu_d + \frac{1}{\sqrt{n_0}}z_0 - \frac{1}{\sqrt{n_1}}z_1}_2^2}_{A'} \times \underbrace{\frac{1}{p} \tr\left(\left(\frac{Z^\top Z}{n-2}\right)^\dagger\right)^2}_{B'}. \label{apeq:q0denominator}
\end{align}
The convergence of term $A'$ follows the same argument of term $A$ in the numerator, and we have
\begin{align}\label{apeq:A'}
A' \overset{\rm a.s.}{\longrightarrow} \Delta^2 + \gamma_1 + \gamma_0.
\end{align}
Conditioned on $Z^\top Z$ being invertible, term $B'$ can be written as
\begin{align*}
B' = \int \frac{1}{s^2} d F_{\frac{Z^\top Z}{n-2}}(s) \longrightarrow \int \frac{1}{s^2} dF_\gamma (s).
\end{align*}
To compute the limiting integral above, we differentiate the Stieljes transformation $m(\zeta)$. By sending $\zeta \rightarrow 0$ again, we can derive
\begin{align}\label{apeq:B'}
B' \overset{\rm a.s.}{\longrightarrow} \lim_{\zeta \rightarrow 0} \frac{d}{d\zeta} m(\zeta) = \frac{1}{(1-\zeta)^3}.
\end{align}
Substituting \eqref{apeq:A'} and \eqref{apeq:B'} into \eqref{apeq:q0denominator} yields
\begin{align}\label{apeq:divider}
\|\hat{\beta}\|_{\Sigma} \overset{\rm a.s.}{\longrightarrow} \sqrt{\Delta^2 + \gamma_1 + \gamma_0} \frac{1}{(1-\gamma)^{3/2}}.
\end{align} 

Combining \eqref{apeq:numerator} and \eqref{apeq:divider}, as well as putting the threshold term $\ln (n_1/n_0)$ back, we obtain
\begin{align}\label{apeq:q0}
q_0 \overset{\rm a.s.}{\longrightarrow} \frac{-\frac{\Delta^2 + \gamma_1 - \gamma_0}{2(1-\gamma)}+\ln \frac{\gamma_0}{\gamma_1}}{\sqrt{\Delta^2 + \gamma_1 + \gamma_0} \frac{1}{(1-\gamma)^{3/2}}}.
\end{align}
The same argument of analyzing $q_0$ applies to $q_1$, and therefore, we have
\begin{align}\label{apeq:q1}
q_1 \overset{\rm a.s.}{\longrightarrow} \frac{-\frac{\Delta^2 + \gamma_0 - \gamma_1}{2(1-\gamma)}+\ln \frac{\gamma_1}{\gamma_0}}{\sqrt{\Delta^2 + \gamma_1 + \gamma_0} \frac{1}{(1-\gamma)^{3/2}}}.
\end{align}
To complete the proof in the case of $0 < \gamma < 1$, we plugging \eqref{apeq:q0} and \eqref{apeq:q1} into \eqref{apeq:risk}.

\noindent $\bullet$ {\bf Case 2. $\gamma > 1$}. The goal is still to find the limits of $q_0$ and $q_1$. Consider $q_0$ first. We observe that both \eqref{apeq:q0numeratorstep2} and \eqref{apeq:q0denominator} are valid for $\gamma > 1$. However, a key difference is that $Z^\top Z$ is rank deficient in the limit considering $\gamma > 1$. To resolve this issue, we observe that $Z^\top Z$ and $ZZ^\top$ share all nonzero eigenvalues. Therefore, we can replace $Z^\top Z$ by $ZZ^\top$ in \eqref{apeq:q0numeratorstep2} and \eqref{apeq:q0denominator} without changing their values. Such a reasoning is justified by Lemma \ref{lemma:trace}.
We can now rewrite \eqref{apeq:q0numeratorstep2} as
\begin{align}
\hat{\beta}^\top (\hat{\alpha} - \mu_0) &\overset{d}{=} \left[ -\frac{1}{2}\|\mu_d -\tfrac{1}{\sqrt{n_1}}z_1\|_2^2 + \frac{1}{2n_0} \|z_0\|_2^2 \right] \times \frac{(n-2)^2}{p^2}\frac{1}{n-2}\tr\left(\frac{Z Z^\top}{p}\right)^{\dagger}. \nonumber\\
&\overset{\rm a.s.}{\longrightarrow} -\frac{1}{2}\left[\Delta^2 + \gamma_1 - \gamma_0 \right] \times \frac{1}{\gamma(\gamma - 1)}.\label{apeq:numeratorresult}
\end{align}
Similarly we can rewrite \eqref{apeq:q0denominator} as 
\begin{align}
    \| \hat \beta \|_\Sigma^2 & \mathop{=}\limits^d \norm{ \mu_d + \frac{1}{\sqrt{n_0}} z_0 - \frac{1}{\sqrt{n_1}} z_1 }_2^2 \times \frac{(n-2)^3}{p^3} \frac{1}{n-2} \tr \left( \left(\frac{Z^\top Z}{n-2}\right)^\dagger \right)^2 \nonumber \\
    \label{apeq:q0denominatorresult}&\overset{\rm a.s.}{\longrightarrow} (\Delta^2 + \gamma_0 + \gamma_1) \times \frac{1}{(\gamma-1)^3}. \\
\end{align}
Combining \eqref{apeq:numeratorresult} and \eqref{apeq:q0denominatorresult}, as well as putting the threshold term $\ln(n_1/n_0)$ back, we obtain
\begin{align}
q_0 & \overset{\rm a.s.}{\longrightarrow} \frac{-\frac{\Delta^2 + \gamma_1 - \gamma_0}{2\gamma(\gamma-1)} + \ln \frac{\gamma_0}{\gamma_1}}{\sqrt{\Delta^2 + \gamma_1 + \gamma_0} \frac{1}{(\gamma-1)^{3/2}}}, \nonumber
\end{align}
The same argument of analyzing $q_0$ applies to $q_1$, and therefore, we have
\begin{align}
q_1 & \overset{\rm a.s.}{\longrightarrow} \frac{-\frac{\Delta^2 + \gamma_0 - \gamma_1}{2\gamma(\gamma-1)} + \ln \frac{\gamma_1}{\gamma_0}}{\sqrt{\Delta^2 + \gamma_1 + \gamma_0} \frac{1}{(\gamma-1)^{3/2}}}. \nonumber
\end{align}
The misclassification error in the case of $\gamma > 1$ follows by substituting $q_0, q_1$ above into \eqref{apeq:risk}. The proof is complete.
\end{proof}



\begin{figure}[!tp]
    \begin{center}
        \includegraphics[width=0.5\textwidth]{figures_final/figure7.pdf}
    \end{center}
    \vspace{-0.2cm}
    \caption{Misclassification error in the highly overparametrized regime as a function of $n_1/n_0$, with $\Delta^2 = 9$.}
    \label{fig:highlyoverparametried}
    \vspace{-0.5cm}
\end{figure}


\save{\section{Conclusion and Discussion}
\label{sec:ConclusionAndDiscussion}}

This paper provides a theoretical analysis on the performance of LDA under label shift, in both the under- and over-parametrized regime. We explicitly quantify the misclassification error in the proportional limit of $n \rightarrow \infty$ and $p/n_{\ell}  \rightarrow \gamma_{\ell}$ for $\ell=0,1$, where $\gamma_{\ell} > 0$ is a constant.
Our theory shows a peaking phenomenon when the sample size is close to the data dimension.
We demonstrate a {\it phase transition} phenomenon about data imbalance: The misclassification error exhibits different behaviors as the two-class ratio $n_1/n_0$ varies, depending on the value of $\gamma_0$. We clearly characterize the three behaviors of the misclassification error in the underparametrized, lightly overparametrized, and overparametrized regions depending on $\gamma_0$. We also investigate the regularized LDA, and show that the peaking and phase transition phenomenons disappear when the regularization becomes strong.

Additionally, in the highly overparametrized regime, the misclassification error under label shift has multiple maxima, as shown in Figure \ref{fig:highlyoverparametried}. We believe that the multi-peaking behavior of the misclassification error reflects the intertwined interaction between the lack of information in data and the model mismatch under label shift.







