
\section{Introduction}
\label{sec:introduction}
\input{sections/10-intro}

\section{Related Work}
\label{sec:related}
\input{sections/20-related}

\section{Estimating Risk under Exposure Bias}
\label{sec:estimating-risk}
\input{sections/30-estimating-risk}

\section{Learning Propensities and Link Probabilities}
\label{sec:learning-estimators}
\input{sections/40-learning-estimators}

\section{Feedback Loops}
\label{sec:feedback-loops}
\input{sections/50-feedback-loops}

\section{Experiments}
\label{sec:experiments}
\input{sections/60-experiments}

\section{Conclusion}
\label{sec:conclusion}
\input{sections/70-conclusion}








\subsection{Semi-Synthetic Dataset} \label{sec:experiments-semi-synthetic}

\begin{table}
\caption{Evaluation metrics 
on the test set of the semi-synthetic data computed against
known ground truth citation links.
}
\label{table:semi-synthetic-eval-metrics}
\begin{center}
\begin{small}
\begin{sc}
\begin{tabular}{lcccccr}
\toprule
Model & Prec. & Rec. & AUC & MAP \\
\midrule
\makecell[l]{No Prop.} & 67.24 & 54.81 & 84.45 & 41.87 \\
\makecell[l]{MLE} & 81.04 & 60.19 & 93.12 & 56.77 \\
\makecell[l]{$\widehat{R}_w$} & \textbf{83.28} & 63.73 & \textbf{96.42} & 56.96 \\
\makecell[l]{$\widehat{R}_\text{PU}$} & 82.16 & 63.07 & 94.28 & 58.01 \\
\makecell[l]{$\widehat{R}_\text{AP}$} & 83.01 & \textbf{65.54} & 95.38 & \textbf{59.90} \\
\bottomrule
\end{tabular}
\end{sc}
\end{small}
\end{center}
\end{table}

Since we do not have ground truth exposure values in the MAG dataset, 
we cannot know whether a paper was not cited
due to a lack of exposure or due to irrelevancy.
As a result, we construct a semi-synthetic dataset 
with simulated propensity scores and link probabilities.
We use a subset of $41{,}600$ papers.
We generate train-test-validation splits 
by taking a topological ordering of the nodes 
and use the subgraph created 
from the first $70\%$ for training, 
next $10\%$ for validation, 
and the remaining $20\%$ for testing.
We use the real text and FOS for each paper.
The simulated propensity matrix $\pi$ 
is a $19 \times 19$ matrix 
with its diagonal and off-diagonal entries initialized
from $U(0.7, 1)$ and $U(0.1, 0.3)$, respectively, 
where $U(.)$ is the uniform distribution.
The link probability is simulated using 
$y_{ij} = \sigma(w^\top (h_i \odot h_j) + b)$,
where $\sigma$ is the sigmoid function, 
$\odot$ is element-wise product, 
and $w, b$ are fixed \emph{known} vectors.

\begin{table}
\caption{RMSE of the estimated risk 
with respect to the true risk 
computed using our proposed estimators.
The first column shows the risk 
used in the loss function in Eq.~\ref{eq:combined-loss-function} 
to learn $\widehat{\pi}$ and $\widehat{y}$.
}
\label{table:semi-synthetic-rmse-metrics}
\begin{center}
\begin{small}
\begin{sc}
\begin{tabular}{l|ccccr}
\toprule
\multirow{2}{*}{\makecell{Trained \\ Using}} & \multicolumn{4}{c}{Estimator Used} \\
& $\widehat{R}_\text{naive}$ & $\widehat{R}_w$ & $\widehat{R}_\text{PU}$ & $\widehat{R}_\text{AP}$ \\
\midrule
No Prop. & 1.50 & - & - & - \\
\makecell[l]{MLE} & 0.67 & \textbf{0.23} & 0.24 & 0.32 \\
\makecell[l]{$\widehat{R}_w$} & 0.43 & \textbf{0.04} & 0.10 & 0.11 \\
\makecell[l]{$\widehat{R}_\text{PU}$} & 0.38 & 0.05 & 0.11 & \textbf{0.04} \\
\makecell[l]{$\widehat{R}_\text{AP}$} & 0.41 & 0.06 & 0.08 & \textbf{0.03} \\
\bottomrule
\end{tabular}
\end{sc}
\end{small}
\end{center}
\end{table}

We show 
the evaluation metrics for five models on the test set 
computed against the simulated true citations 
(not the observed citations) (Table \ref{table:semi-synthetic-eval-metrics}).
\textit{No Prop} is the model trained naively on the observed data using only
the output model in Eq. \ref{eq:experiments-output-model}. 
\textit{MLE} is the model trained using the loss function 
in Eq. \ref{eq:combined-loss-function} with $\lambda_R = 0$.
The remaining three are models trained using
$\widehat{R}_w, \widehat{R}_{\text{PU}}$, 
and $\widehat{R}_{\text{AP}}$ with $\lambda_R = 10$ and $\lambda_L = 1$.
We see that all other estimators 
significantly outperform \textit{No Prop}.
Additionally, our proposed estimators 
lead to improved performance over the MLE. 
We emphasize that these metrics 
are computed against true citations 
and thus are a measure of true risk 
which is the appropriate metric 
for evaluating an RS's performance.
This shows the utility of accounting for exposure bias 
and learning using our proposed loss function.

\begin{figure}[t]
\begin{center}
\begin{subfigure}{0.47\columnwidth}
\includegraphics[width=0.9\columnwidth]{figures/ground_truth_propensities.pdf}
\caption{Ground Truth}
\end{subfigure}
\begin{subfigure}{0.47\columnwidth}
\includegraphics[width=0.9\columnwidth]{figures/estimated_propensities.pdf}
\caption{Estimated}
\end{subfigure}
\caption{The estimated propensities propensities are close to the true simulated values when learned using $\widehat{R}_w$.}
\label{fig:ground-truth-and-estimated-propensities}
\end{center}
\end{figure}

In this work, we tackle two separate (but related) challenges.
The first challenge is \textit{learning} link probabilities in such a way that they are not underestimated due to exposure bias.
The second challenge is \textit{evaluating} a RS given learned link probabilities and propensity scores, i.e., computing a good estimate of the true risk.
We demonstrate the efficacy of our methods for the second challenge and show
that our proposed weighting schemes lead to good estimates of the true risk (Table~\ref{table:semi-synthetic-rmse-metrics}).
We show 
the RMSE of the risk estimated 
using the proposed estimators
with respect to the true risk. 
The first column denotes risk function 
used to train the model 
(as described in Section~\ref{sec:learning-estimators}).
The rest of the columns denote 
the estimators used to estimate true risk 
using the learned propensities and link probabilities 
from the trained model in the first column.
We trained each model $10$ times to compute the RMSE.
The RMSE estimated using $\widehat{R}_{\text{naive}}$ 
is always greater than that of the other estimators, 
which shows that leveraging the learned propensities 
leads to substantially better estimates of the true risk 
(and thus more accurately evaluates the RS). 
The RMSE when trained using the \textit{MLE}
is higher than when trained using our proposed estimators, 
showing the benefit of our proposed estimators over the \textit{MLE}.
This also qualitatively validates the generalization bound
proved in Section~\ref{sec:learning-estimators} 
by showing that minimizing 
$\widehat{R} \in \{ \widehat{R}_w, \widehat{R}_{\text{PU}}, \widehat{R}_{\text{AP}} \}$
also leads to small values of the true risk.

The heatmap of simulated propensities 
and estimated propensities when using $\widehat{R}_w$
shows that the estimated propensities 
are close to the true propensities (Figure \ref{fig:ground-truth-and-estimated-propensities}). 
The mean relative error between the true 
and estimated propensities is $19.47\%$, 
demonstrating that the training procedure recovers the propensities.
Together, these results show that our methods 
successfully mitigate exposure bias in this dataset.

\subsection{Real-World Datasets}\label{sec:real-world-datasets}

\begin{table}
\caption{Evaluation metrics for various models 
computed on the test sets 
of the two real-world citation datasets.}
\label{table:real-data-standard-metrics}
\begin{center}
\begin{small}
\begin{sc}
\begin{tabular}{lccccccr}
\toprule
Model & Prec. & Rec. & F1 & AUC & \makecell[l]{MAP}\\
\midrule
\multicolumn{6}{l}{\textbf{Dataset 1}} \\
\midrule
\makecell[l]{No Prop.} & 29.45 & 78.30 & 42.81 & 84.44 & 24.10 \\
\makecell[l]{MLE} & 30.24 & 77.84 & 43.56 & 84.41 & 24.60 \\
\makecell[l]{$\widehat{R}_w$} & 31.46 & 78.02 & 44.84 & 84.74 & 25.60  \\
\makecell[l]{$\widehat{R}_\text{PU}$} & 30.98 & \textbf{78.94} & 44.49 & \textbf{85.24} & 25.11 \\
\makecell[l]{$\widehat{R}_\text{AP}$} & \textbf{36.07} & 76.08 & \textbf{48.94} & 84.67 & \textbf{28.58} \\
\midrule
\multicolumn{6}{l}{\textbf{Dataset 2}} \\
\midrule
\makecell[l]{No Prop.} & 44.86 & 70.85 & 54.94 & 83.22 & 33.19 \\
\makecell[l]{MLE} & 44.43 & 74.66 & 55.71 & 84.97 & 34.39 \\
\makecell[l]{$\widehat{R}_w$} & \textbf{48.70} & 71.62 & \textbf{57.98} & 83.90 & \textbf{36.25} \\
\makecell[l]{$\widehat{R}_\text{PU}$} & 42.17 & \textbf{76.15} & 54.28 & \textbf{85.43} & 33.26 \\
\makecell[l]{$\widehat{R}_\text{AP}$} & 47.22 & 71.84 & 56.98 & 83.89 & 35.27 \\
\bottomrule
\end{tabular}
\end{sc}
\end{small}
\end{center}
\end{table}
We now evaluate our proposed method 
on a real-world citation network. 
We construct two datasets
by using disjoint subgraphs of the MAG. 
The first generated dataset has
$2{,}442{,}008$ papers and $7{,}577{,}886$ edges. 
The second dataset has $1{,}328{,}664$ papers and $1{,}469{,}899$ edges. 
Thus the second graph is sparser than the first one.
The FOS distribution is also different in both datasets 
(see details in Appendix \ref{appendix:experiments}). 
We use 70-10-20\% train-validation-test splits 
generated similarly to the semi-synthetic dataset.
We do not have access to true exposure values
and thus we evaluate our methods against the observed citation links.

We show 
the evaluation metrics for the proposed estimators (Table \ref{table:real-data-standard-metrics}).
Since we do not have access 
to the true citation links in the real dataset,
we compute these metrics over the observed links. 
In other words, this is a measure of the naive risk.
We see that our proposed estimators 
achieve comparable metrics 
to \textit{No Prop}.
For both datasets, the best numbers for each metric 
are achieved by estimators other than \textit{No Prop}.
Moreover, the models using the weighted estimators 
outperform the \textit{MLE} estimator in both datasets.
These results show that our proposed estimators 
achieve comparable performance 
even when evaluated on the observed citation data.
Similarly, Table \ref{table:real-data-link-pred-metrics} shows 
link prediction metrics
for various models computed against \textit{observed} citations.
\textit{Recall@100} refers to the recall in the top 100 recommendations 
averaged across all papers in the test set. 
\textit{Mean Rank} is the mean rank 
of the cited papers
averaged across all the papers.
\textit{Entropy@100 of True Positives} is the entropy in the FOS 
of the true positives in the top 100 recommendations for each paper; 
we use it to measure the diversity in the FOS of the recommendations.
Our proposed estimators achieve comparable \textit{Recall@100} and \textit{Mean Rank} to \textit{No Prop}
for both datasets.
As expected, propensity based estimators have higher FOS entropy scores than \textit{No Prop}, with $\widehat{R}_w$ achieving the highest FOS entropy in both datasets. Thus our proposed estimators \emph{recommend more relevant papers from different FOS and still maintain comparable performance} to \textit{No Prop}.

At first blush, the comparable performance to \textit{No Prop} may not seem compelling. However, this is a strong result.
Our goal is to correct exposure bias and minimize \textit{true risk}, not observed (or naive) risk.
Since Tables \ref{table:real-data-standard-metrics} and \ref{table:real-data-link-pred-metrics} are computed against observed citations, our proposed methods should not be expected to outperform \textit{No Prop} as they are not trying to optimize metrics against the observed links.
In Section \ref{sec:experiments-semi-synthetic}, we showed that our methods correct exposure bias and achieve lower true risk. Coupled with those results, our goal in this section was to show that our methods do not lower performance even if evaluated against the standard evaluation metrics.
\textcolor{black}{We suspect that the negative log-likelihood term $\mathcal{L}(o|\widehat{\pi}, \widehat{y})$ in Eq.~\ref{eq:combined-loss-function} is likely responsible for the comparable performance against the observed risk. This is because, as seen from the results, the MLE also performs well as compared to \textit{No Prop}.}.

\begin{table}
\caption{Link prediction metrics for various models when evaluated on the test set of the real datasets.}
\label{table:real-data-link-pred-metrics}
\begin{center}
\begin{small}
\begin{sc}
\begin{tabular}{lcccccr}
\toprule
Model & \makecell{Recall\\@100} & \makecell{Mean \\ Rank} & \makecell{Entropy@100 \\ True positives}\\
\midrule
\multicolumn{4}{l}{\textbf{Dataset 1}} \\
\midrule
\makecell[l]{No Prop.} & 24.39 & \textbf{2247.27} & 1.65 \\
\makecell[l]{MLE} & 24.70 & 2891.40 & 1.73 \\
\makecell[l]{$\widehat{R}_w$} & 25.03 & 2836.73 & \textbf{1.74} \\
\makecell[l]{$\widehat{R}_\text{PU}$} & 24.61 & 2875.13 & 1.73 \\
\makecell[l]{$\widehat{R}_\text{AP}$} & \textbf{26.66} & 2425.51 & 1.71 \\
\midrule
\multicolumn{4}{l}{\textbf{Dataset 2}} \\
\midrule
\makecell[l]{No Prop.} & \textbf{6.32} & \textbf{10170.26} & 1.06 \\
\makecell[l]{MLE} & 6.07 & 10731.88 & 1.08 \\
\makecell[l]{$\widehat{R}_w$} & 6.30 & 10873.19 & \textbf{1.12} \\
\makecell[l]{$\widehat{R}_\text{PU}$} & 5.92 & 10717.06 & 1.08 \\
\makecell[l]{$\widehat{R}_\text{AP}$} & 5.99 & 10801.11 & 1.10 \\
\bottomrule
\end{tabular}
\end{sc}
\end{small}
\end{center}
\end{table}

\subsection{Feedback Loops}

\begin{figure}[t]
\begin{center}
\begin{subfigure}{0.49\columnwidth}
\includegraphics[width=1\columnwidth]{figures/feedback_loops_no_propensity.pdf}
\caption{No Propensity}
\label{fig:feedback-loops-fos-no-propensity}
\end{subfigure}
\begin{subfigure}{0.49\columnwidth}
\includegraphics[width=1\columnwidth]{figures/feedback_loops_with_propensity.pdf}
\caption{With Propensity ($\widehat{R}_w$)}
\label{fig:feedback-loops-fos-with-propensity}
\end{subfigure}
\caption{The fraction of recommended papers from the same FOS over time.}
\label{fig:feedback-loops-fos}
\end{center}
\end{figure}

We run simulations to examine what happens 
when a citation recommender is trained repeatedly 
on data collected from users
interacting with its recommendations. 
We use the iterative training procedure 
described in Section~\ref{sec:feedback-loops}. 
We construct a training set of $410$ papers from the MAG,
with their corresponding real FOS and text embeddings. 
For ease of exposition, we use an arbitrary 
but fixed mapping to map the 19 FOS to two FOS. 
The synthetic propensities and link probabilities
are simulated similarly to Sec~\ref{sec:experiments-semi-synthetic}.
In the first iteration, the models are trained 
using the observed citation network. 
For subsequent training iterations, 
the training data is generated as follows.
For each paper, we recommend 20 papers. 
The probability of recommending a paper is proportional 
to its estimated citation probability.
We then simulate the user's interaction with the recommendations 
according to the known simulated propensities and citation probabilities.
This generates the training set for the subsequent iteration.
We then repeat this process.

We show how the fraction of recommended papers 
from the same FOS changes over multiple training iterations 
for models trained without propensity, i.e., \textit{No Prop} 
and the model trained using $\widehat{R}_w$ 
(Figure \ref{fig:feedback-loops-fos}).
We plot this time series for both FOS in our dataset.
For \textit{No Prop}, the fraction of papers recommended 
from the same FOS increases over time for both FOS (Figure~\ref{fig:feedback-loops-fos-no-propensity}).
This demonstrates the existence of a feedback loop 
that worsens exposure bias and reduces 
the number of papers recommended 
from a different FOS over time.
On the other hand, when we train our models using $\widehat{R}_w$, 
the feedback loop no longer exists 
and the fraction of papers recommended 
from a different FOS 
remains stable over time
(Figure~\ref{fig:feedback-loops-fos-with-propensity}). 
This shows that our proposed estimator continues 
to recommend relevant papers from a different FOS 
and corrects the feedback loop.

\section{Bias and Variance} \label{appendix:bias-variance}

\begin{lemma}\label{lemma:variance-bernoulli}
Let $X \sim \text{Bernoulli}(\theta)$ and $Y = aX + b(1-X)$, where $a$ and $b$ are some constants. Then
\begin{align*}
    \text{Var}(Y) = \theta (1 - \theta) (a - b)^2.
\end{align*}
\end{lemma}

\textbf{Proof of Lemma \ref{lemma:bias-var-naive}}

\begin{lemma*}
The bias and variance of $\widehat{R}_{\text{naive}}(\widehat{o})$ are
\begin{align*}
    \text{B}(\widehat{R}_{\text{naive}}) &= \left|{\bf E}[\widehat{R}_{\text{naive}}] - R(\widehat{o})\right| \\
        &= \frac{\Delta}{|\mathcal{U}|} \left| \sum_{(i,j) \in \mathcal{U}} y_{ij} (1 - \pi_{ij}) (1 - 2 \widehat{o}_{ij}) \right|,\\
    \text{Var}(\widehat{R}_{\text{naive}}) &= \frac{\Delta^2}{|\mathcal{U}|^2} \sum_{(i,j) \in \mathcal{U}} y_{ij}\pi_{ij} (1 - y_{ij}\pi_{ij}).
\end{align*}
\end{lemma*}
\begin{proof}
We have
\begin{align*}
    \widehat{R}_{\text{naive}}(o, \widehat{y}) &= \frac{1}{|\mathcal{U}|} \sum_{(i,j) \in \mathcal{U}} \delta(o_{ij}, \widehat{o}_{ij}) \\
        &= \frac{1}{|\mathcal{U}|} \sum_{(i,j) \in \mathcal{U}} \left[ o_{ij} \delta(1, \widehat{o}_{ij}) + (1-o_{ij}) \delta(0, \widehat{o}_{ij}) \right] \\
    \therefore\,\, {\bf E}_o[\widehat{R}_{\text{naive}}(o, \widehat{y})] &= \frac{1}{|\mathcal{U}|} \sum_{(i,j) \in \mathcal{U}} \left[ y_{ij} \pi_{ij} \delta(1, \widehat{o}_{ij}) + (1-y_{ij} \pi_{ij}) \delta(0, \widehat{o}_{ij}) \right] \\
    &= \frac{1}{|\mathcal{U}|} \sum_{(i,j) \in \mathcal{U}} \left[ y_{ij} \pi_{ij} (1 - \widehat{o}_{ij})\delta(1, 0) + (1-y_{ij} \pi_{ij}) \widehat{o}_{ij} \delta(0, 1) \right] \\
    &= \frac{\Delta}{|\mathcal{U}|} \sum_{(i,j) \in \mathcal{U}} \left[ y_{ij} \pi_{ij} (1 - \widehat{o}_{ij}) + (1-y_{ij} \pi_{ij}) \widehat{o}_{ij} \right].
\end{align*}
The true risk is
\begin{align*}
    R(\widehat{y}) &= \frac{1}{|\mathcal{U}|} \sum_{(i,j) \in \mathcal{U}} \left[y_{ij} \delta(1, \widehat{o}_{ij}) + (1 - y_{ij}) \delta(0, \widehat{o}_{ij}) \right] \\
    &= \frac{1}{|\mathcal{U}|} \sum_{(i,j) \in \mathcal{U}} \left[y_{ij} (1-\widehat{o}_{ij}) \delta(1, 0) + (1 - y_{ij}) \widehat{o}_{ij} \delta(0, 1) \right] \\
    &= \frac{\Delta}{|\mathcal{U}|} \sum_{(i,j) \in \mathcal{U}} \left[y_{ij} (1-\widehat{o}_{ij}) + (1 - y_{ij}) \widehat{o}_{ij} \right].
\end{align*}
Thus the bias is
\begin{align*}
    \text{B}(\widehat{R}_{\text{naive}}) &= \left|{\bf E}[\widehat{R}_{\text{naive}}] - R(\widehat{o})\right| \\
    &= \left| \frac{\Delta}{|\mathcal{U}|} \sum_{(i,j) \in \mathcal{U}} \left[ y_{ij} \pi_{ij} (1 - \widehat{o}_{ij}) + (1-y_{ij} \pi_{ij}) \widehat{o}_{ij} - y_{ij} (1-\widehat{o}_{ij}) - (1 - y_{ij}) \widehat{o}_{ij} \right]  \right| \\
    &= \frac{\Delta}{|\mathcal{U}|} \left| \sum_{(i,j) \in \mathcal{U}} y_{ij} (1 - \pi_{ij}) (1 - 2 \widehat{o}_{ij}) \right|.
\end{align*}
The variance is
\begin{align*}
    \text{Var}(\widehat{R}_{\text{naive}}) &= \text{Var}\left( \frac{1}{|\mathcal{U}|} \sum_{(i,j) \in \mathcal{U}} \left[ o_{ij} \delta(1, \widehat{o}_{ij}) + (1-o_{ij}) \delta(0, \widehat{o}_{ij}) \right] \right) \\
    &= \frac{1}{|\mathcal{U}|^2} \sum_{(i,j) \in \mathcal{U}} \text{Var}\left(  o_{ij} \delta(1, \widehat{o}_{ij}) + (1-o_{ij}) \delta(0, \widehat{o}_{ij}) \right) \\
    &= \frac{1}{|\mathcal{U}|^2} \sum_{(i,j) \in \mathcal{U}} y_{ij} \pi_{ij} (1 - y_{ij} \pi_{ij}) \left(\delta(1, \widehat{o}_{ij}) - \delta(0, \widehat{o}_{ij}) \right)^2 \,\,\, \text{(using Lemma \ref{lemma:variance-bernoulli})} \\
    &= \frac{\Delta^2}{|\mathcal{U}|^2} \sum_{(i,j) \in \mathcal{U}} y_{ij} \pi_{ij} (1 - y_{ij} \pi_{ij}).
\end{align*}
\end{proof}
Lemmas \ref{lemma:rw_bias_var}, \ref{lemma:bias-var-pu}, and \ref{lemma:bias-var-ap} can be proved similarly.
\newline

\textbf{Proof of Theorem \ref{thm:comp_var}}

\begin{theorem*}[\textbf{Comparison of Variances}]
For all values of $\widehat{\pi}, \widehat{y}$, we have $\text{Var}(\widehat{R}_{\text{AP}}) < \text{Var}(\widehat{R}_{\text{naive}}), \,\, \text{and} \\
\text{Var}(\widehat{R}_{\text{AP}}) < \text{Var}(\widehat{R}_w) < \text{Var}(\widehat{R}_{\text{PU}})$
\end{theorem*}
\begin{proof}
First we show that $\text{Var}(\widehat{R}_{\text{AP}}) < \text{Var}(\widehat{R}_{\text{naive}})$. We have
\begin{align*}
    \text{Var}(\widehat{R}_{\text{AP}}) &= \frac{\Delta^2}{|\mathcal{U}|^2} \sum_{(i,j) \in \mathcal{U}} y_{ij}\pi_{ij} (1 - y_{ij}\pi_{ij}) \psi^2_{ij}, \\
    \text{where} \,\, & \psi_{ij} = \frac{1-\widehat{y}_{ij}}{1-\widehat{\pi}_{ij}\widehat{y}_{ij}} < 1, \\
    \text{Var}(\widehat{R}_{\text{naive}}) &= \frac{\Delta^2}{|\mathcal{U}|^2} \sum_{(i,j) \in \mathcal{U}} y_{ij}\pi_{ij} (1 - y_{ij}\pi_{ij}).
\end{align*}
Using the fact that $\psi_{ij}^2 < 1 \,\, \forall (i,j) \in \mathcal{U}$, we get $\text{Var}(\widehat{R}_{\text{AP}}) < \text{Var}(\widehat{R}_{\text{naive}})$.

Next, we show that $\text{Var}(\widehat{R}_w) < \text{Var}(\widehat{R}_{\text{PU}})$:
\begin{align*}
    \text{Var}(\widehat{R}_w) &= \frac{\Delta^2}{|\mathcal{U}|^2} \sum_{(i,j) \in \mathcal{U}} y_{ij}\pi_{ij} (1 - y_{ij}\pi_{ij}) \left( \frac{1 - \widehat{o}_{ij}}{\widehat{\pi}^2_{ij}} + \widehat{o}_{ij} \psi^2_{ij} \right) \\
    \text{Var}(\widehat{R}_{\text{PU}}) &= \frac{\Delta^2}{|\mathcal{U}|^2} \sum_{(i,j) \in \mathcal{U}} \frac{y_{ij}\pi_{ij} (1 - y_{ij}\pi_{ij})}{\widehat{\pi}^2_{ij}} \\
        &= \frac{\Delta^2}{|\mathcal{U}|^2} \sum_{(i,j) \in \mathcal{U}} y_{ij}\pi_{ij} (1 - y_{ij}\pi_{ij}) \left( \frac{1 - \widehat{o}_{ij}}{\widehat{\pi}^2_{ij}} + \frac{\widehat{o}_{ij}}{\widehat{\pi}^2_{ij}}  \right) \\
    \therefore\,\, \text{Var}(\widehat{R}_w) - \text{Var}(\widehat{R}_{\text{PU}}) &= \frac{\Delta^2}{|\mathcal{U}|^2} \sum_{(i,j) \in \mathcal{U}} y_{ij}\pi_{ij} (1 - y_{ij}\pi_{ij}) \widehat{o}_{ij} \left( \psi^2_{ij} - \frac{1}{\widehat{\pi}^2_{ij}} \right) \\
    \therefore\,\, \text{Var}(\widehat{R}_w) - \text{Var}(\widehat{R}_{\text{PU}}) &< 0 \,\,\, \left(\text{because $\psi_{ij} < 1$ and $\frac{1}{\widehat{\pi}_{ij}} > 1$}\right) \\
    \therefore\,\, \text{Var}(\widehat{R}_w) &< \text{Var}(\widehat{R}_{\text{PU}}).
\end{align*}
Next, we show that $\text{Var}(\widehat{R}_{\text{AP}}) < \text{Var}(\widehat{R}_w)$:
\begin{align*}
    \text{Var}(\widehat{R}_{\text{AP}}) - \text{Var}(\widehat{R}_w) &= \frac{\Delta^2}{|\mathcal{U}|^2} \sum_{(i,j) \in \mathcal{U}} y_{ij}\pi_{ij} (1 - y_{ij}\pi_{ij}) (1 - \widehat{o}_{ij}) \left( \psi^2_{ij} - \frac{1}{\widehat{\pi}^2_{ij}} \right) \\
    \therefore\,\, \text{Var}(\widehat{R}_{\text{AP}}) - \text{Var}(\widehat{R}_w) &< 0 \,\,\, \left(\text{because $\psi_{ij} < 1$ and $\frac{1}{\widehat{\pi}_{ij}} > 1$}\right) \\
    \therefore\,\, \text{Var}(\widehat{R}_{\text{AP}}) &< \text{Var}(\widehat{R}_w).
\end{align*}
\end{proof}

\textbf{Proof of Theorem \ref{thm:bias-comparisons}}

\begin{theorem*}[\textbf{Comparison of Biases}]
Under the bias approximations, a sufficient condition 
for $\text{B}(\widehat{R}_w) = \text{B}(\widehat{R}_{\text{PU}}) < \text{B}(\widehat{R}_{\text{naive}})$ is
\begin{align*}
    \frac{\pi_{ij}}{2 - \pi_{ij}} < \widehat{\pi}_{ij} < 1, \,\, \forall (i, j) \in \mathcal{U},
\end{align*}
and for $\text{B}(\widehat{R}_{\text{AP}}) < \text{B}(\widehat{R}_{\text{naive}})$ is
\begin{align*}
    & \frac{\pi_{ij}}{2 - \pi_{ij}} < \widehat{\pi}_{ij} < 1 \,\, \text{and} \,\, 0 < \widehat{y}_{ij} < c y_{ij}, \,\, \forall (i,j) \in \mathcal{U} \\
    & \text{where}\,\, c = \frac{2 (1 - \pi_{ij})}{1 - \widehat{\pi}_{ij} - \pi_{ij} y_{ij} + (2 - \pi_{ij}) \widehat{\pi}_{ij} y_{ij}} \geq 1.
\end{align*}
\end{theorem*}
\begin{proof}
We first derive the sufficient condition for $\text{B}(\widehat{R}_w) = \text{B}(\widehat{R}_{\text{PU}}) < \text{B}(\widehat{R}_{\text{naive}})$. We have
\begin{align*}
    \text{B}(\widehat{R}_{\text{naive}}) &\approx \frac{\Delta}{|\mathcal{U}|} \sum_{(i,j) \in \mathcal{U}'} y_{ij} (1 - \pi_{ij}), \\
    \text{B}(\widehat{R}_w) \approx \text{B}(\widehat{R}_{\text{PU}}) &\approx \begin{aligned}[t]
        & \frac{\Delta}{|\mathcal{U}|} \left| \sum_{(i,j) \in \mathcal{U}'} y_{ij}\left(1 - \frac{\pi_{ij}}{\widehat{\pi}_{ij}} \right) \right|.
        \end{aligned}
\end{align*}
If $1 > \widehat{\pi}_{ij} > \pi_{ij} \forall (i,j) \in \mathcal{U}$, we have
\begin{align*}
    & \left(1 - \frac{\pi_{ij}}{\widehat{\pi}_{ij}} \right) > 0 \,\,\, \forall (i,j) \in \mathcal{U} \\
    \therefore\,\, & \text{B}(\widehat{R}_w) \approx \text{B}(\widehat{R}_{\text{PU}}) \approx \begin{aligned}[t]
        & \frac{\Delta}{|\mathcal{U}|} \sum_{(i,j) \in \mathcal{U}'} y_{ij}\left(1 - \frac{\pi_{ij}}{\widehat{\pi}_{ij}} \right).
        \end{aligned} \\
    \therefore\,\, & \begin{aligned}[t]
    \text{B}(\widehat{R}_w) - \text{B}(\widehat{R}_{\text{naive}}) &= \frac{\Delta}{|\mathcal{U}|} \sum_{(i,j) \in \mathcal{U}'} y_{ij} \left( \pi_{ij} - \frac{\pi_{ij}}{\widehat{\pi}_{ij}} \right) \\
        &= \frac{\Delta}{|\mathcal{U}|} \sum_{(i,j) \in \mathcal{U}'} y_{ij} \pi_{ij} \left( 1 - \frac{1}{\widehat{\pi}_{ij}} \right) \\
        &< 0 \,\,\, \text{(because $\widehat{\pi}_{ij} < 1$)}
    \end{aligned} \\
    \therefore\,\, & \text{B}(\widehat{R}_w) < \text{B}(\widehat{R}_{\text{naive}}).
\end{align*}
If $0 < \widehat{\pi}_{ij} \leq \pi_{ij} \forall (i,j) \in \mathcal{U}$, we have
\begin{align*}
    & \left(1 - \frac{\pi_{ij}}{\widehat{\pi}_{ij}} \right) \leq 0 \,\,\, \forall (i,j) \in \mathcal{U} \\
    \therefore\,\, & \text{B}(\widehat{R}_w) \approx \text{B}(\widehat{R}_{\text{PU}}) \approx \begin{aligned}[t]
        & \frac{\Delta}{|\mathcal{U}|} \sum_{(i,j) \in \mathcal{U}'} y_{ij}\left( \frac{\pi_{ij}}{\widehat{\pi}_{ij}} - 1 \right).
        \end{aligned}.
\end{align*}
Then, a sufficient condition for $\text{B}(\widehat{R}_w) = \text{B}(\widehat{R}_{\text{PU}}) < \text{B}(\widehat{R}_{\text{naive}})$ is
\begin{align*}
    y_{ij}\left( \frac{\pi_{ij}}{\widehat{\pi}_{ij}} - 1 \right) &< y_{ij} (1 - \pi_{ij}) \,\,\, \forall (i,j) \in \mathcal{U} \\
    \therefore\,\, y_{ij}\left( \frac{\pi_{ij}}{\widehat{\pi}_{ij}} - 1 \right) &< y_{ij} (1 - \pi_{ij}) \,\,\, \forall (i,j) \in \mathcal{U} \\
    \therefore\,\, \widehat{\pi}_{ij} &> \frac{2}{2 - \pi_{ij}} \,\,\, \forall (i,j) \in \mathcal{U}.
\end{align*}
Next, we derive the sufficient condition for $\widehat{R}_{\text{AP}} < \widehat{R}_{\text{naive}}$. Observe that
\begin{align*}
    & \frac{\pi_{ij}}{2 - \pi_{ij}} < \widehat{\pi}_{ij} < 1 \,\,\, \forall (i, j) \in \mathcal{U} \\
    \therefore\,\, & (1 -\pi_{ij}) y_{ij} - (1-\pi_{ij}y_{ij}) \tau_{ij} \geq 0 \,\,\, \forall (i, j) \in \mathcal{U} \\
    \therefore\,\, & \widehat{R}_{\text{AP}} \approx \frac{\Delta}{|\mathcal{U}|} \sum_{(i,j) \in \mathcal{U}'} \left[ (1-\pi_{ij}) y_{ij} - (1-\pi_{ij}y_{ij}) \tau_{ij} \right], \,\, \text{where} \,\, \tau_{ij} = \left( \frac{\widehat{y}_{ij}(1 - \widehat{\pi}_{ij})}{1 - \widehat{\pi}_{ij} \widehat{y}_{ij}} \right).
\end{align*}
Therefore, when $\frac{\pi_{ij}}{2 - \pi_{ij}} < \widehat{\pi}_{ij} < 1 \,\,\, \forall (i, j) \in \mathcal{U}$, a sufficient condition for $\widehat{R}_{\text{AP}} < \widehat{R}_{\text{naive}}$ is
\begin{align*}
    (1 -\pi_{ij}) y_{ij} - (1-\pi_{ij}y_{ij}) \tau_{ij} < y_{ij} (1 - \pi_{ij})  \,\,\, \forall (i, j) \in \mathcal{U} \\
    \therefore\,\, 0 < \widehat{y}_{ij} < \left( \frac{2 (1 - \pi_{ij})}{1 - \widehat{\pi}_{ij} - \pi_{ij} y_{ij} + (2 - \pi_{ij}) \widehat{\pi}_{ij} y_{ij}} \right) y_{ij} \,\,\, \forall (i, j) \in \mathcal{U}.
\end{align*}
\end{proof}

\section{Generalization Bound} \label{appendix:generalization-bound}

\textbf{Proof of Theorem \ref{thm:generalization-bound}}

\begin{theorem*}[\textbf{Generalization Bound}]
Let $\mathcal{F}$ be a class of functions 
$(\widehat{\pi}, \widehat{y})$.
Let $\delta(o_{ij}, \widehat{y}_{ij}) \leq \eta\,\,\forall (i,j) \in \mathcal{U}$ 
and $\widehat{\pi}_{ij} \geq \epsilon > 0\,\,\forall (i,j) \in \mathcal{U}$.
Then, for $\widehat{R} \in \left\{ \widehat{R}_w,  \widehat{R}_{\text{PU}}, \widehat{R}_{\text{AP}} \right\}$, with probability at least $1 - \delta$, we have
\begin{align}
    R(\widehat{y}) &\leq \widehat{R}(\widehat{y}, \widehat{\pi}) + B(\widehat{R}) + 2 \mathcal{G}(\mathcal{F}, \widehat{R}) + M \label{eq:rademacher-bound} \\
        &\leq \widehat{R}(\widehat{y}, \widehat{\pi}) + B(\widehat{R}_w) + 2 \widehat{\mathcal{G}}(\mathcal{F}, \widehat{R}_w) + 3 M, \label{eq:rademacher-empirical-bound}
\end{align}
where $M = \sqrt{\frac{4 \eta^2}{\epsilon^2 |\mathcal{U}|} \log(\frac{2}{\delta})}$ and $B(\widehat{R})$ is the bias of $\widehat{R}$ derived in Section~\ref{sec:estimating-risk}.
\end{theorem*}
\begin{proof}
We proceed similarly to the standard Rademacher complexity generalization bound
proof \citep{shalev2014understanding}[Ch. 26]. Observe that
\begin{align}
    R(\widehat{y}) &= R(\widehat{y}) - {\bf E}_o[\widehat{R}(o, \widehat{y}, \widehat{\pi})] + {\bf E}_o[\widehat{R}(o, \widehat{y}, \widehat{\pi})] \nonumber \\
        &\leq B(\widehat{R}) + {\bf E}_o[\widehat{R}(o, \widehat{y}, \widehat{\pi})]. \label{eq:rademacher-proof-1}
\end{align}
Let $\Phi(o) = \sup_{(\widehat{\pi}, \widehat{y}) \in \mathcal{F}} \left[ {\bf E}_o[\widehat{R}(o, \widehat{y}, \widehat{\pi})] - \widehat{R}(o, \widehat{y}, \widehat{\pi}) \right]$. Then
\begin{align}\label{eq:rademacher-proof-2}
     {\bf E}_o[\widehat{R}(o, \widehat{y}, \widehat{\pi})]| \leq \widehat{R}(o, \widehat{y}, \widehat{\pi}) + \Phi(o).
\end{align}
Now we upper bound $\Phi(o)$.
Since $\delta(o_{ij}, \widehat{y}_{ij}) \leq \eta\,\,\forall (i,j)$ 
and $\widehat{\pi}_{ij} \geq \epsilon > 0,\,\,\forall (i,j)$ and $\forall \, \widehat{R} \in \left\{ \widehat{R}_w,  \widehat{R}_{\text{PU}}, \widehat{R}_{\text{AP}} \right\}$,
we have
\begin{align*}
    |\Phi(o) - \Phi(\Tilde{o})| \leq \frac{2 \eta}{\epsilon},
\end{align*}
if $o$ and $\Tilde{o}$ differ in only one coordinate, i.e., 
$o_{ij} \neq \Tilde{o}_{ij}$ for some $(i,j) \in \mathcal{U}$ 
and $o_{lm} = \Tilde{o}_{lm} \forall (l,m) \in \mathcal{U}$ 
s.t. $(i, j) \neq (l,m)$. 
Using McDiarmid’s Inequality, with probability at least $1 - \delta$, 
we have
\begin{align} \label{eq:rademacher-proof-3}
    \Phi(o) \leq {\bf E}[\Phi(o)] + C.
\end{align}
Next, we upper bound ${\bf E}[\Phi(o)]$.
Let $\bar{o}$ be a ghost sample independently drawn having the same distribution as $o$. We have
\begin{align}
    {\bf E}[\Phi(o)] &= {\bf E}_o\left[ \sup_{(\widehat{\pi}, \widehat{y}) \in \mathcal{F}} \left[ {\bf E}_o[\widehat{R}(o, \widehat{y}, \widehat{\pi})] - \widehat{R}(o, \widehat{y}, \widehat{\pi}) \right] \right] \nonumber \\
    &= {\bf E}_o\left[ \sup_{(\widehat{\pi}, \widehat{y}) \in \mathcal{F}} {\bf E}_{\bar{o}}\left[\widehat{R}(\bar{o}, \widehat{y}, \widehat{\pi}) - \widehat{R}(o, \widehat{y}, \widehat{\pi}) \,\big| \, o \right] \right] \nonumber \\
    &= {\bf E}_o\left[ \sup_{(\widehat{\pi}, \widehat{y}) \in \mathcal{F}} {\bf E}_{\bar{o}}\left[ \frac{1}{|\mathcal{U}|} \sum_{(i,j) \in \mathcal{U}} r(\bar{o}_{ij}, \widehat{\pi}_{ij}, \widehat{y}_{ij}) - \frac{1}{|\mathcal{U}|} \sum_{(i,j) \in \mathcal{U}} r(o_{ij}, \widehat{\pi}_{ij}, \widehat{y}_{ij}) \,\Big| \, o \right] \right] \nonumber \\
    &\leq {\bf E}_{o, \bar{o}} \left[ \sup_{(\widehat{\pi}, \widehat{y}) \in \mathcal{F}} \left[ \frac{1}{|\mathcal{U}|} \sum_{(i,j) \in \mathcal{U}} r(\bar{o}_{ij}, \widehat{\pi}_{ij}, \widehat{y}_{ij}) - \frac{1}{|\mathcal{U}|} \sum_{(i,j) \in \mathcal{U}} r(o_{ij}, \widehat{\pi}_{ij}, \widehat{y}_{ij}) \right] \right] \,\, \text{(Jensen's Inequality)} \nonumber \\
    &= {\bf E}_{o, \bar{o}, \sigma} \left[ \sup_{(\widehat{\pi}, \widehat{y}) \in \mathcal{F}} \left[ \frac{1}{|\mathcal{U}|} \sum_{(i,j) \in \mathcal{U}} \sigma_{ij} r(\bar{o}_{ij}, \widehat{\pi}_{ij}, \widehat{y}_{ij}) - \frac{1}{|\mathcal{U}|} \sum_{(i,j) \in \mathcal{U}} \sigma_{ij} r(o_{ij}, \widehat{\pi}_{ij}, \widehat{y}_{ij}) \right] \right] \nonumber \\
    &= {\bf E}_{o, \bar{o}, \sigma} \left[ \sup_{(\widehat{\pi}, \widehat{y}) \in \mathcal{F}} \left[ \frac{1}{|\mathcal{U}|} \sum_{(i,j) \in \mathcal{U}} \sigma_{ij} r(\bar{o}_{ij}, \widehat{\pi}_{ij}, \widehat{y}_{ij}) + \frac{1}{|\mathcal{U}|} \sum_{(i,j) \in \mathcal{U}} \sigma_{ij} r(o_{ij}, \widehat{\pi}_{ij}, \widehat{y}_{ij}) \right] \right] \nonumber \\
    &\leq {\bf E}_{o, \bar{o}, \sigma} \left[ \sup_{(\widehat{\pi}, \widehat{y}) \in \mathcal{F}} \left[ \frac{1}{|\mathcal{U}|} \sum_{(i,j) \in \mathcal{U}} \sigma_{ij} r(\bar{o}_{ij}, \widehat{\pi}_{ij}, \widehat{y}_{ij}) \right] + \sup_{(\widehat{\pi}, \widehat{y}) \in \mathcal{F}} \left[ \frac{1}{|\mathcal{U}|} \sum_{(i,j) \in \mathcal{U}} \sigma_{ij} r(o_{ij}, \widehat{\pi}_{ij}, \widehat{y}_{ij}) \right] \right] \nonumber \\
    &= 2 \mathcal{G}(\mathcal{F}, \widehat{R}). \label{eq:rademacher-proof-4}
\end{align}
Combining Eqs. \ref{eq:rademacher-proof-1},
\ref{eq:rademacher-proof-2}, \ref{eq:rademacher-proof-3},
and \ref{eq:rademacher-proof-4}, we get Eq. \ref{eq:rademacher-bound}.
Another application of McDiarmid’s Inequality 
allows us to obtain Eq. \ref{eq:rademacher-empirical-bound}
from Eq. \ref{eq:rademacher-bound}.
\end{proof}

\section{Feedback Loops} \label{appendix:feedback-loops}

\begin{lemma}[\textbf{Binomial Tail Bound}]\label{lemma:appendix-binomial-tail}
If the random variable $X_n \sim \frac{1}{n} \text{Binomial}(n, \theta)$, then for $\epsilon > 0$, we have
\begin{align*}
    {\bf P}(|X_n - \theta| > \epsilon) \leq 2 \exp\left( -2n \epsilon^2 \right).
\end{align*}
\begin{proof}
Observe that $X_n \in [0, 1]$. Applying Hoeffding's inequality gives us the desired result.
\end{proof}
\end{lemma}

\begin{lemma}\label{lemma:appendix-thm-proof-lemma-1}
Let $n \in \mathbb{N}$ and $\kappa$ be a fixed $C - 1$ simplex such that $\kappa_v n \in \mathbb{N} \,\,\, \forall v \in [C]$. The random variable $\Tilde{q}_v \sim \frac{1}{\kappa_v n} \text{Binomial}(\kappa_v n, q_v)$, where $q_v \in (0, 1)$. Assume that $q_v > q_w$ if $v > w$. We denote as $\widehat{e}$ the following $C-1$ simplex:
\begin{align*}
    \widehat{e} = \frac{1}{Z} \left[ \kappa_1 \Tilde{q}_1, \kappa_2 \Tilde{q}_2, \hdots, \kappa_C \Tilde{q}_C \right], \,\,\, \text{where} \,\, Z = \sum_{i \in [C]} \kappa_i \Tilde{q}_i.
\end{align*}
Let $\widehat{e}_{vw} = \frac{\widehat{e}_v}{\widehat{e}_v + \widehat{e}_w} = \frac{\kappa_v \Tilde{q}_v}{\kappa_w \Tilde{q}_w + \kappa_w \Tilde{q}_w}$ and $\kappa_{vw} = \frac{\kappa_v}{\kappa_v + \kappa_w}$. Then for a constant $\rho_{vw}$ such that
\begin{align*}
    0 < \rho_{vw} < \frac{\kappa_v \kappa_w (q_v - q_w)}{q_v \kappa^2_v + (q_v+q_w)\kappa_v \kappa_w + q_w \kappa^2_w},
\end{align*}
we have
\begin{align*}
    & |\Tilde{q}_v - q_v| < \epsilon_{vw}, \,\, |\Tilde{q}_w - q_w| < \epsilon_{vw} \implies \widehat{e}_{vw} - \kappa_{vw} > \rho_{vw}, \\
    & \text{for some constant $\epsilon_{vw}$ s.t.} \,\, 0 < \epsilon_{vw} < \frac{\rho_{vw} q_v \kappa^2_v - \kappa_v \kappa_w (q_w - q_v) + q_w \rho_{vw} \kappa_v (\kappa_v - \kappa_w)}{\rho_{vw} (\kappa^2_v - \kappa^2_w ) - 2 \kappa_v \kappa_w }.
\end{align*}
This is saying that, for $(v, w)$ s.t. $v > w$ the simplex $\widehat{e}$ will be more skewed towards $v$ than the simplex $\kappa$ if the sampled $\Tilde{q}_v$ and $\Tilde{q}_w$ are close to their mean values $q_v$ and $q_w$, respectively.
\end{lemma}
\begin{proof}
Observe that if $|\Tilde{q}_v - q_v| < \epsilon_{vw}$ and $|\Tilde{q}_w - q_w| < \epsilon_{vw}$, then the lowest value that $\widehat{e}_{vw}$ can take is
\begin{align*}
    & \widehat{e}^{(\text{min})}_{vw} = \frac{\kappa_v (q_v - \epsilon_{vw})}{\kappa_v (q_v - \epsilon_{vw}) + \kappa_w (q_w + \epsilon_{vw})}, \,\,\, \text{and} \\
    & \widehat{e}^{(\text{min})}_{vw} - \kappa_{vw} > \rho_{vw} \implies \widehat{e}_{vw} - \kappa_{vw} > \rho_{vw}.
\end{align*}
Therefore, we have
\begin{align*}
    & \widehat{e}^{(\text{min})}_{vw} - \kappa_{vw} > \rho_{vw} \,\,\, \text{and} \,\,\, \epsilon_{vw} < q_w \\
    & \impliedby \underbrace{\frac{\kappa_v (q_v - \epsilon_{vw})}{\kappa_v (q_v - \epsilon_{vw}) + \kappa_w (q_w + \epsilon_{vw})} - \frac{\kappa_v}{\kappa_v + \kappa_w} > \rho_{vw}}_{(1)} \,\,\, \text{and} \,\,\, \rho_{vw} < \frac{\kappa_v \kappa_w (q_v - q_w)}{q_v \kappa^2_v + (q_v+q_w)\kappa_v \kappa_w + q_w \kappa^2_w}.
\end{align*}
The inequality (1) above can further be simplified as
\begin{align*}
    & \frac{\kappa_v (q_v - \epsilon_{vw})}{\kappa_v (q_v - \epsilon_{vw}) + \kappa_w (q_w + \epsilon_{vw})} - \frac{\kappa_v}{\kappa_v + \kappa_w} > \rho_{vw} \\
    \impliedby & \epsilon_{vw} < \frac{\rho_{vw} q_v \kappa^2_v - \kappa_v \kappa_w (q_w - q_v) + q_w \rho_{vw} \kappa_v (\kappa_v - \kappa_w)}{\rho_{vw} (\kappa^2_v - \kappa^2_w ) - 2 \kappa_v \kappa_w }.
\end{align*}
This completes the proof.
\end{proof}

\begin{lemma}\label{lemma:appendix-thm-proof-lemma-2}
Let $\alpha$ be a fixed $C-1$ simplex and $\widehat{e}$ be the following $G-1$ simplex, $\widehat{e} = \frac{1}{Z} [\alpha_1 \Tilde{q}_1, \alpha_2 \Tilde{q}_2, \hdots, \alpha_C \Tilde{q}_C]$, where $Z = \sum_{z \in [C]} \alpha_z \Tilde{q}_z$ and the vector $\kappa \sim \frac{1}{n} \text{Multinomial}(n, \widehat{e})$. Let $\widehat{e}_{vw} = \frac{\widehat{e}_v}{\widehat{e}_v + \widehat{e}_w} = \frac{\Tilde{q}_v}{\Tilde{q}_w + \Tilde{q}_w}$ and $\kappa_{vw} = \frac{\kappa_v}{\kappa_v + \kappa_w}$.

Assume that $|\Tilde{q}_z - q_z| < \epsilon \,\,\, \forall z \in [C]$ where $q_z \in (0, 1)$ are fixed. If $|\kappa_v - \widehat{e}_v < \frac{\eta_{nw}}{C}|$ and $|\kappa_w - \widehat{e}_w < \frac{\eta_{nw}}{C}|$, then for some constant $\rho$, we have
\begin{align*}
    \widehat{e}_{vw} - \kappa_{vw} < \rho, \,\,\,
    \text{when} \,\,\, \eta_{vw} < \rho \left( \frac{q_v + q_w}{\max_{z \in [C]}q_z + \epsilon} \right).
\end{align*}
\begin{proof}
If $|\kappa_v - \widehat{e}_v < \frac{\eta_{nw}}{C}|$ and $|\kappa_w - \widehat{e}_w < \frac{\eta_{nw}}{C}|$, then the smallest value that $\kappa_{vw}$ can achieve is
\begin{align*}
    \kappa^{\text{(min)}}_{vw} = \frac{\widehat{e}_v - \frac{\eta_{nw}}{C}}{\widehat{e}_v + \widehat{e}_w}.
\end{align*}
This means that
\begin{align*}
    & \widehat{e}_{vw} - \kappa_{vw} < \rho \\
    \impliedby & \widehat{e}_{vw} - \kappa^{\text{(min)}}_{vw} < \rho \\
    \iff & \frac{\frac{\eta_{nw}}{C}}{\widehat{e}_v + \widehat{e}_w} < \rho \\
    \iff & \frac{\eta_{nw}}{C} < \rho (\widehat{e}_v + \widehat{e}_w).
\end{align*}
Since $|\Tilde{q}_z - q_z| < \epsilon \,\,\, \forall z \in [C]$, we have
\begin{align*}
    \widehat{e}_v &= \frac{\alpha_v \Tilde{q}_v}{\sum_{z \in [C]} \alpha_z \Tilde{q}_z} \\
        &> \frac{\alpha_v (q_v - \epsilon)}{\alpha_v (q_v - \epsilon) + \sum_{z \in [C], z \neq v} \alpha_z (q_z + \epsilon)} \\
        &> \frac{\alpha_v (q_v - \epsilon)}{C \max_{z \in [C]} \alpha_z (q_z + \epsilon)}, \\
    \text{and similarly} \,\,\, \widehat{e}_w &> \frac{\alpha_w (q_w - \epsilon)}{C \max_{z \in [C]} \alpha_w (q_z + \epsilon)} \\
    \therefore\,\, \widehat{e}_v + \widehat{e}_w &> \frac{ \alpha_v (q_v - \epsilon) + \alpha_w (q_w - \epsilon) }{C \max_{z \in [C]} (q_z + \epsilon)}.
\end{align*}
Therefore, we can set $\eta_{vw}$ such that
\begin{align*}
    \frac{\eta_{nw}}{C} &< \rho \left(\frac{\alpha_v (q_v - \epsilon) + \alpha_w (q_w - \epsilon) }{C \max_{z \in [C]} (q_z + \epsilon)} \right) \\
    \therefore\,\, \eta_{nw} &< \rho \left( \frac{\alpha_v (q_v - \epsilon) + \alpha_w (q_w - \epsilon) }{\max_{z \in [C]} q_z + \epsilon} \right).
\end{align*}
\end{proof}
\end{lemma}

\textbf{Proof of Theorem \ref{thm:feedback-loop-finite-sample}}

\begin{theorem*}
Suppose that $q_v > q_w$ if $v > w$. 
Let $\kappa^{(t)}_{vw} = \frac{\kappa^{(t)}_v}{\kappa^{(t)}_v + \kappa^{(t)}_w}$.
Let $A^{(t)}_{vw}$ represent the event 
that relative fraction of recommendations 
from $c_v$ to that from $c_w$ increases at time $t$,
i.e., $\kappa^{(t+1)}_{vw} > \kappa^{(t)}_{vw}$. 
Let $A^{(t)}$ be the event
that all relative fractions get skewed
towards $c_v$ from $c_w$ if $q_v > q_w$, 
i.e., $A^{(t)} = \bigcap_{(v,w) \in \mathcal{S}} A^{(t)}_{vw}$,
where $\mathcal{S} = \{ (v, w): v \in [C], w \in [C], v > w \}$.
Then, for constants $\epsilon, \eta > 0$ 
that only depend on $\kappa^{(t)}$ and $q$, we have
\begin{align*}
    & {\bf P}(A^{(t)} |\kappa^{(t)}) \begin{aligned}[t]
        & \geq 1 - 2C \exp\left( -2 n \left[ \epsilon^2 + \frac{\eta^2}{C^2} \right]  \right) \\
        & \geq 1 - 2C \exp\left( - \mathcal{O}\left(\frac{n}{C^2}\right) \right).
    \end{aligned}
\end{align*}
\end{theorem*}
\begin{proof}
We know that the estimated probabilities $\widehat{q}^{(t)}_v$ have distribution $\widehat{q}^{(t)}_v | \kappa^{(t)} \sim \frac{1}{n} \text{Binomial}(n \kappa^{(t)}_v, q_v)$. The simplex with normalized probabilities is $\widehat{e}^{(t+1)} = \frac{1}{Z} [\widehat{q}^{(t)}_1, \widehat{q}^{(t)}_2, \hdots, \widehat{q}^{(t)}_C]$, where $Z = \sum_{z \in [C]} \widehat{q}^{(t)}_z$.

Let $\Tilde{q}^{(t)}_v = \frac{\widehat{q}^{(t)}_v}{\kappa^{(t)}_v}$.
Observe that $\Tilde{q}^{(t)}_v | \kappa^{(t)} \sim \frac{1}{n \kappa^{(t)}_v} \text{Binomial}(n \kappa^{(t)}_v, q_v)$. We denote by $\widehat{e}^{(t+1)}_{vw}$,
\begin{align*}
    \widehat{e}^{(t+1)}_{vw} = \frac{\widehat{e}^{(t+1)}_v}{\widehat{e}^{(t+1)}_v + \widehat{e}^{(t+1)}_w} = \frac{\kappa^{(t)}_v \Tilde{q}^{(t)}_v}{\kappa^{(t)}_v \Tilde{q}^{(t)}_v + \kappa^{(t)}_w \Tilde{q}^{(t)}_w}.
\end{align*}
There are two main parts to the proof. 
First, we show that, with high probability, 
$\widehat{e}^{(t+1)}_{vw} - \kappa^{(t)}_{vw} > \rho \,\, \forall (v,w) \in \mathcal{S}$ for some constant $\rho$. 
Then, we show that, with high probability, 
$\widehat{e}^{(t+1)}_{vw} - \kappa^{(t+1)}_{vw} < \rho \,\, \forall (v,w) \in \mathcal{S}$. We combine these two results to show that, with high probability, $\kappa^{(t+1)}_{vw} > \kappa^{(t)}_{vw} \,\, \forall (v,w) \in \mathcal{S}$.

Using Lemma \ref{lemma:appendix-thm-proof-lemma-1}, for some $(v, w) \in \mathcal{S}$, we know that for some constant $\rho_{vw}$ such that
\begin{align*}
    0 < \rho_{vw} < \frac{\kappa^{(t)}_v \kappa^{(t)}_w (q_v - q_w)}{q_v (\kappa^{(t)}_v)^2 + (q_v+q_w)\kappa^{(t)}_v \kappa^{(t)}_w + q_w (\kappa^{(t)}_w)^2},
\end{align*}
we have
\begin{align*}
    & |\Tilde{q}^{(t)}_v - q_v| \leq \epsilon_{vw} \,\, \textit{and} \,\, |\Tilde{q}^{(t)}_w - q_w| \leq \epsilon_{vw} \implies \widehat{e}^{(t+1)}_{vw} - \kappa^{(t)}_{vw} \geq \rho_{vw}, \\
    & \text{for a constant $\epsilon_{vw}$ s.t.} \,\, 0 < \epsilon_{vw} < \frac{\rho_{vw} q_v (\kappa^{(t)}_v)^2 - \kappa^{(t)}_v \kappa^{(t)}_w (q_w - q_v) + q_w \rho_{vw} \kappa^{(t)}_v (\kappa^{(t)}_v - \kappa^{(t)}_w)}{\rho_{vw} ((\kappa^{(t)}_v)^2 - (\kappa^{(t)}_w)^2 ) - 2 \kappa^{(t)}_v \kappa^{(t)}_w } \\
    \implies & {\bf P}\left(\widehat{e}^{(t+1)}_{vw} - \kappa^{(t)}_{vw} \geq \rho_{vw}\right) \geq {\bf P}\left(|\Tilde{q}^{(t)}_v - q_v| \leq \epsilon_{vw}, |\Tilde{q}^{(t)}_w - q_w| \leq \epsilon_{vw}\right).
\end{align*}
Intuitively, this is saying 
that $\widehat{e}^{(t+1)}_{vw} - \kappa^{(t)}_{vw} > \rho_{vw}$ if $\Tilde{q}^{(t)}_v$ and $\Tilde{q}^{(t)}_w$ are close to $q_v$ and $q_w$, respectively.
Let $\rho = \min_{(v,w) \in \mathcal{S}} \rho_{vw}$ and $\epsilon = \min_{(v,w) \in \mathcal{S}} \epsilon_{vw}$. Then we have
\begin{align}
    {\bf P} \left( \bigcap_{(v,w) \in \mathcal{S}} \widehat{e}^{(t+1)}_{vw} - \kappa^{(t)}_{vw} \geq \rho \right) &\geq {\bf P} \left( \bigcap_{z \in [C]} |\Tilde{q}^{(t)}_z - q_z| \leq \epsilon \right) \\
    &= 1 - {\bf P} \left( \bigcup_{z \in [C]} |\Tilde{q}^{(t)}_z - q_z| \geq \epsilon \right) \\
    &\geq 1 - \sum_{z=1}^{C} {\bf P} \left( |\Tilde{q}^{(t)}_z - q_z| \geq \epsilon \right) \,\,\, \text{(Union Bound)} \\
    &\geq 1 - \sum_{z=1}^{C} 2 \exp\left( -2n \epsilon^2 \right)
    \,\,\, \text{(using Lemma \ref{lemma:appendix-binomial-tail})} \nonumber \\
    &= 1 - 2C \exp\left( -2n \epsilon^2 \right). \label{eq:feedback-finite-sample-proof-1}
\end{align}
Now, we show that $\widehat{e}^{(t+1)}_{vw}$ 
is close to $\kappa^{(t+1)}_{vw}$.
We know that $\kappa^{(t+1)} \sim \frac{1}{n} \text{Multinomial}(n, \widehat{e}^{(t+1)})$.
Let the event $Q^{(t)} = \bigcap_{z \in [C]} |\Tilde{q}^{(t)}_z - q_z| \leq \epsilon$.
Using Lemma \ref{lemma:appendix-thm-proof-lemma-2}, we know that, under $Q^{(t)}$, for some constant $\eta_{vw}$, we have 
\begin{align*}
    & \left|\widehat{e}^{(t+1)}_v - \kappa^{(t+1)}_v \right| < \frac{\eta_{vw}}{C} \,\,\, \text{and} \,\,\, \left|\widehat{e}^{(t+1)}_w - \kappa^{(t+1)}_w \right| < \frac{\eta_{vw}}{C} \implies \widehat{e}^{(t+1)}_{vw} - \kappa^{(t+1)}_{vw} < \rho, \\
    & \text{where} \,\, 0 < \eta_{vw} < \frac{ \kappa^{(t+1)}_v (q_v - \epsilon) + \kappa^{(t+1)}_w (q_w - \epsilon)}{\max_{z \in [C]} \kappa^{(t+1)}_z (q_z + \epsilon)} \\
    \implies & {\bf P} \left( \widehat{e}^{(t+1)}_{vw} - \kappa^{(t+1)}_{vw} < \rho \, \big| \, Q^{(t)} \right) \geq {\bf P}\left(\left|\widehat{e}^{(t+1)}_v - \kappa^{(t+1)}_v \right| < \frac{\eta_{vw}}{C}, \,\,\, \left|\widehat{e}^{(t+1)}_w - \kappa^{(t+1)}_w \right| < \frac{\eta_{vw}}{C} \right).
\end{align*}
Intuitively, this is saying that $\widehat{e}^{(t+1)}_{vw} - \kappa^{(t+1)}_{vw} < \rho$ if $\kappa^{(t+1)}_v$ and $\kappa^{(t+1)}_w$ are close to $\widehat{e}^{(t+1)}_v$ and $\widehat{e}^{(t+1)}_w$, respectively. Thus, for $\eta = \min_{(v,w) \in \mathcal{S}} \eta_{vw}$, we have
\begin{align}
    {\bf P} \left( \bigcap_{(v,w) \in \mathcal{S}} \widehat{e}^{(t+1)}_{vw} - \kappa^{(t+1)}_{vw} \leq \rho \, \Big| \, Q^{(t)} \right) &\geq {\bf P}\left( \bigcap_{z \in [C]} |\widehat{e}^{(t+1)}_z - \kappa^{(t+1)}_z| \leq \frac{\eta}{C} \right) \nonumber \\
    &= 1 - {\bf P}\left( \bigcup_{z \in [C]} |\widehat{e}^{(t+1)}_z - \kappa^{(t+1)}_z| \geq \frac{\eta}{C} \right) \nonumber \\  
    &\geq 1 - \sum_{z=1}^{C} {\bf P}\left(|\widehat{e}^{(t+1)}_z - \kappa^{(t+1)}_z| > \frac{\eta}{C} \right) \,\,\, \text{(Union Bound)} \nonumber \\
    &\geq 1 - 2C \exp\left( -\frac{2n \eta^2}{C^2} \right) \,\,\, \text{(using Lemma \ref{lemma:appendix-binomial-tail})}. \label{eq:feedback-finite-sample-proof-2}
\end{align}
Combining Eq. \ref{eq:feedback-finite-sample-proof-1} and \ref{eq:feedback-finite-sample-proof-2}, 
we get the desired result as follows:
\begin{align*}
    {\bf P} \left( \bigcap_{(v,w) \in \mathcal{S}} A^{(t)}_{vw} \right) &= {\bf P} \left( \bigcap_{(v,w) \in \mathcal{S}} \kappa^{(t+1)}_{vw} > \kappa^{(t)}_{vw} \right) \\
    &\geq {\bf P} \left( \bigcap_{(v,w) \in \mathcal{S}} \left( \widehat{e}^{(t+1)}_{vw} - \kappa^{(t)}_{vw} \geq \rho, \,\,\, \widehat{e}^{(t+1)}_{vw} - \kappa^{(t+1)}_{vw} \leq \rho  \right) \right) \\
    &\geq {\bf P}\left( \bigcap_{z \in [C]} \left( |\widehat{e}^{(t+1)}_z - \kappa^{(t+1)}_z| \leq \frac{\eta}{C}, \,\,\, |\Tilde{q}^{(t)}_z - q_z| \leq \epsilon \right) \right) \\
    &= {\bf P}\left( \bigcap_{z \in [C]} |\widehat{e}^{(t+1)}_z - \kappa^{(t+1)}_z| \leq \frac{\eta}{C} \, \Big| \, Q^{(t)} \right) {\bf P}\left( \bigcap_{z \in [C]}  |\Tilde{q}^{(t)}_z - q_z| \leq \epsilon \right) \\
    &\geq \left( 1 - 2C \exp\left( -\frac{2n \eta^2}{C^2} \right) \right) \left( 1 - 2C \exp\left( -2n \epsilon^2 \right) \right) \\
    &\geq 1 - 2C \left[ \exp\left( -2n \epsilon^2 \right) + \exp\left( -\frac{2n \eta^2}{C^2} \right) \right] \\
    &\geq 1 - 2C \exp\left( - \mathcal{O}\left(\frac{n}{C^2}\right) \right).
\end{align*}
\end{proof}

\textbf{Proof of Theorem \ref{thm:worsen_rate}}

\begin{lemma}[Convergence in Probability]\label{lemma:appendix-convergence-in-prob}
Let $X_n, Y_n$, and $Z$ be random variables such that $X_n \overset{p}{\to} Y_n$ and $Y_n \overset{p}{\to} Z$, then $X_n \overset{p}{\to} Z$.
\begin{proof}
For any $\epsilon > 0$, we have
\begin{align*}
    {\bf P}(|X_n - Z| \geq \epsilon) &= {\bf P}(|X_n - Y_n + Y_n - Z| \geq \epsilon) \\
    &\leq {\bf P}(|X_n - Y_n| + |Y_n-Z| \geq \epsilon) \\
    &\leq {\bf P}\left(|X_n - Y_n| \geq \frac{\epsilon}{2}\right) + {\bf P}\left(|X_n - Y_n| \geq \frac{\epsilon}{2}\right) \\
    &= 0.
\end{align*}
Therefore, $X_n \overset{p}{\to} Z$.
\end{proof}
\end{lemma}

\begin{theorem*}
Suppose that $q_v > q_w$. As $n \rightarrow \infty$, $\kappa^{(t)}_{vw} \overset{p}{\rightarrow} 1 - \frac{1}{1 + c^t}$, where $c = \frac{q_v}{q_w}$.
\end{theorem*}
\begin{proof}
At time step $t$, the fraction of recommendations 
from each group is $\kappa_t$. 
From group $g_v$, the user cites papers 
according to probability $q_v$. 
Therefore, $\widehat{q}^{(t)}_v \overset{p}{\rightarrow} \kappa^{(t)}_v q_v$. 
And the normalized estimate is
$\widehat{e}^{(t+1)} = \frac{1}{S} [\kappa^{(t)}_1 q_1, \hdots, \kappa^{(t)}_C q_C]$, where $S = \sum_{z \in [C]} \kappa^{(t)}_z q_z$. 
Since $\kappa^{(t+1)} \sim \frac{1}{n} \text{Multinomial}(n, \widehat{e}^{(t+1)})$,
we have
\begin{align}
    \kappa^{(t+1)} &\overset{p}{\rightarrow} \widehat{e}^{(t+1)} \nonumber \\
    \frac{\kappa^{(t+1)}_v}{\kappa^{(t+1)}_w} &\overset{p}{\rightarrow} \frac{q_v \kappa^{(t)}_v}{q_w \kappa^{(t)}_w} \nonumber \\
    &= c \frac{ \kappa^{(t)}_v}{ \kappa^{(t)}_w}. \label{eq:appendix-feedback-asymptotic-proof}
\end{align}
We know that $\frac{ \kappa^{(1)}_v}{ \kappa^{(1)}_w} \overset{p}{\rightarrow} c$. 
Combining this with Eq. \ref{eq:appendix-feedback-asymptotic-proof} and using Lemma \ref{lemma:appendix-convergence-in-prob} recursively, we get
\begin{align*}
    \frac{\kappa^{(t)}_v}{\kappa^{(t)}_w} &\overset{p}{\rightarrow} c^t \\
    \therefore\,\, 1 - \frac{1}{1 + \frac{\kappa^{(t)}_v}{\kappa^{(t)}_w}} &\overset{p}{\rightarrow} 1 - \frac{1}{1 + c^t} \,\,\, \text{(Continuous mapping theorem)} \\
    \therefore\,\, \frac{\kappa^{(t)}_v}{\kappa^{(t)}_v + \kappa^{(t)}_w} &\overset{p}{\rightarrow} 1 - \frac{1}{1 + c^t} \\
    \therefore\,\, \kappa^{(t)}_{vw} &\overset{p}{\rightarrow} 1 - \frac{1}{1 + c^t}.
\end{align*}
\end{proof}

\section{Experiments} \label{appendix:experiments}

\begin{table}
\caption{The distribution of the FOS in the two real-world datasets.
}
\label{table:fos-distributions-in-dataset}
\begin{center}
\begin{sc}
\begin{tabular}{lcc}
\toprule
FOS & Dataset 1 & Dataset 2 \\
\midrule
Art & $0.03\%$ & $0.08\%$ \\
Biology & $26.48\%$ & $23.43\%$ \\
Business & $0.38\%$ & $0.10\%$ \\
Chemistry & $10.11\%$ & $15.67\%$ \\
Computer Science & $9.40\%$ & $3.42\%$ \\
Economics & $2.51\%$ & $0.03\%$ \\
Engineering & $6.24\%$ & $17.98\%$ \\
Environmental Science & $0.13\%$ & $0.03\%$ \\
Geography & $0.48\%$ & $0.40\%$ \\
Geology & $1.45\%$ & $0.46\%$ \\
History & $0.04\%$ & $0.03\%$ \\
Materials Science & $3.06\%$ & $19.09\%$ \\
Mathematics & $7.17\%$ & $1.03\%$ \\
Medicine & $21.28\%$ & $13.90\%$ \\
Philosophy & $0.03\%$ & $0.01\%$ \\
Physics & $2.99\%$ & $3.14\%$ \\
Political Science & $0.18\%$ & $0.01\%$ \\
Psychology & $7.49\%$ & $1.14\%$ \\
Sociology & $0.55\%$ & $0.05\%$ \\
\bottomrule
\end{tabular}
\end{sc}
\end{center}
\end{table}

Table \ref{table:fos-distributions-in-dataset} provides the distribution of the various FOS in both the datasets used for the real-world dataset experiments (Section~\ref{sec:real-world-datasets}). We can see that the FOS distributions are different. For example, Dataset 2 has substantially more \textit{Materials Science} and \textit{Engineering} papers.