Near-optimal inference in adaptive linear regression
Near-optimal inference in adaptive linear regression

A]Koulik<EMAIL_ADDRESS>B]Yash<EMAIL_ADDRESS>
C]Tor<EMAIL_ADDRESS>D]Lester<EMAIL_ADDRESS>E]Martin J<EMAIL_ADDRESS>
[A]Department of Statistics,
Rutgers University,

Voleon Group, e2

[C]DeepMind, e3
[D]Microsoft Research, e4

[E]Lab for Information and Decision Systems,
Statistics and Data Science Center,
Massachusetts Institute of Technology

Departments of Statistics and EECS,
UC Berkeley

When data is collected in an adaptive manner, even simple methods like
ordinary least squares can exhibit non-normal asymptotic behavior.  As
an undesirable consequence, hypothesis tests and confidence intervals
based on asymptotic normality can lead to erroneous results.  We
propose a family of online debiasing estimators to correct these
distributional anomalies in least squares estimation.  Our proposed
methods take advantage of the covariance structure present in the
dataset and provide sharper estimates in directions for which more
information has accrued.  We establish an asymptotic normality
property for our proposed online debiasing estimators under mild
conditions on the data collection process and provide asymptotically
exact confidence intervals. We additionally prove a minimax lower
bound for the adaptive linear regression problem, thereby providing a
baseline by which to compare estimators. There are various conditions
under which our proposed estimators achieve the minimax lower bound.
We demonstrate the usefulness of our theory via applications to
multi-armed bandit, autoregressive time series estimation, and active
learning with exploration.

[Primary ]62E20
[; secondary ]6208

Adaptive linear regression
Online debiasing
Minimax lower bound
Asymptotic normality
Multi-armed bandits
Autoregressive time series

§ INTRODUCTION

Consider a prediction problem in which we observe $\numobs$ datapoints
of the form $(\x_i, \y_i) \in \real^\Dim \times \real$ with covariate
vector $\x_i$ and response $\y_i$ linked via the linear model

\begin{align}
\label{eqn:linear-reg-model}
\y_i = \myinprod{\x_i}{\thetastar} + \error_i \quad \mbox{for $i = 1,
\ldots, \numobs$.}
\end{align}

Here the vector $\thetastar \in \real^\Dim$ is an unknown parameter of
interest, and $\error_i$ is additive noise.  When the datapoints are
generated via some i.i.d. sampling process, this model, and in
particular the behavior of the ordinary least squares (OLS) estimate
$\thetaLS$, is very well-understood.  The focus of this paper is the
more challenging setting in which the covariate vectors $\Nseq{\x}{i}$
have been adaptively collected, meaning that the choice of
$\x_i$ can depend on the entire set of previous observations $\{\x_j,
\y_j\}_{j =1}^{i-1}$.

More precisely, given a filtration $\Nseq{\filtration}{i}$, assume
that $\x_i$ is $\filtration_{i - 1}$-measurable and
that the additive error $\Nseq{\error}{i}$ is a martingale difference
sequence with respect to $\Nseq{\filtration}{i}$, with

\begin{align}
\Exs[\error_i \mid \filtration_{i - 1} ] = 0, \quad \mbox{and} \quad
\Exs[\error_i ^2 \lvert \filtration_{i-1}] = \sigma^2,
\end{align}

for some non-random scalar $\sigma^2 > 0$.  We refer to the
combination of the linear observation
model (<ref>) with such (potentially) adaptive
collection procedures as the adaptive linear regression model.
Instances of adaptive linear regression arise in a variety of
applications, including multi-armed
bandits [20], active
learning [9], times series
modeling [4], stochastic
control [1], and adaptive stochastic
approximation schemes [18, 6].

Let us discuss some known results for the OLS estimate $\thetaLS$.
The estimate can be expanded in the form

\begin{align}
\label{EqnLSDecompose}
\thetaLS = \SigmaMat_\numobs^{-1}\X_n^\top \y_n = \thetastar
+ \SigmaMat_\numobs^{-1} \sum_{i=11}^\numobs \x_i \error_i, \qquad
\text{where} \quad \SigmaMat_\numobs \mydefn \sum_{i=1}^\numobs \x_i \x_i^\top.
\end{align}

This decomposition reveals that the statistical properties of the OLS
estimate depend on the martingale transform $\Nsum{i} \x_i \error_i$,
along with the random matrix $\SigmaMat_\numobs$.  There is a lengthy
literature on conditions under which the OLS estimate is
consistent [18, 1, 4, 10, 16, 17].
Notably, <cit.> show that the OLS estimate is
strongly consistent, meaning that $\thetaLS \almostSurely \thetastar$,

\begin{align}
\label{eqn:strong-consistency-cond}\lambda_{\min}(\SigmaMat_\numobs) \almostSurely \infty \quad \text{and}
\quad
\frac{\log  \lambda_{\max}(\SigmaMat_\numobs)}{\lambda_{\min}(\SigmaMat_\numobs)}
\almostSurely 0.
\end{align}

Arguably, these conditions for consistency are quite mild.  In
contrast, <cit.> also show that asymptotic
normality of the least squares estimator in the adaptive linear
regression model holds under a stability condition that is
substantially more restrictive—namely, the existence of a sequence
$\{\SigmaLim_\numobs\}_{\numobs \geq 1}$ of non-random strictly
positive definite matrices such that

\begin{align}
\label{eqn:Cond-Cov-stability} \SigmaLim_\numobs^{-1} \SigmaMat_\numobs \stackrel{\smallop{p}}{\longrightarrow} \MyIdDim.
\end{align}

Moreover, <cit.> demonstrate through the
example of a unit root autoregressive model that the OLS estimator
fails to be asymptotically normal in absence of the stability
property (<ref>).  In such cases, confidence
intervals and other forms of inference performed using Gaussian limit
theory are no longer valid.

§.§ Contributions

In this paper, we propose and analyze a new family of estimators for
the parameter vector $\thetastar$ (or linear functionals thereof)
based on online debiasing techniques.  We show that, under mild
conditions, our proposed estimators are both asymptotically unbiased
and asymptotically normal. The underlying assumptions are less
stringent than the stability condition (<ref>)
and are satisfied by a large class of models for data generation and
protocols for choosing covariate vectors.  We provide a detailed
discussion of three such example classes in <ref>.
By deriving minimax lower bounds on the performance of any estimator,
we show that our estimators are minimax optimal.  We also show that
the asymptotic performance of these estimators are near-optimal in an
instance dependent sense, in that they match the performance of the
best problem-specific behavior up to a logarithmic factor.

§.§ Related work

The broader literature on bandit algorithms and experimentation
focuses mostly on a single statistical objective, with standard
examples being minimizing regret or selecting an optimal arm with high
probability. In the papers [29, 26],
the authors empirically observed that bandit algorithms induce bias,
which can be problematic for ex-post inference.  Later
works [21, 22, 23]
characterize the sign and bound the magnitude of this bias. In the
paper [11], the authors develop estimators that
use propensity scores for the multi-armed bandit setting, a special
case of the stochastic regression model (<ref>)
in which the covariate vectors $\x_i$ are restricted to standard basis
vectors. However, it is not clear how to extend this approach to
general designs. Also in the bandit setting, Zhang et
al. [30] develop a least squares estimator that
exploits an assumed batch structure, meaning that only a fixed, finite
number of adaptive decisions are made. This approach, however, does
not apply to more general schemes that make adaptive decisions at each
round. Recently, Zhang et al. [31] proposed a
weighted  for contextual bandit problems where the
bandit algorithm is known. It is also not clear how to generalize this
approach to a more general data collection scheme or to the case when
the data collection algorithm is not completely known.

There is also a parallel line of work that exploits concentration of
measure results (e.g., see the papers [3, 27]) to
develop confidence regions that are valid uniformly in time.  This
approach has its roots in the bandits
literature [2, 13] and has been
refined in more recent
work [12, 14].  An advantage of
this approach is that it yields bounds that are uniform in time. On
the flip side, it requires very strong exponential tail conditions on
the error sequence in contrast to the relatively mild moment conditions
that we impose.  Overall, we view this line of work as being
complementary to our goal of developing corrected estimators that obey
asymptotic normality.

This paper builds upon and extends past work, due to a subset of the
current authors [6, 5],
using online debiasing techniques. In <ref>,
we prove a lower bound that shows how the matrix sequence
$\SigmaMat_\numobs^{-1}$ controls the fundamental difficulty of the
problem, and this lower bound also motivates the particular form of
debiasing proposed in this paper.  The construction used in past
work [6, 5] is based on a
non-adaptive upper bound of the form $\lambda_* \MyIdDim$, where the
scalar $\lambda_*$ is chosen to be much larger than
$\lambda_{\max}(\SigmaMat_\numobs^{-1})$ with high probability.  By
sharp contrast, our analysis instead makes use of an adaptive upper
bound that simultaneously respects the structure of
$\SigmaMat_\numobs^{-1}$ and leads to a stable martingale transform;
this particular construction and our analysis thereof allows us to
obtain sharper guarantees than past
work [6, 5].

§.§ Notation
Let us summarize some notation used throughout
the remainder of the paper.  For a positive integer $n$, we make use
of the convenient shorthand $[n] \mydefn \{1, 2, \ldots, n\}$.  We use
$\myvector{e}_j$ to denote the $j$th standard basis vector in
$\real^{\Dim}$.  For a matrix $\mymatrix{M}$, we use the notation
$\opnorm{\mymatrix{M}}$ and $\frobnorm{\mymatrix{M}}$ to denote the
operator norm (maximum singular value) and the Frobenius norm of the
matrix $\mymatrix{M}$, respectively; similarly, we use the notation
$\maxnorm{\mymatrix{M}}$ to denote the maximum entry in absolute
value.  For a square matrix $\mymatrix{S}$, the quantities
$\lambda_{\max}(\mymatrix{S})$ and $\lambda_{\min}(\mymatrix{S})$
respectively denote the maximum and minimum eigenvalue of the matrix
$\mymatrix{S}$.  The quantity $\tr(\mymatrix{S})$ denotes the sum of
diagonal entries of the square matrix $\mymatrix{S}$. For a pair of
squares matrices $(\mymatrix{A}, \mymatrix{B})$ of compatible
dimensions, we use the notation $\mymatrix{A} \succcurlyeq
\mymatrix{B}$ to indicate that the difference matrix $\mymatrix{A} -
\mymatrix{B}$ is positive semidefinite; we use the notation
$\mymatrix{A} \succ \mymatrix{B}$ when the difference matrix
$\mymatrix{A} - \mymatrix{B}$ is positive definite. The relations
$\mymatrix{A} \preccurlyeq \mymatrix{B}$ and $\mymatrix{A} \prec
\mymatrix{B}$ are defined analogously. For a symmetric positive
semidefinite matrix $\mymatrix{S}$, we use $\matsqrt{\mymatrix{S}}$ to
denote its symmetric matrix square root.

For a sequence of random variables $\{\RV_\numobs\}_{\numobs \geq 1}$
and a random variable $\RV$, we write $\RV_\numobs
\stackrel{\smallop{p}}{\longrightarrow} \RV$ to mean that the sequence
of random variables $\{\RV_\numobs\}_{\numobs \geq 1}$ converges to
$\RV$ in probability; the notation $\RV_\numobs \indistrb \RV$
indicates convergence in distribution.  For a sequence of real-valued
random variables $\{\RV_\numobs\}_{\numobs \geq 1}$ and a sequence of
non-zero real numbers $\{\realseq_\numobs\}_{\numobs \geq 1}$, we
write $\RV_\numobs = \smalloP(\realseq_\numobs)$ to mean that the
ratio $\tfrac{\RV_\numobs}{\realseq_\numobs}
\stackrel{\smallop{p}}{\longrightarrow} 0$.  We write $\RV_\numobs =
\BigoP(\realseq_\numobs)$ to mean that the ratio
$\RV_\numobs/\realseq_\numobs$ is stochastically bounded.  More
precisely, for every scalar $\epsilon >0$, there exits a positive real
number $C_{\epsilon}$ such that

§ FROM ORDINARY LEAST SQUARES TO ONLINE DEBIASING

In this section, we begin by motivating the work by discussing how
classical theory about ordinary least squares estimate can break down
when data is collected in an adaptive manner.  We then introduce an
online debiasing approach to computing alternative estimates.

§.§ Breakdown of the ordinary least squares estimator

Distribution of $\frac{\thetahat_2 -
\theta_2^*}{\sqrt{(\SigmaMat_n^{-1})_{22}} }$

Distribution of $\left(\SigmaMat_n/n\right)_{22}$

Quantitative behavior of the first coordinate $\thetahat_1$
of the ordinary least squares estimator $\thetaLS := (\thetahat_1,
\thetahat_2)$ on a dataset drawn from the $\greedy$-greedy two-armed
bandit model of <ref>. The results are
obtained with a dataset of size $n = 1000$ and $5000$ independent
replications. (a) The distribution of the rescaled difference $\frac{\thetahat_2 -
\theta_2^*}{\sqrt{(\SigmaMat_n^{-1})_{22}} }$ is far from standard
Gaussian.  (b) The bimodal distribution of $\left(\SigmaMat_n/n
\right)_{22}$ suggests that the scaled covariance matrix
$\SigmaMat_n/n$ does not converge to a deterministic matrix.

Let us begin by considering the behavior of the OLS estimator
$\thetaLS$ from equation (<ref>).  When the covariates
$\{\x_i\}_{i \geq 1}$ are either fixed or independently sampled from a
fixed distribution, it has several optimality properties.
Accordingly, it is natural to ask what the performance of the estimator is when the covariates $\{\x_i\}$ are drawn in an adaptive

In order to fix ideas, let us consider a two-armed bandit
problem [20], a special case of the linear
regression model (<ref>) with each $\x_i$ chosen
to be either $(1,0)^\top$ or $(0, 1)^\top$ based on the prior data
$\{\x_j, \y_j \mid j \leq i - 1 \}$.  In order to generate the
covariates $\{\x_i\}_{i = 1}^{n}$, suppose that we apply the
$\greedy$-greedy selection algorithm, a popular choice for tackling
bandit problems [20].

A simple simulation reveals some interesting phenomena.  We generated
linear regression data using $\greedy$-greedy selection algorithm with
the choices $\greedy = 0.1$, $\thetastar = (0.3, 0.3)^\top$, and noise
variables $\error_i \stackrel{i.i.d.}{\sim} \Ncal(0, 1)$.  Let
$\thetahat_2$ denote the first coordinate of the estimator fit to
the bandit data.  Figure <ref>(a) demonstrates that the
distribution of $\thetahat_2$, even after proper re-centering and
scaling, does not converge to a standard normal distribution.  As an
undesirable consequence, the confidence intervals for $\thetastarScalar_2$,
usually constructed using the quantiles of a standard normal random
variable, are not valid.

Let us try to understand why the estimate
fails to be asymptotically normal.
Figure <ref>(b) plots a histogram
of the $(2,2)$ entry of the scaled sample covariance
matrix $\SigmaMat_n/n$.
The bimodal behavior suggests that $\SigmaMat_n/n$ fails to converge to a
non-random matrix $\SigmaLim$ and indicates that the stability
condition (<ref>)
is not satisfied. Indeed, in a recent paper [30],
the authors show that when $\theta_1^* = \theta_2^*$ as in our example, the estimator,
after proper centering and scaling, converges to
a distribution which is not a standard Gaussian distribution.

It turns out that this distributional anomaly of the estimator is
neither specific to the two-armed bandit
problem [20] nor to the $\greedy$-greedy
algorithm used to simulate the data for
Figure <ref>. The same phenomenon was documented in the
time-series and forecasting literature half a century ago, dating back
to the works of 28, 7
and 18. More recent
work [6, 30] has highlighted
that a similar phenomenon commonly occurs in multi-armed bandit
problems when using popular selection algorithms, including Thompson
sampling and the upper confidence bound (UCB)
algorithm [20].

In Section <ref>, we rectify the distributional anomaly
of the estimator by proposing an estimator based on the online
debiasing principles of 6 and show that
our online debiasing estimator exhibits asymptotic normality even in
the absence of the stability condition (<ref>).
In Section <ref>, we demonstrate the usefulness of
our theory via applications to the multi-armed bandit problems,
autoregressive time series, and active learning problems with

§.§ Online debiasing estimator

In this section, we propose and analyze an estimator based on an
online debiasing technique motivated by the
work of 6. At a high-level, the estimator
involves a specific perturbation of the ordinary least squares
estimator $\thetaLS$.  This perturbation is constructed via a linear
combination of the prediction errors $\NseqNosub{\y_i -
\myinprod{\x_i}{\thetaLS}}{i}$ along with a carefully chosen
sequence of weight vectors $\Nseq{\w}{i}$.  The key property ensured by the
construction is that the weight vector $\w_i$ is $\filtration_{i - 1}$
measurable for each $i \in [\numobs]$.

Concretely, for weight vectors $\Nseq{\w}{i}$, we
compute the

\begin{align}
\label{eqn:ThetaDecorr-defn}
\thetaDecorr & \mydefn \thetaLS + \SigmaMat_\numobs^{-\frac{1}{2}}
\sum_{i = 1}^{\numobs} \w_i (\y_i - \myinprod{\x_i}{\thetaLS}).
\end{align}

Here the reader should recall our earlier definition
$\SigmaMat_\numobs \mydefn \sum_{i=1}^\numobs \x_i \x_i^\top$, and
throughout, we assume that the sample covariance $\SigmaMat_\numobs$
is invertible. The matrix $\SigmaMat_\numobs^{-\frac{1}{2}}$ denotes a
symmetric matrix square root of $\SigmaMat_\numobs^{-1}$.

Of course, there is an infinite family of estimators of the
form (<ref>), and the key question is how to
define the weight vectors.  In this paper, we propose an estimator in
which the sequence $\Nseq{\w}{i}$ is obtained by solving an
optimization problem that takes three inputs:

(i) the original data $\NseqNosub{(\x_i, \y_i)}{i}$,

(ii) a non-random scalar $\tuneParScaled_{\numobs} \in (0, 1]$,

(iii) a sequence of symmetric positive semidefinite matrices
$\Nseq{\XscaleMat}{i}$ such that $\XscaleMat_{i} \in \filtration_{i -
1}$ for each $i \in [\numobs] \mydefn \{1, \ldots, \numobs \}$.

In order to simplify notation, we adopt the shorthand $\z_i \mydefn
\XscaleMat_i^{-\frac{1}{2}} \x_i$. Moreover, for each index  we define the matrices

\begin{align*}
\Z_{i}^\top \mydefn \begin{bmatrix} \z_1 & \z_2 & \cdots & \z_i
\end{bmatrix}, \quad
\mbox{and} \quad \W_{i} \mydefn \begin{bmatrix} \w_{1} & \w_2 & \cdots
& \w_i
\end{bmatrix}
\end{align*}

We also define $\W_{0} = 0$ and $\Z_0 = 0$.  With these definitions,
the vectors $\Nseq{\w}{i}$ are obtained recursively by solving the
following convex program

\begin{align}
\label{eqn:Wn-expression-original}
\w_{i} & \mydefn \arg \min _{\w \in \real^\Dim} \left \{
\frobnorm{ \MyIdDim - \W_{i - 1} \Z_{i - 1} - \w \z_i^\top }^{2} +
\frac{\tuneParScaled_\numobs}{2} \|\w\|_2^2 \right \} \quad
\mbox{for $i = 1, 2, \ldots, \numobs$.}
\end{align}

Conveniently, this optimization problem has the explicit solution

\begin{align}
\label{eqn:Wn-expression}
\w_i & = \frac{(\MyIdDim - \W_{i - 1} \Z_{i -
1})\z_i}{(\tuneParScaled_\numobs/2) + \| \z_i \|_2^2}.
\end{align}

§ MAIN RESULTS

Having motivated and introduced the online debiasing approach, we now
turn to some theoretical guarantees that can be given for these
methods.  We begin in <Ref> by providing sufficient
conditions for the online debiasing estimator of <Ref>
to exhibit asymptotically Gaussian behavior
(<Ref>). In <Ref>, we
provide an asymptotically exact confidence region for $\thetastar$ as
well as an asymptotically exact confidence interval
(Proposition <ref>) for $\dir^\top \thetastar$, where
$\dir$ is an arbitrary fixed direction $\dir \in \real^\dims$.
In <Ref>—in particular,
see <Ref>—we complement these results by
providing minimax lower bounds on a family of Mahalanobis errors and
the length of confidence intervals.  These lower bounds apply to any
estimator for the stochastic regression model which does not know the
true value of the target parameter $\thetastar$ but may have the full
knowledge of how the data was collected. Finally,
in <Ref> we provide general strategies which can
be used to verify the conditions of <Ref>.  All
of our asymptotic statements assume that the dimension $\Dim$ is fixed
(constant) while the sample size $\numobs$ grows.

§.§ Asymptotic normality guarantees

The main result of this section is an asymptotic normality guarantee
for the proposed estimator (<ref>), where the
weight vectors are defined via the
recursion (<ref>).

We begin by stating our assumptions and providing some intuition about
their role in the theorem.

§.§.§ Assumption A

*  There are positive scalars $\sigma$ and
$\mompar$ such that the noise sequence $\Nseq{\error}{i}$ satisfies
the conditions $\Exs[\error_i \mid \filtration_{i - 1} ]
= 0$ and $\Exs[\error_i^2 \mid \filtration_{i - 1}] = \noisesd^2$ for all
$i \in [\numobs]$ and moreover

\begin{align*}
\max_{i \in [\numobs]} \; \Exs[\error_i^{2+\mompar} \mid
\filtration_{i - 1}] < \infty.
\end{align*}

The sequence of matrices $\{\SigmaMat_\numobs\}_{\numobs \geq1}$
satisfy the conditions  and

For each $\numobs$, the scalar $\tuneParScaled_\numobs > 0$ and
positive semidefinite matrices $\Nseq{\XscaleMat}{i}$ with
are chosen such that:

\begin{align}
\tag*{(a) Asymptotic negligibility:} \max \limits_{i \in [\numobs]}
\;\; \left \{ \frac{1}{\tuneParScaled_\numobs}
\myinprod{\x_i}{\XscaleMat_i^{-1} \x_i} \right \} &
\stackrel{\smallop{p}}{\longrightarrow} 0, \\
\tag*{(b) Vanishing bias:} \sqrt{\tuneParScaled_\numobs \log
\lambda_{\max}(\SigmaMat_\numobs) } \cdot \opnorm{\MyIdDim -
\W_\numobs \X_\numobs \SigmaMat_\numobs^{-\frac{1}{2}}} & \stackrel{\smallop{p}}{\longrightarrow}
0, \quad \mbox{and} \\
\tag*{(c) Variance stability:} \opnorm{\MyIdDim - \Nsum{i} \w_i
\x_i^\top \XscaleMat_i^{-\frac{1}{2}}} &
\stackrel{\smallop{p}}{\longrightarrow} 0.
\end{align}

Let us provide some intuition for the role of each of these
assumptions in the theorem.  First, Assumption <ref> is
quite simple: it imposes relatively mild moment conditions on the
noise variables.  Second, as discussed in the introduction,
Assumption <ref> is standard in guaranteeing the
consistency of the least squares estimate.  Both
Assumptions <ref> and <ref> are viewed as
mild conditions in the stochastic linear regression literature and are
satisfied by many practical models including those studied in the
papers [16, 18, 15, 6].
Note that Assumptions <ref> and <ref>
concern the regression model itself as opposed to the method: in
particular, they do not depend on the algorithm parameters
$\tuneParScaled_\numobs$ and $\Nseq{\XscaleMat}{i}$.

The more subtle requirements for our theorem to apply, which do
depend on the algorithm parameters, are stated in
Assumption <ref>.  We discuss the technical role of
these conditions in the comments after <Ref>, to
be stated momentarily.  In <Ref> to follow, we
provide concrete choices of the algorithm parameters
$\tuneParScaled_\numobs$ and $\Nseq{\XscaleMat}{i}$ that ensure that
Assumption <ref> holds.

Finally, it should be noted that Assumption <ref> is
weaker than the stability condition (<ref>).
Indeed, if the stability condition (<ref>) and
the growth condition <ref> are satisfied, then we may
take $\gamma_n = 1$, $\XscaleMat_i = \SigmaLim_\numobs$ and $\w_i =
\SigmaLim_\numobs^{-\frac{1}{2}} \x_i$. With these choices, the
conditions <ref> are automatically satisfied; moreover,
in this particular case, the online debiased estimator reduces to
ordinary least squares; see <Ref> for

With these preliminaries in
place, we are now equipped to state our main theorem on the online
debiasing estimator $\thetaDecorr$:

Let $\widehat{\noisesd}^2$ be any consistent estimator of
$\noisesd^2$.  Then under
, we have

\begin{align}
\label{eqn:asymp-normality}
\sqrt{\frac{\tuneParScaled_\numobs}{\widehat{\noisesd}^2}} \cdot
\SigmaMat_\numobs^{\frac{1}{2}}(\thetaDecorr - \thetastar) \indistrb
\Ncal(0, \MyIdDim).
\end{align}

We prove this theorem in <Ref>.

A few comments on this theorem are in order.  First, needing a
consistent estimate of the error variance $\noisesd^2$ is a mild
requirement.  For instance, under our conditions, the estimator

\begin{align*}
\widehat{\noisesd}^2 = \frac{1}{\numobs} \sum_{i = 1}^\numobs
\left(\y_i - \x_i^\top \thetaLS \right)^2
\end{align*}

is strongly consistent; see Lemma 3 in the paper [18]
for details.

A second important fact is that Assumption <ref> is
considerably weaker than the stability
condition (<ref>) required for asymptotic
normality of the OLS estimate.  To reinforce this point,
<ref> provides a detailed discussion of three
classes of problems for which OLS fails to be asymptotically normal
but the guarantee (<ref>) still holds for the
online debiasing estimator.

Of all the conditions of <ref>, verifying the
variance stability condition in part (c) of
Assumption <ref> is the most challenging, and our
arguments for doing so vary from problem to problem.  In
Corollaries <ref> and <ref>, we
verify the variance stability condition for multi-armed bandit
problems and autoregressive time series models, respectively. In
Corollary <ref>, we verify this condition for a
large class of problems satisfying a sufficient exploration
condition. In Section <ref> to follow, we argue
that when $\lambda_{\min}(\SigmaMat_\numobs) \geq \log^2(n)$, we can
always find choices of the tuning parameters $\gamma_n$ and
$\{\Gamma_i\}$ such that Assumption <ref> is
satisfied. Additionally, in absence of such lower bound on
$\lambda_{\min}(\SigmaMat_\numobs)$, we propose a data-augmentation
strategy that gets rids of this growth condition at the cost of
collecting additional $\log^2(n)$ many data-points.

Let us now discuss how Assumption <ref> enters the proof of
Theorem <ref>.  Our argument is based on the

\begin{align}
\sqrt{\tuneParScaled_\numobs} \cdot
\SigmaMat_\numobs^{\frac{1}{2}}(\thetaDecorr - \thetastar) & =
\Bias_\numobs + \ZeroMartin_\numobs, \quad \mbox{where} \\
\Bias_\numobs & \mydefn \sqrt{\tuneParScaled_\numobs} \cdot \left(
\MyIdDim - \W_\numobs \X_\numobs \SigmaMat_\numobs^{-\frac{1}{2}}
\right)(\thetaLS - \thetastar) \quad \mbox{and} \\
\ZeroMartin_\numobs & \mydefn \sqrt{\tuneParScaled_\numobs} \cdot
\sum_{i = 1}^\numobs \w_i \error_i.
\end{align}

By construction (and suggested by our notation), the term
$\Bias_\numobs$ corresponds to the bias in our estimate, a quantity
that must be shown to vanish in order for our claim to hold. In order
to do so, we first derive an upper bound on the norm
$\|\Bias_\numobs\|_2$.  The “vanishing bias” condition stated in
Assumption <ref>(b) enters in showing that, via our
choices of the tuning parameters $\tuneParScaled_\numobs$ and
$\Nseq{\XscaleMat}{i}$, this upper bound converges to zero in

The random vector $\ZeroMartin_\numobs$ defines a zero-mean
martingale, and our proof controls its behavior via a standard
martingale central limit theorem.  Doing so requires a Lindeberg type
condition on the weight vectors $\Nseq{\w}{i}$, as given in part (a)
of Assumption <ref>.  Moreover, it requires that the
conditional covariance of the martingale behave suitably, in which
context part (c) of Assumption <ref> enters.

§.§ Obtaining confidence regions and intervals

In this section, we use the online
debiasing procedure to obtain asymptotically exact confidence regions
and intervals.

§.§.§ Confidence region for $\thetastar$

First, for some user-defined level $\alpha \in (0,1)$, consider the
problem of finding a confidence region for $\thetastar$—that is, a
(random) set $\CIset_{1 - \alpha}$ that contains $\thetastar$ with
probability at least $1 - \alpha$.  We would like a set that is as
small as possible, asymptotically exact in the sense that its coverage
converges to $1-\alpha$.

Theorem <ref> allows us to construct such a set in
the following straightforward way.  For any $\alpha \in (0,1)$,
consider the subset of $\real^\dims$ given by

\begin{align*}
\CIset_{1 - \alpha} = \left\{ \thetaVec\in \real^{\dims} \mid
\frac{\tuneParScaled_\numobs}{\sigmahat^2} \cdot (\thetaDecorr -
\theta)^\top \SigmaMat_\numobs (\thetaDecorr - \theta) \leq
\chisq_{\dims, 1 - \alpha} \right\}
\end{align*}

where $\chisq_{\dims, 1 - \alpha}$ denotes the $(1 - \alpha$)-quantile
for a standard chi-squared distribution with degrees of freedom $d$.
From the result  of Theorem <ref>, we have the

\begin{align*}
\lim_{\numobs \rightarrow \infty} \Prob(\thetastar\in\CIset_{1 - \alpha}) = 1 - \alpha,
\end{align*}

In many applications, however, instead of a confidence region for the
full vector $\thetastar$, we are instead interested in obtaining a
confidence interval for the scalar quantity $\dir^\top \thetastar$,
where $\dir \in \real^\dims$ is a fixed direction.  It turns out that
Theorem <ref> no longer provides a straightforward
answer to this question.  In order to understand why, it is useful to
begin by following a naive line of reasoning that is incorrect
and then show how it can be fixed.

§.§.§ An incorrect argument

In order to obtain a confidence interval for $\dir^\top \thetastar$,
it might be tempting to “directly invert” the distributional
property (<ref>).  In particular, letting
denote the $1 - \tfrac{\alpha}{2}$ quantile of
the standard Gaussian distribution, we might claim that the interval

\begin{align}
\label{eqn:wrong-CI}
\left[ \myinprod{\dir}{\thetaDecorr} -
\tfrac{\sigmahat}{\sqrt{\tuneParScaled_\numobs}}
(\myinprod{\dir}{\SigmaMat_\numobs^{-1} \dir})^{\frac{1}{2}}
\NormalQuantile{1 - \alpha/2}, \qquad
\myinprod{\dir}{\thetaDecorr} +
\tfrac{\sigmahat}{\sqrt{\tuneParScaled_\numobs}} (\myinprod{\dir}
\NormalQuantile{1 - \alpha/2} \right],
\end{align}

is an asymptotically exact $1 - \alpha$ confidence interval for
$\dir^\top \thetastar$.

Unfortunately, the conclusion (<ref>) is based on faulty
logic, namely the assertion that the asymptotic
guarantee (<ref>) implies that

\begin{align}
\label{eqn:wrong-ASN}
\frac{\sqrt{\tuneParScaled_\numobs}}{\sigmahat} \cdot (\thetaDecorr
- \thetastar) \; - \; \Ncal(0, \SigmaMat_\numobs^{-1}) \indistrb 0.
\end{align}

It is now interesting to understand when the implication above follows from Theorem <ref>.

Under the stability condition $\mathbf{C}_n^{-1} \mathbf{S}_n \indistrb \Id$ for a sequence of deterministic matrices $\mathbf{C}_n$, the conclusion (<ref>) follows from the Theorem <ref> by Slutsky's theorem.
However, in the absence of this stability condition, the sample covariance matrix remains random and may depend on $\thetaDecorr$.

The following counterexample shows that when a sequence of random vectors $\mathbf{b}_n$ and matrices $\mathbf{A}_n$ are dependent, $\mathbf{A}_n \mathbf{b}_n \indistrb \mathcal{N}(0,\Id)$ does not imply $\mathbf{b}_n - \mathcal{N}(0,(\mathbf{A}_n^\top \mathbf{A}_n)^{-1}) \indistrb 0$ in general.

Let $X, Y$ be independent standard Gaussian random variables. Consider the vector $\mathbf{b} = (\textrm{sign}(X), \textrm{sign}(Y))^\top$, and a $2\times2$  matrix $\mathbf{A}$ with $\mathbf{A}_{11} = |X|/\sqrt{2}, \mathbf{A}_{12} = |Y|/\sqrt{2}, A_{21} = - |X|/\sqrt{2}, \mathbf{A}_{22} = |Y|/\sqrt{2}$. Simple calculations yields that $\mathbf{Ab} = (X + Y,\; X - Y)/\sqrt{2} \sim \mathcal{N}\left((0, 0)^\top, \Id \right)$. Note that each entry of $\mathbf{b}$ is $\pm1$ with probability $1/2$. Additionally,
the variable $\mathcal{N}( (0, 0)^\top (\mathbf{A}^\top \mathbf{A})^{-1})$ is a scale mixture of Gaussians — a continuous distribution —- and therefore does not match the distribution of $\mathbf{b}$.

In summary, we conclude that additional justification is needed to guarantee
the asymptotic validity of the confidence interval (<ref>) for the functional $\dir^\top \thetastar$ in general.

Nonetheless, there are certain special cases in which the
interval (<ref>) is a valid CI.  Concretely, suppose that
$\dir = \e_j$ is one of the standard coordinate basis vectors and that
$\SigmaMat_\numobs$ is diagonal as in the multi-armed bandit setting
studied in <Ref>.  In this case, the calculations
of <Ref> show that the interval (<ref>)
is valid.  More generally, given an arbitrary direction $\dir$, our
strategy will be to run a variant of online debiasing that effectively
reduces the problem to this favorable case.

§.§.§ Correct fixed-direction confidence intervals

Let us now describe the variant of online debiasing that can be used
to obtain asymptotically correct confidence intervals for fixed
directions.  Let $\dir \in \real^\dims$ be the direction of interest;
without loss of generality, we assume that $\enorm{\dir} = 1$.  We now
form an orthonormal basis of $\real^\dims$ with $\dir$ as its first
element—that is, a collection of orthonormal vectors $\{\dir_1 =
\dir, \dir_2, \ldots, \dir_\dims \}$.  Let $\Basis$ be the matrix with
$\dir_j^\top$ as its $j^{th}$ row.  Note that we have $\Basis \Basis^\top
= \MyIdDim$ and $\Basis^\top \e_1 = \dir$ by construction.  Using these
two properties, we can rewrite our model as

\begin{align*}
\dir^\top \thetastar = \e_1^\top \Basis \thetastar \;\;
\text{and} \;\; \y_i = \inprod{\Basis \x_i}{ \Basis \thetastar} +
\error_i \quad \text{for all} \;\; i = 1, \ldots, \numobs.
\end{align*}

Consequently, in this new basis, estimating the scalar $\dir^\top
\thetastar$ is same as estimating the first coordinate of transformed
vector $\Basis \thetastar$.

This fact allows us to define a variant of online debiasing that
supports asymptotically exact confidence intervals for $\dir^\top
\thetastar$.  In particular, let us introduce the notation

\begin{gather*}
\x_{\dir, i} = \Basis \x_i, \quad \X_{\dir, \numobs} = \X_\numobs
\Basis^\top, \quad \SigmaMat_{\dir, \numobs}
= (\Basis \SigmaMat_\numobs \Basis^\top)^{-1}, \quad \text{and} \quad
\InvSigmaV_{\dir, \numobs} = \SigmaMat_{\dir, \numobs}^{-1}.
\end{gather*}

Define the block diagonal matrix $\diagCov_{\dir,\numobs}$ as

\begin{align}
\label{eqn:Blk-diag-scaling-defn}
\diagCov_{\dir, \numobs}^{\frac{1}{2}} = \begin{pmatrix} \omega_{11}^{-\frac{1}{2}} & \mathbf{0}^\top \\
\mathbf{0} & \InvSigmaV_{22}^{-\frac{1}{2}} \end{pmatrix} \qquad  \text{where} \qquad
\InvSigmaV_{\dir, \numobs} = \begin{pmatrix} \omega_{11} &\InvSigmaV_{12} \\ \InvSigmaV_{21} & \InvSigmaV_{22} \end{pmatrix},
\end{align}

where the matrix $\InvSigmaV_{22}^{-\frac{1}{2}}$ denotes the symmetric square root of the matrix $\InvSigmaV_{22}^{-1}$.  Now consider the estimator

\begin{align}
\label{eqn:thetaBasis-defn}
\thetaBasis & \mydefn \thetaBasisLS + \OpTuning_\numobs \cdot \diagCov_{\dir,
\numobs}^{-\frac{1}{2}} \sum_{i = 1}^{\numobs} \w_i (\y_i -
\x_{\dir, i}^\top  \thetaBasisLS),
\end{align}

where $\OpTuning_\numobs \mydefn
\opnorm{\diagCov_{\dir,\numobs}^{\frac{1}{2}} \SigmaMat_{\dir,
\numobs}^{-\frac{1}{2}}}$ is a scalar which is at least one by
definition, and  is the OLS estimator
using the data $\{\y_i, \Basis \x_i \}_{i = 1}^\numobs$. We analyze
the behavior of $\thetaBasis$ under the following variant of
Assumption <ref>.

§.§.§ Assumption (A3)$^\prime$

*  For each $\numobs$, the scalar
$\tuneParScaled_\numobs > 0$ and positive semidefinite matrices
$\Nseq{\XscaleMat}{i}$ with  are
chosen such that:

\begin{align}
\tag*{(a) Asymptotic negligibility:} \max \limits_{i \in [\numobs]}
\;\; \left \{ \frac{1}{\tuneParScaled_\numobs}
\myinprod{\x_{\dir, i}}{\XscaleMat_i^{-1} \x_{\dir, i}} \right \} & \stackrel{\smallop{p}}{\longrightarrow} 0, \\
\tag*{(b) Vanishing bias:} \sqrt{\tuneParScaled_\numobs \log
\lambda_{\max}(\SigmaMat_{\dir,\numobs}) } \cdot \opnorm{ \tfrac{1}{\OpTuning_\numobs} \cdot \diagCov_{\dir, \numobs}^{\frac{1}{2}} \SigmaMat_{\dir, \numobs}^{-\frac{1}{2}}  -
\W_\numobs \X_{\dir, \numobs} \SigmaMat_{\dir, \numobs}^{-\frac{1}{2}}} & \stackrel{\smallop{p}}{\longrightarrow}
0, \quad \mbox{and} \\
\tag*{(c) Variance stability:} \opnorm{\MyIdDim - \Nsum{i} \w_i
\x_{\dir, i}^\top \XscaleMat_i^{-\frac{1}{2}}}  & \stackrel{\smallop{p}}{\longrightarrow} 0.
\end{align}

Under , given any consistent estimator
$\widehat{\noisesd}^2$ of $\noisesd^2$, the following interval is an
asymptotically exact $1 - \alpha$ confidence interval for $\dir^\top
\thetastar$

\begin{align}
\label{eqn:general-CI}
\left[ \e_1^\top \thetaBasis -
\tfrac{ \OpTuning_\numobs \sigmahat }{\sqrt{\tuneParScaled_\numobs}}
(\myinprod{\dir}{\SigmaMat_\numobs^{-1} \dir})^{\frac{1}{2}}
\NormalQuantile{1 - \alpha/2}, \  \e_1^\top \thetaBasis +
\tfrac{ \OpTuning_\numobs \sigmahat}{\sqrt{\tuneParScaled_\numobs}}
(\myinprod{\dir} {\SigmaMat_\numobs^{-1}\dir})^{\frac{1}{2}}
\NormalQuantile{1 - \alpha/2} \right].
\end{align}

See Section <ref> for the proof of this

A few comments regarding Proposition <ref> are in order.
Observe that
the length of the confidence intervals (<ref>)
matches the length of the confidence interval (<ref>)
up to a multiplicative factor $\OpTuning_\numobs$. Thus, it is interesting to understand
the value of the scalar $\OpTuning_\numobs$.  Note that $\OpTuning_\numobs^2 =
\lambda_{\max}(\diagCov_{\dir, \numobs}^{\frac{1}{2}} \SigmaMat_{\dir, \numobs}^{-1} \diagCov_{\dir, \numobs}^{\frac{1}{2}})$, and a little calculation yields

\begin{align*}
\diagCov_{\dir, \numobs}^{\frac{1}{2}} \SigmaMat_{\dir, \numobs}^{-1} \diagCov_{\dir, \numobs}^{\frac{1}{2}}
= \Id_{\dims} +   \begin{pmatrix}
0  & (\InvSigmaV_{22}^{-\frac{1}{2}} \InvSigmaV_{21} \omega_{11}^{-\frac{1}{2}})^\top \\
\InvSigmaV_{22}^{-\frac{1}{2}} \InvSigmaV_{21} \omega_{11}^{-\frac{1}{2}}  & {\mathbf{0}}_{\dims - 1}
\end{pmatrix}.
\end{align*}

Thus, we have

\begin{align}
\label{eqn:OpTuning-upperbound}
1 \leq \OpTuning_\numobs^2 = 1 + 2 \cdot\enorm{\InvSigmaV_{22}^{-\frac{1}{2}} \InvSigmaV_{21} \omega_{11}^{-\frac{1}{2}}}^2,
\end{align}

which yields $\OpTuning_\numobs \approx 1$ when the vector
$\InvSigmaV_{21} \approx 0$.  To gain further intuition on when
$\OpTuning_\numobs \approx 1$, let us assume $\dir = \e_1$, i.e., we
are interested in obtaining a confidence interval for the coordinate
$\thetastarScalar_1$.  In this case, a natural choice of the basis matrix is
$\Basis = \MyIdDim$, and as a result, we have $\InvSigmaV_{21} =
(\SigmaMat_\numobs^{-1})_{21} = (\SigmaMat^{-1}_{21}, \ldots,
\SigmaMat^{-1}_{\dims1})$. Recall that for $j \neq 1$, the entry
$\SigmaMat^{-1}_{j1}$ is proportional to the (empirical) partial
correlation coefficient between the first and the $j^{th}$ coordinate,
conditioned the remaining $\dims - 2$ coordinates of $\x$; meaning
that $\InvSigmaV_{21} \approx 0$ when the first coordinate of $\x$ has
small correlation with all linear functions of the other $\dims - 1$
coordinates of $\x$.

Finally, we point out that the assumption <ref>
is not significantly stronger than the original
assumption <ref>. To fix ideas, we again assume $\dir = \e_1$
and $\Basis = \Id$.
In that case, assumption <ref> and <ref>
only differ in the vanishing bias condition (b). Assuming
and condition <ref>(b) holds, we have

\begin{align}
&\phantom{=} \sqrt{\tuneParScaled_\numobs \cdot \log
\lambda_{\max}(\SigmaMat_{\e_1,\numobs})} \cdot \opnorm{ \tfrac{1}{\OpTuning_\numobs} \cdot
\diagCov_\numobs^{\frac{1}{2}} \SigmaMat_\numobs^{-\frac{1}{2}}-
\W_\numobs \X_\numobs \SigmaMat_\numobs^{-\frac{1}{2}} } \notag \\
\sqrt{\tuneParScaled_\numobs \cdot \log
\lambda_{\max}(\SigmaMat_\numobs) } \cdot \opnorm{ \tfrac{1}{\OpTuning_\numobs}
\cdot\MyIdDim-
\W_\numobs \X_\numobs \SigmaMat_\numobs^{-\frac{1}{2}} } \notag \\
& \quad \quad \quad \quad +
\sqrt{\tuneParScaled_\numobs \cdot \log
\lambda_{\max}(\SigmaMat_\numobs) } \cdot (1 +
\tfrac{1}{\OpTuning_\numobs} \cdot \opnorm{\diagCov_\numobs^{\frac{1}{2}}
\SigmaMat_\numobs^{-\frac{1}{2}}}) \notag \\
&= \smalloP(1) + \smalloP(1) \stackrel{\smallop{p}}{\longrightarrow} 0
\label{eqn:A3_to_A3prime} .
\end{align}

The first step  uses the fact $\lambda_{\max}(\SigmaMat_\numobs) = \lambda_{\max}(\SigmaMat_{\dir, \numobs})$ for any basis matrix $\Basis$.
The second step uses the vanishing bias condition <ref>(b),
the fact that the dimension $\dims$ is
fixed, and the upper
bound .

§.§.§ Comparison with least squares confidence intervals

An interesting consequence of the results so far, is that we can construct a confidence intervals based on any consistent estimator $\thetahat$ of $\thetastar$. To illustrate this, we focus on constructing a confidence interval for the mean of the first arm ($\theta_1$) of a multiarmed bandit problem. Theorem <ref> and Proposition <ref> ensures that

$$[\widehat{\theta}^{OD}_{1} - \sqrt{1/\gamma_n} \cdot \frac{z_{1 - \alpha/2}}{\sqrt{\SigmaMat_{11}}}, \;\;\widehat{\theta}^{OD}_1 + \sqrt{1/\gamma_n} \cdot \frac{z_{1 - \alpha/2}}{\sqrt{\SigmaMat_{11}}}]$$

is an asymptotically exact $1 - \alpha$ confidence interval for $\theta^\star$. As an immediate consequence, we have that for any weakly consistent estimator $\widehat{\theta}_1$ for $\theta^\star_1$ the confidence interval

$$[\widehat{\theta}_{1} -  \sqrt{1/\gamma_n} \cdot \frac{z_{1 - \alpha/2}}{\sqrt{\SigmaMat_{11}}}, \;\;\widehat{\theta}_1 + \sqrt{1/\gamma_n} \cdot \frac{z_{1 - \alpha/2}}{\sqrt{\SigmaMat_{11}}}]$$

is also an asymptotically exact $1 - \alpha$ confidence interval. Accordingly, when the OLS estimator $\widehat{\theta}^{OLS}$ is consistent, the above argument also allows us to construct
an asymptotically valid confidence interval for the least squares estimator.

Interestingly, while the width of a valid OLS interval is tightly constrained by its own bias, the width of a valid online debiasing interval can be significantly smaller. In particular, one can construct confidence intervals which are $\sqrt{\log(n)}$ times smaller in width than the ones obtained from lest square estimator; see Appendices <ref> and <ref> for details.

§.§ Minimax lower bounds

Thus far, we have derived two guarantees for online debiasing
procedures: asymptotic normality in Theorem <ref>
along with confidence intervals in Proposition <ref>.
It is natural to wonder in what sense these guarantees are optimal.
Accordingly, this section is devoted to lower bounds that apply to the
performance of any estimator $\thetahat$.  These bounds are
derived within the classical minimax framework and cover two
particular risk measures.

Our first risk measure involves the Mahalanobis pseudometric:
given an arbitrary positive semi-definite matrix $\mahal$, possibly
random, this pseudometric[We parameterize the Mahalnobis
pseudometric slightly differently than standard definitions, using
$\mahal$ as opposed to its inverse for the quadratic form.  This is
only for notational ease when $\mahal$ has a non-trivial null space.]
is given by

\begin{align}
\label{EqnMahal}
\lVert\thetahat - \thetastar\rVert_\mahal \mydefn \lVert
\mahal^{\frac{1}{2}}(\thetahat - \thetastar) \rVert_2,
\end{align}

and we provide lower bounds on the squared form of this pseudometric
in part (a) of Theorem <ref>, below.  Notably,
our analysis allows for the matrix $\mahal$ to also depend on the
dataset $\{ \x_i, \y_i \}_{i = 1}^\numobs$ itself, so that for example, setting $\mahal =
\SigmaMat_\numobs$ is a valid choice.

Our second risk measure corresponds to the length of a two-sided
confidence interval. For a given vector $\dir \in \real^\Dim$ and
significance level $\alpha \in (0,1)$, let  be any level
$\alpha$ confidence interval for the scalar
$\myinprod{\dir}{\thetastar}$, so that by definition, we have

\begin{align}
\label{eqn:valid-CI}
\Prob_{\thetastar} \big[\myinprod{\dir}{\thetastar} \in \CI_{\alpha, \dir}
\big] & \geq 1 -  \alpha \qquad \text{for
all } \;\; \thetastar \in \real^\dims.
\end{align}

We are interested in finding the smallest such confidence interval,
and part (b) of Theorem <ref> provides a lower
bound on its length ${|\CI_{\alpha, \dir}| \mydefn u_\alpha -
\ell_\alpha}$.

Our bounds apply to any estimator $\thetahat$, meaning a measurable
function of the data as well as the data collection process.  The data
collection process is summarized by a collection of (potentially
randomized) selection algorithms, each of the form
which take the observed data $\{ (\x_j, \y_j)
\}_{j=1}^{i-1}$ up to time $i$ and output a new observation
$\x_{i}$. With a slight abuse of notation, we refer to

as the selection algorithm of the data collection process.

Fix any selection algorithm $\Selections_\numobs$.  Under the linear
model (<ref>) with i.i.d. Gaussian noise
and data collected using
$\Selections_\numobs$, the following claims hold:

* For any (possibly random) matrix $\mahal$ such that $\Exs[\tr(
\SigmaMat_\numobs^{-1} \mahal)]$ is finite, we have

\begin{align}
\label{EqnMinimaxMSE}
\inf_{\thetahat} \sup_{\thetastar \in \real^\dims}
\E \lVert
\thetahat - \thetastar \rVert^{2}_{\mahal} \geq \sigma^2 \Exs[\tr(
\SigmaMat_\numobs^{-1} \mahal) \big],
\end{align}

where the infimum is taken over any estimator $\thetahat$, potentially
depending on $\Selections_\numobs$, in addition to the data.

* There is a universal constant $C > 0$ such that for any pair
$(\numobs, \Dim)$ with $n \geq d^3/C$ and $d \geq 2$, there exists a selection
algorithm $\Psi_\numobs$ and direction $\dir \in \reals^d$ such that

\begin{align*}
\inf_{\cE} \sup_{\thetastar} \E\left[\frac{(\ip{\dir, \thetastar} -
\cE)^2}{\norm{\dir}^2_{\SigmaMat_n^{-1}}}\right] \geq C \cdot d
\sigma^2 \log\left(n\right)\,,
\end{align*}

where the infimum is taken over all $\filtration_n$-measurable
estimates $\cE$ of the scalar $\ip{\dir,\thetastar}$

* Suppose that $\Exs[\SigmaMat_n]$ exists and is invertible.  Then
for any direction $\dir \in \real^\Dim$ and scalar  we have

\begin{align}
\label{EqnMinimaxCI}
\inf_{\CI_{\alpha, \dir}} \sup_{\thetastar \in \real^\dims} \E \big[
\lvert \CI_{\alpha, \dir} \rvert \big] \geq 2\noisesd \cdot z_{1-
\alpha/2} \cdot\Exs\left\{\left(\dir^\top \SigmaMat_n^{-1}
\dir\right)^\frac{1}{2}\right\},
\end{align}

where the infimum is taken over any procedure, potentially depending
on $\Selections_\numobs$ in addition to the data, that returns valid
level $1- \alpha$ confidence intervals $\CI_{\alpha, \dir}$ for
$\ip{\dir,\thetastar}$ (cf. definition (<ref>)).

We provide the proofs of parts (a), (b) and (c)
of <Ref>
in <Ref>,
and <Ref>, respectively.

§.§ Comments on part (a): instance-dependent lower bound on MSE

In order to gain intuition for the MSE bound in part (a), it is
helpful to begin with the simplest case—that is, the non-adaptive
setting.  Consider the classical problem of fixed design linear
regression, in which the covariates (and hence $\SigmaMat_\numobs$)
are viewed as fixed, and the additive noise is zero-mean Gaussian with
variance $\noisesd^2$.  In this case, the standard OLS estimate
$\thetaLS$ follows the Gaussian distribution $\normal(\thetastar,
\noisesd^2 \SigmaMat_\numobs^{-1})$ for any sample size $\numobs$.
Consequently, for any fixed matrix $\mahal$, we have the equality

\begin{align}
\E \lVert \thetahat - \thetastar \rVert^{2}_{\mahal} = \sigma^2
\Exs[\tr(\SigmaMat_\numobs^{-1} \mahal)].
\end{align}

This simple calculation shows that the lower
bound (<ref>) is unimprovable in general.

Of course, the more substantive content of
Theorem <ref>(a) lies in the fact that it
allows for adaptive data collection, along with potentially random
choices of $\mahal$.  One interesting choice is the random matrix
$\mahal = \SigmaMat_\numobs$, for which the
bound (<ref>) guarantees that
It is worth comparing this lower bound to <Ref>.
From the arguments used to prove this theorem, and under a
mildly stronger version of
Assumption <ref>—which the convergence in distribution
conditions are replaced by convergence in $L_1$—it can be shown that

\begin{align}
\label{EqnL1-convergence}
\lim_{\numobs \rightarrow +\infty} \tuneParScaled_\numobs
\|\thetaDecorr - \thetastar\|_{\SigmaMat_\numobs}^2 & = \noisesd^2
\dims.
\end{align}

See <Ref> for the details of this argument.

As discussed in <Ref> to follow, in many practical
problems of interest, the tuning parameter $\tuneParScaled_\numobs$
typically scales logarithmically in the sample size $\numobs$, and
also our choice of the tuning parameters ensure that the
aforementioned stronger version of assumption <ref> is
satisfied. Consequently, the result (<ref>), when
combined with the lower bound (<ref>), shows that the
online debiasing procedure is instance-optimal up to logarithmic

§.§.§ Comments on part (b): minimax lower bound on MSE

Part (b) of Theorem <ref> provides a minimax
lower bound on the MSE. It shows that in the adaptive setting, the
scaling of $\gamma_n$ needs to at least $\log(\numobs)$, and moreover,
this logarithmic scaling is unavoidable when dimension $d \geq 2$. As
discussed in the last comment, for many problems, we can take
$\gamma_n$ arbtitarily close to $\log(\numobs)$; thus, we conclude
that the online debiased estimator is miniamx optimal when $d \geq
2$. Finally, we point out that the lower bound result is not
true when $d = 1$. In this case the logarithmic dependence on $n$
becomes doubly logarithmic, as is consistent with the law of the
iterated logarithm that underlies this behavior in that case.

§.§ Comments on part (c): bounds on lengths of CIs

<Ref>(c) provides a lower bound on the width
of any confidence interval for the scalar $\dir^\top \thetastar$ that
is valid when the data set is collected in an adaptive manner. To the
best of our knowledge, this is the first result providing a lower
bound on the width of confidence intervals in an adaptive setting.

§.§ Choices of the tuning parameters

Let us now return to the practical issue of choosing the tuning
parameters $\tuneParScaled_\numobs$ and $\Nseq{\XscaleMat}{i}$ of our
debiasing procedures.  In particular, these parameters must be chosen
appropriately so as to ensure that either
Assumption <ref>, or its variant in
Assumption <ref>, is satisfied.

We analyze practical default choices that are based upon on a
deterministic matrix $\SigmaMatLB_\numobs$ that acts as a lower bound
on the sample covariance matrix $\SigmaMat_\numobs$.  Let $\{
\SigmaMatLB_\numobs \}_{\numobs \geq 1}$ be a sequence of
$\Dim$-dimensional diagonal matrices with nonnegative entries
such that

\begin{align}
\label{eqn:SigmaLB-conditions}
\opnorm{ \SigmaMatLB_\numobs^{\frac{1}{2}}
\diag(\SigmaMat_\numobs^{-1}) \SigmaMatLB_\numobs^{\frac{1}{2}} } =
\BigoP(1) \quad \text{and} \quad \lambda_{\min}(\SigmaMatLB_\numobs)
\almostSurely \infty.
\end{align}

For a given $\numobs$, we define a collection of (diagonal) scaling

\begin{align}
\label{eqn:scaling_mat-choice}
\XscaleMat_{i, \numobs} & \mydefn \max \left \lbrace \diag\left(
\SigmaMat_i^{-1} \right)^{-1}, \SigmaMatLB_\numobs \right \rbrace,
\end{align}

where $\max\{\cdot, \cdot\}$ denotes the element-wise maximum
operator.[The choice (<ref>) of
scaling matrix is especially easy to understand for multi-armed bandit
problems, where the scaling matrix $\XscaleMat_{i, \numobs}$ can be
written as $\XscaleMat_{i, \numobs} = \max \left\lbrace \SigmaMat_i,
\;\; \SigmaMatLB_\numobs \right\rbrace$.  Assuming that $\SigmaMat_i =
\max\{\SigmaMatLB_\numobs, \SigmaMat_i \}$ for large value of $i$, we
see that the tuning parameter $\XscaleMat_{i, \numobs}$ is the sample
covariance matrix up to time $i$.  This assumption indeed holds for
Corollaries <ref>– <ref> to be
presented in the sequel.]

We point out that it is relatively straightforward to find a diagonal
matrix $\SigmaMatLB_\numobs$ satisfying the
condition (<ref>) as long as
$\lambda_{\min}(\SigmaMat_\numobs) \rightarrow \infty$ almost
surely. For simplicity, let us assume that the covariates $\x_i$ are
uniformly bounded and that the minimum eigenvalue of the matrix
$\SigmaMat_\numobs$ is lower bounded as  with high probability
(see <Ref> for one sufficient condition
for this bound to hold).  Then, with our recommended default choices

and $\gamma_n = \frac{1}{\log(n)\log\log(n)}$ the
condition (<ref>) is satisfied.

Consider the solutions $\Nseq{\w}{i}$ obtained from the optimization
problem (<ref>) using parameters
$\NseqNosub{\XscaleMat_{i, \numobs}}{i}$
defined in equation (<ref>) and (non-random) $\tuneParScaled_\numobs > 0$. Then we have the operator norm bound

\begin{align}
\label{EqnVanishBiasProp}
\opnorm{ \Id - \W_\numobs \X_\numobs
\SigmaMat_\numobs^{-\frac{1}{2}} } = \BigoP(\Dim^2).
\end{align}

In particular, if $\tuneParScaled_\numobs = \smalloP
\left(\Dim^4\log \lambda_{\max}(\SigmaMat_\numobs) \right)$,
the vanishing bias and asymptotic negligibility conditions in
Assumption <ref> are satisfied.

See <ref> for the proof of this

§.§.§ Sharper bound for multi-armed bandits
The dimension
dependence of the upper bound (<ref>) can be removed
in many concrete applications in which we have additional information
about the data generating process.

As one concrete example, in the multi-armed bandit model of the sequel
(<ref>), the upper bound can be sharpened to

\begin{align}
\label{eqn:bias-bound-improved}
\opnorm{\Id - \W_\numobs \X_\numobs \SigmaMat_\numobs^{-\frac{1}{2}} }
= \BigoP(1).
\end{align}

See the end of Appendix <ref> for a proof of
this claim, and see the proofs of the
Corollaries <ref>, <ref>,
and <ref> in the sequel for more details.

§.§.§ Verifying the growth condition on $\lambda_{\min}(\SigmaMat_\numobs)$

In the absence of any additional assumption on the data collection
method, it may be difficult to verify the lower bound

A simple fix to this problem is to collect $\log^2(\numobs)$ many
additional data points with $\x_i$ chosen uniformly at random from a
$d$-dimensional unit sphere and append them to the original
dataset. Then the condition
$\lambda_{\min}(\SigmaMat^{\textrm{new}}_\numobs) \geq
\log^2(\numobs)$ is automatically satisfied for
$\SigmaMat^{\textrm{new}}_\numobs$ — the sample covariance matrix of
the new augmented dataset.

§ APPLICATIONS

We next illustrate the concrete consequences of our results.
Sections <ref> and <ref> are devoted to
multi-armed bandit problems and autoregressive time series models,
respectively, while <ref> discusses active
learning with exploration.  We end each section with an empirical
evaluation of online debiasing.  Specifically, we compare the
confidence interval (CI) coverage and width of four methods: our
online debiasing estimator (<ref>),
OLS (<ref>) with standard but potentially invalid
Gaussian intervals, the
of 6, and a valid CI based on the
concentration inequality of 2.  We highlight
that the CIs for OLS are based on the distributional assumption
. This property, while true when the covariates $\{\x_i
\}_{i = 1}^n$ are selected in a non-adaptive manner, need not hold for
adaptively collected
covariates [6, 30], and as a
consequence, the corresponding CIs need not give the correct coverage.
Meanwhile, the valid concentration inequality-based
intervals [2] are guaranteed to provide at least
the nominal coverage but are often unnecessarily wide.

§.§ Multi-armed bandits

Consider a multi-armed bandit with $\Dim$ arms indexed by the set
At each time $i \in
[\numobs]$, a bandit algorithm selects an arm $k_i \in [\Dim]$
and observes the reward

\begin{align}
\label{model:bandits}
\y_i = \myinprod{\basis_{k_i}}{\armMeans} + \epsilon_i,
\end{align}

where $\basis_{k_i}$ is the $k_i^{th}$ basis vector in dimension
$\Dim$ and $\armMeans \in \real^\Dim$ is the vector containing the
mean rewards of $\Dim$ arms.  We assume that the noise sequence
$\Nseq{\error}{i}$ satisfies Assumption <ref>.  Notably,
the multi-armed bandit model (<ref>) is a special case
of the adaptive linear regression model (<ref>)
with $\x_i = \basis_{k_i}$ for each $i \in [\numobs]$.

Since the bandit observation model (<ref>) has a simple
linear form, the OLS solution $\thetaLS$ is a standard estimate of the
reward vector $\thetastar \in \real^\Dim$.  As we mentioned earlier,
the behavior of the OLS estimate depends on the stability of the
matrix $\SigmaMat_\numobs$; see the covariance stability
condition (<ref>).  In the
paper [6], the authors conjectured based on
empirical evidence that for various popular data selection algorithms,
including the Upper Confidence Bound (UCB), Thompson Sampling, and
$\greedyEps$-greedy algorithms (see the
book [20]), the stability
condition (<ref>) is not satisfied when there
are multiple optimal arms. In recent work, Zhang et
al. [30] established the validity of this
conjecture for the two-armed bandit problem: when the two means are
equal, then the OLS estimate fails to have a Gaussian limiting

In sharp contrast to these negative results for ,
<ref> to follow guarantees that the online debiasing
estimator (<ref>) is asymptotically normal under
a mild assumption on the minimum number of times that each arm is
pulled.  More precisely, for each arm $k \in [\Dim]$ and round $i \in
[\numobs]$, let $\armCount_{k, i}$ denote the number of times $k$ is
pulled in the first $i$ rounds, and define the minimum
and maximum  arm counts.  Then the scaled
sample covariance is a $\Dim \times \Dim$ diagonal matrix, in which
the $k^{th}$ diagonal entry corresponds to the number of times that
arm $k$ is pulled within the first $i$ rounds:

\begin{align}
\label{eqn:bandit-notations}
\armCountMat_i & = \diag \left( \armCount_{1, i}, \ldots,
\armCount_{\Dim, i} \right).
\end{align}

We assume a lower bound on the minimum number of times that each arm
is pulled—namely,

\begin{align}
\label{eqn:bandits-minimum-arm-LB}
\armCount_{\min} \geq ( \log \numobs)^2.
\end{align}

Moreover, we implement the debiasing
estimate (<ref>) with the choice of tuning

\begin{align}
\label{eqn:bandits-tuning}
\tuneParScaled_\numobs = \frac{1}{(\log \numobs) \cdot \log
\log(\numobs)} \quad \text{and} \quad \XscaleMat_{i, \numobs} =
\max \left\lbrace \armCountMat_i, \;\; (\log \numobs)^2 \cdot
\Id_{\Dim} \right \rbrace,
\end{align}

where $\max \left\lbrace \cdot, \cdot \right\rbrace$ denotes the
element-wise maximum operator.

Suppose the minimum arm pull
condition (<ref>) and the moment
condition <ref> are valid. Then, given any consistent
estimate $\sigmahat^2$ of the error variance $\errorstd^2$, the
estimate $\thetaDecorr$ obtained using the tuning parameter
choices (<ref>) satisfies

\begin{align}
\label{EqnBanditsCorr}
(\widehat{\noisesd}^2 \cdot (\log \numobs) \cdot \log\log(\numobs) \big)^{-1/2}
\cdot \armCountMat_\numobs^{\frac{1}{2}} \left( \thetaDecorr -
\armMeans \right) \indistrb \Ncal(0, \MyIdDim).
\end{align}

See <ref> for the proof of this claim.
<ref> also enables us to construct asymptotically exact
confidence regions for $\thetastar$. Moreover, the sample covariance matrix $\SigmaMat_n$
is diagonal, and as a result, we can also construct confidence intervals of the coordinates
$\theta_i^*$; see the proof of Proposition <ref> for details. Finally,
for a direction $\dir$ which is not a standard basis direction,
we can obtain an asymptotically exact $1 - \alpha$ confidence interval
of $\dir^\top \thetastar$ using Proposition <ref>; see the comments following
Corollary <ref> for further details.

§.§.§ Numerical experiment

illustrates the performance of online debiasing with bandit tuning (<ref>).
we consider a two-armed bandit problem (<ref>) with arm-mean vector
$\thetastar = (0.3, 0.3)^\top$ and i.i.d. standard normal error
$\{ \error_i \}_{i = 1}^\numobs$.
The covariates $\{ \x_i \}_{i = 1}^\numobs$  were
generated using the Thompson sampling algorithm [24], and we consider confidence intervals (CIs) for $\theta_1^*$.

We observe first that online debiasing provides appropriate coverage for all confidence levels.
Meanwhile, the lower tail interval severely undercovers, and  W-decorrelation undercovers for both tails despite having larger widths than online debiasing.
Finally, the concentration CI provides 100% coverage for all  confidence levels but yields intervals uniformly larger than the online debiasing CIs.
In Appendix <ref>, we present analogous results for two other popular multi-armed bandit algorithms, the upper confidence bound (UCB) and $\greedy$-greedy algorithms.

Average coverage and width of confidence intervals for $\theta_1^*$ across 1000 independent replications  of a multi-armed bandit experiment (<ref>) with $\thetastar \equiv (\theta_1^*, \theta_2^*) =  (0.3, 0.3)^\top$. The covariates $\{\x_i\}_{i = 1}^{1000}$
were selected using the UCB algorithm [20], and the error bars represent $\pm 1$ standard error.
Left and Center: Coverage of one-sided $1 - \alpha$
intervals for $\theta_1^*$. Right: Width of two-sided $1 - \alpha$ intervals for $\theta_1^*$.
The first row corresponds to data generated following the UCB algorithm, the second row corresponds to performances when additional data of size $\log^2(n)$ is added to the data; see Section <ref>. We observe that the online debiased estimator provide qualitatively similar performances in the both cases.
See Section <ref> for additional simulations.

§.§ Autoregressive time series model

Our next example involves estimating the parameters of an autoregressive
time series model.  It is well-known that the estimate can exhibit non-Gaussian limit behavior for
versions of such processes that are unstable [18].  In order to focus attention on the key
issues, we restrict ourselves here to the simple case of a scalar
autoregressive process.

More precisely, given the initial point $\y_0 = 0$ and an unknown
scalar $\thetastarScalar \in (-1, 1]$, consider a stochastic process
generated by the first-order autoregression

\begin{align}
\label{eqn:autoreg-model}
\y_i = \thetastarScalar \y_{i - 1} + \error_i \quad \mbox{for $i =
1, \ldots, \numobs$.}
\end{align}

We assume that the noise sequence $\Nseq{\error}{i}$ consists of
i.i.d. standard normal random variables.  Note that the
autoregression (<ref>) is a special case of the
stochastic linear regression model (<ref>), in
particular one with $x_i = \y_{i - 1}$ for all $i \in [\numobs]$. An
especially interesting instantiation of the
autoregression (<ref>) is obtained by setting
$\thetastarScalar = 1$. Such a process is a special case of a
unit root autoregression, a class of models that play an
important role in econometric time series analysis [4].

With the choice $\thetastarScalar = 1$, the
process (<ref>) is a random walk and so has
a variance that grows linearly with time.

Moreover, by an application of Donsker's theorem (cf. Example
3 in the paper [18]), we have

\begin{align}
\label{eqn:Var-scaling-AR1}
\frac{1}{\numobs^2} \Nsum{i} x_i^2 \mydefn \frac{1}{\numobs^2}
\Nsum{i} \y_{i - 1}^2 &\indistrb \int_{0}^1 \wiener^2(t) dt, \qquad
\mbox{and} \\
\sqrt{\Nsum{i} y_{i - 1}^2} \cdot (\thetaLSScalar - \thetastarScalar)
& \indistrb \frac{\wiener^2(1) - 1}{2 \int_{0}^1 \wiener^2(t) dt},
\notag
\end{align}

where $\wiener(t)$ denotes the standard Wiener process (see the
paper [28] for details).  Put simply, in the autoregressive time series
model (<ref>) with $\thetastar = 1$ the
stability condition (<ref>) is not satisfied,
and the distribution of the OLS estimate $\thetaLSScalar$ is not
asymptotically normal.

In contrast to this negative result for the OLS estimate, we can show
that the debiasing estimate $\thetaDecorrScalar$, after suitable
centering and scaling, does indeed converge in distribution to a
standard Gaussian.  Our result is based on the tuning parameters
and scaling matrices chosen as

\begin{align}
\label{eqn:autoregressive-tuning}
\tuneParScaled_\numobs = \frac{1}{(\log \numobs)\cdot \log\log(\numobs)}, \quad
\mbox{and} \quad {\XscaleMat_{i, \numobs} = \max \left \{ (\log \numobs)^2 y_{i - 1}^2, \;\; \sum_{j =1}^{i-1} y_j^2 \right \}}.
\end{align}

Given a sequence $\Nseq{\y}{i}$ generated from the autoregressive
model (<ref>), the estimate
$\thetaDecorrScalar$ (<ref>) obtained with the
tuning parameters (<ref>) satisfies

\begin{align}
\label{EqnAutoRegClima}
\sqrt{\frac{\Nsum{i} \y_{i - 1}^2}{(\log \numobs) \cdot \log\log(\numobs)}}
\cdot ( \thetaDecorrScalar - \thetastarScalar ) \indistrb \Ncal(0,
\end{align}

See <ref> for the proof of this
<ref> enables us to construct asymptotically exact
confidence intervals for $\thetastarScalar$. We also reiterate that the above result holds for any ${\thetastar \in (-1, 1]}$.

§.§.§ Numerical experiment

<ref> illustrates the performance of online
debiasing with autoregression
tuning (<ref>).  Here our data is generated
from the time series model (<ref>) with $\thetastar
= 1$.  We again find that online debiasing provides appropriate
coverage for all confidence levels.  Meanwhile, the lower tail
interval exhibits
severe undercoverage, and W-decorrelation exhibits ranges of undercoverage for both tails.  Finally, the concentration-based CI again
provides 100% coverage for all confidence levels, at the expense of
interval lengths that are uniformly longer than the online debiasing

Average coverage and width of confidence intervals for $\thetastar=1$ across 1000 independent replications of an autoregressive time series experiment (<ref>). The error bars represent $\pm1$ standard error.
Left and Center: Coverage of one-sided $1 - \alpha$
intervals for $\thetastar$. Right: Width of two-sided $1 - \alpha$ intervals for $\thetastar$.
See Section <ref> for details.

§.§ Active learning with exploration

In our third example, we focus on the case where the covariates
$\Nseq{\x}{i}$ are generated using any algorithm satisfying a sufficient exploration property.

We say that a selection algorithm $\Selections_\numobs =
\Nseq{\selection}{i}$ admits a $\greedy$-exploration property

\begin{align}
\label{eqn:selection-with-exploration}
\x_i \mydefn
\begin{cases}
\algoDir_i & \text{with probability $1 - \greedyEps_i$ for some $\algoDir_i \in
\filtration_{i - 1}$, and} \\
\randDir_i & \mbox{with probability $\greedyEps_i$ for some $\randDir_i$
independent of $\filtration_{i - 1}$.}
\end{cases}
\end{align}

Here the exploration probability sequence $\Nseq{\greedyEps}{i}$
consists of nonnegative scalars in the interval $(0,1)$, and the
vectors $\Nseq{v}{i}$ are i.i.d. random vectors such that

\begin{align}
\label{eqn:exploration-lower-bound}
\Exs[\randDir_i \randDir_i^\top] \succeq \ExsXsqLB
\qquad \text{where}  \quad \ExsXsqLB \succeq \mathbf{0}
\end{align}

In words, the selection algorithm $\selection_i$ behaves as follows:
with probability $1 - \greedyEps_i$, it chooses vector $\algoDir_i$
based on the previous data points $\{(\x_j, \y_j)\}_{j =1}^{i-1}$, and
with probability $\greedyEps_i$, it chooses a random direction
$\randDir_i$, independent of the previous data points.

§.§.§ Example: $\greedy$-greedy linear bandits

Let us briefly consider a concrete instance of a selection algorithm
$\Nseq{\selection}{i}$ that is of the $\greedyEps$-greedy type.  In
the linearly parameterized bandit problem, at each time $i \in
[\numobs]$, an algorithm $\selection_i$ chooses an action vector
$\x_i$, usually lying within some bounded set $\actionset_i$, and
obtains a reward .
A popular and simple strategy for regret minimization is a special
case of the $\greedyEps$-greedy selection
algorithm [20].  For linearly parameterized
bandits, the selection algorithm $\selection_i$ chooses

\begin{align}
\label{eqn:greedy-selection}
\x_i \;\; &
\begin{cases}
\in \arg \max \limits_{\x\in \actionset_i} \quad
\myinprod{\x_i}{\thetaRidge^{(i - 1)}} & \mbox{with probability
$1 - \greedyEps_i$, and} \\
\sim \Unif(\actionset_i) & \mbox{with probability
\end{cases}
\end{align}

where $\thetaRidge^{(i - 1)}$ denotes the ridge regression estimator
based on all data observed up to stage $i - 1$, i.e., the collection
of covariate-response pairs $\{ (\x_j, \y_j)\}_{j=1}^{i-1}$.  Put
simply, with probability $1-\greedyEps_i$, the selection algorithm
$\selection_i$ chooses an optimal action given data collected so far
(exploitation), and with probability $\greedyEps_i$, the algorithm
randomizes uniformly amongst its choices (exploration).  In the more
general setting (<ref>) considered
here, it is not necessary to select the optimal action in the
exploitation step.  Rather, our result holds also when an arbitrary,
$\filtration_{i-1}$-measurable choice is made in the first part
of <ref>, as in the EXP3 or UCB algorithms with
exploration.  See the book [20] for more
details.  $\diamondsuit$

Returning to our general
setting (<ref>), we now state a
guarantee for selection algorithms with $\greedyEps$-exploration.  As
is standard in the bandit literature, we assume that the covariates
are uniformly bounded, so that there exists a scalar $K$ satisfying

\begin{align}
\label{eqn:bounded-covariates}
\| \x_i \|_2 \leq K \quad \mbox{for all $i \in [\numobs]$.}
\end{align}

See our discussion following the corollary for how this condition can
be relaxed. In addition, we impose a sufficient exploration
condition, meaning a lower bound on the magnitude of the exploration
probabilities, of the form

\begin{align}
\label{eqn:suff-exploration}
\sum_{i = 1}^\numobs \greedyEps_i \geq \frac{ \Exs[ \max_{i \in
[\numobs]} \| \x_i \|^2_2] }{\lambda_{\min}(\ExsXsqLB)} (\log
\numobs)^2 % \quad \mbox{for some $\pardelta > 0$,}
\end{align}

where the reader should recall that the matrix $\ExsXsqLB$ was defined
in <ref>.  We implement the debiasing
estimate (<ref>) with the choice of tuning

\begin{align}
\label{EqnTuningGreedy}
\XscaleMat_{i, \numobs} = \sum_{j = 1}^\numobs \greedyEps_j \ExsXsqLB
\quad \qtext{and} \quad \tuneParScaled_\numobs =
\frac{1}{(\log\numobs) \cdot \log\log(\numobs)}.
\end{align}

Suppose that Assumptions <ref> and <ref>
hold, the covariates satisfy the bound (<ref>),
and the exploration conditions (<ref>)
and (<ref>) both hold. Then given any consistent
estimator $\sigmahat^2$ of the error variance $\errorstd^2$, the
estimator $\thetaDecorr$ with tuning
parameters (<ref>) satisfies

\begin{align}
\label{eqn:gen-case-asymp-normality}
\big(\log\numobs \cdot \log\log(\numobs) \cdot
\sigmahat^2\big)^{-1/2} \cdot \SigmaMat_\numobs^{\frac{1}{2}} \left(
\thetaDecorr - \thetastar \right) \indistrb \Ncal(0, \MyIdDim).
\end{align}

Moreover, for any $\dir \in \real^\dims$, the following is an
asymptotically exact $1 - \alpha$ confidence intervals for $ \dir^\top
\thetastar$

\begin{align}
\label{eqn:gen-case-simple-CI}
\left[ \e_1^\top \thetaBasis - \tfrac{\OpTuning_\numobs
\sigmahat}{\sqrt{\tuneParScaled_\numobs}}
(\myinprod{\dir}{\SigmaMat_\numobs^{-1} \dir})^{\frac{1}{2}}
\NormalQuantile{1 - \alpha/2}, \qquad \e_1^\top \thetaBasis
\tfrac{\OpTuning_\numobs\sigmahat}{\sqrt{\tuneParScaled_\numobs}}
(\myinprod{\dir} {\SigmaMat_\numobs^{-1}\dir})^{\frac{1}{2}}
\NormalQuantile{1 - \alpha/2} \right],
\end{align}

where $\OpTuning_\numobs = \opnorm{\diagCov_{\dir,
\numobs}^{-\frac{1}{2}} \SigmaMat_{\dir, \numobs}^{\frac{1}{2}}}$
and the estimator $\thetaBasis$ was calculated
using (<ref>) and with the tuning
parameters (<ref>).

See <ref> for the proof.

It is worth noting that the bounded covariate
condition (<ref>) can be relaxed.  For
instance, in absence of the condition (<ref>),
one may obtain a result similar to the part (a) of
Corollary <ref> under the following assumptions:

\begin{align*}
\tuneParScaled_\numobs = \smalloP(\log \lambda_{\max}
(\SigmaMat_\numobs)), \quad \text{and} \quad \sum_{i = 1}^\numobs
\greedyEps_i = \frac{ \max_{i \in [\numobs]} \Exs[\| \x_i \|^2_2]
}{\lambda_{\min}(\ExsXsqLB) \cdot \smalloP(\tuneParScaled_\numobs)}.
\end{align*}

Finally, as a special case, Corollary <ref>
allows us to construct confidence interval for $\dir^\top \thetastar$
for multi-armed bandit problems that we discussed in
Section <ref>.  The
condition (<ref>) is readily satisfied for
multi-armed bandit problems, but the
conditions (<ref>)
and (<ref>) are mildly stronger than the
analogous condition (<ref>).

§.§.§ Ensuring sufficient exploration via data augmentation

It is natural to ask if we can obtain online debiased method when the
sufficient exploration condition is either difficult to verify or is
not satisfied. In such settings, one simple fix is the following is
based on the data-augmentation technique disscued in
Section <ref>. Indeed, if we may collect
$\log^2(n)$ many data points with $\x_i$ chosen uniformly at random
froma $d$-dimensional unit sphere and append it to the new data. The
new data-set has $n + \log^2(n)$ many data points, and the new dataset
satisfy the exploration
condition (<ref>) with

\begin{align*}
\greedyEps_i = 0 \;\; \text{for all} \;\; 1 \leq i \leq \numobs \qquad
\greedyEps_i = 1 \;\; \text{for} \;\; i = \numobs + 1, \ldots, \numobs
+ \log^2(\numobs).
\end{align*}

In summary, the growth condiiton (<ref>) is
satisfied in this case up to a factor $d$.  Hence, the result from
Corollary <ref> holds true in this case.

§.§.§ Numerical simulation

Average coverage and width of confidence intervals for
$\theta_1^*$ and $\theta_2^*$ across 1000 independent
replications of a linear bandits
experiment (<ref>) with $\thetastar \equiv
(\theta_1^*, \theta_2^*) = (0.3, 0.3)^\top$. The covariates
$\{\x_i\}_{i = 1}^{1000}$ were selected using the
$\greedyEps$-greedy linear bandits
algorithm (<ref>), and the error bars
represent $\pm 1$ standard error.  Left and
Center: Coverage of one-sided $1 - \alpha$ intervals
for $\theta_1^*$ and $\theta_2^*$. Right: Width of
two-sided $1 - \alpha$ intervals for $\theta_1^*$ and
$\theta_2^*$.  See <ref> for details.

<ref> illustrates the performance of online
debiasing with the active learning tuning (<ref>).
Here we consider a linear bandits problem with $\thetastar = (0.3,
0.3)^\top$ and i.i.d. standard normal error $\{ \error_i \}_{i =
1}^\numobs$.  The covariates $\{ \x_i \}_{i = 1}^\numobs$ were
generated using the $\greedyEps$-greedy linear bandits
algorithm (<ref>), where, for each stage, the
context set $\actionset_i$ consisted of the same $50$ vectors drawn
and uniformly from the unit sphere in dimension $2$.  For this
problem, the exploration lower bound
<ref> is satisfied with $\mathbf{G} =
\frac{1}{|\actionset|} \cdot \sum_{\action_i \in \actionset} \action_i
\action_i^\top$.  In this setting, 2 only
provide concentration-based CIs based on ridge regression estimators,
rather than .  Here we report the CIs from ridge regression with
regularization parameter $\lambda_{\text{Ridge}} = 0.1$ (which closely
approximates the solution) and display analogous results for
alternative regularization parameters in
Appendix <ref>. We computed the
confidence intervals for $\theta_1^*$ and $\theta_2^*$ using
Corollary <ref>.

We observe once more that online debiasing provides appropriate
coverage for all confidence levels, while the lower tail interval
consistently undercovers. Meanwhile, the concentration CI provides high
coverage for all confidence levels but yields intervals typically
larger than the online debiasing CIs.

§ PROOFS OF THE THEOREMS

In this section, we provide the proofs of our two main results. We
prove Theorem <ref> in
Section <ref>, and
Theorem <ref> in
Section <ref>.

§.§ Proof of Theorem <ref>

Using the condition $\lambda_{\min}(\SigmaMat_\numobs) \almostSurely \infty$
from assumption <ref>, thus we may assume without loss of generality that $\SigmaMat_\numobs$
is invertible. We claim that it suffices to show that
converges in distribution to $\normal(0,
\sigma^2 \MyIdDim)$.  Indeed, when this claim holds, then since
by assumption, Slutsky's theorem
implies the claim of the theorem.

Recall from equation (<ref>) that the random
vector $\sqrt{\tuneParScaled_\numobs} \cdot
\SigmaMat_\numobs^{\frac{1}{2}} (\thetaDecorr - \thetastar)$ can be
decomposed into the sum $\Bias_\numobs + \ZeroMartin_\numobs$. Based
on this decomposition, we see that it is sufficient to prove that
and   The remainder of our proof is devoted to
establishing these two claims.

Analysis of $\Bias_\numobs$

By definition of the operator norm, we have the upper bound

\begin{align}
\label{EqnFirstTriangle}
\| \Bias_\numobs \|_2 & \leq \sqrt{\tuneParScaled_\numobs} \, \opnorm{
\MyIdDim - \W_\numobs \X_\numobs \SigmaMat_\numobs^{-\frac{1}{2}} }
\big \| \SigmaMat_\numobs^{\frac{1}{2}}(\thetaLS - \thetastar)
\big\|_2.
\end{align}

Lemma 1 from the paper [18] guarantees that

\begin{align}
\label{eqn:ls-error}
\big\| \SigmaMat_\numobs^{\frac{1}{2}}(\thetaLS - \thetastar)
\big\|_2 & = \Bigo \left(\sqrt{ \log
\lambda_{\max}(\SigmaMat_\numobs) } \right) \quad \text{almost
\end{align}

On the other hand, the vanishing bias condition from
Assumption <ref>(b) guarantees that

\begin{align}
\label{eqn:biasterm}
\sqrt{\tuneParScaled_\numobs \log \lambda_{\max}(\SigmaMat_\numobs) }
\cdot \opnorm{ \MyIdDim - \W_\numobs \X_\numobs
\SigmaMat_\numobs^{-\frac{1}{2}} } \stackrel{\smallop{p}}{\longrightarrow} 0.
\end{align}

Applying the bounds (<ref>) and (<ref>) to
the right-hand side of the inequality (<ref>)
shows that $\|\Bias_\numobs\|_2 \stackrel{\smallop{p}}{\longrightarrow} 0$.

Analysis of $\ZeroMartin_\numobs$

In order to control the second term, we seek to apply a classical
martingale central limit theorem (cf. Theorem 2.2 in the
paper [8]).  We begin by observing that
$\NseqNosub{\sqrt{\tuneParScaled_\numobs} \w_i \err_i}{i}$ is a
martingale difference sequence with respect to the sigma-field
$\Nseq{\filtration}{i}$. Noting that the tuning parameter
$\tuneParScaled_\numobs$ is non-random, it follows that the sum
$\Nsum{i}\sqrt{\tuneParScaled_\numobs} \w_i \err_i$ has zero mean and
moreover that

\begin{align}
\Nsum{i} \Cov \left[\sqrt{\tuneParScaled_\numobs} \w_i \err_i \mid
\filtration_{i - 1}\right] &= \tuneParScaled_\numobs \Nsum{i} \w_i
\w_i^\top.
\end{align}

Consequently, in order to apply the martingale CLT so as to obtain the
stated claim, we need to show that

\begin{align*} \tuneParScaled_\numobs \Nsum{i} \w_i
\w_i^\top \stackrel{\smallop{p}}{\longrightarrow} \noisesd^2 \Id_{\Dim}.
\end{align*}

Doing so requires the following auxiliary lemma, which characterizes
the behavior of the weight vector sequence $\Nseq{\w}{i}$ constructed
in equation (<ref>).

Under the Assumption <ref> parts (a) and (c), the
sequence of vectors $\Nseq{\w}{i}$ obtained from
equation (<ref>) has the following properties:

\begin{align*}
\text{(Stability:)} &\qquad \tuneParScaled_\numobs \sum_{i = 1
}^{\numobs} \w_i \w_i^\top \stackrel{\smallop{p}}{\longrightarrow} \Id_p, \quad \text{and}
\\ \text{(Vanishing norm:)} &\qquad \max_{i \in [\numobs]}
\sqrt{\tuneParScaled_\numobs} \|\w_i\|_2 \stackrel{\smallop{p}}{\longrightarrow} 0.
\end{align*}

See Appendix <ref> for the proof of
this lemma.

With the above lemma in hand, we now apply a standard martingale
central limit theorem[Concretely, by applying Theorem 2.2
from the paper [8], we first show that
for any unit vector $u$, the inner product $\tfrac{1}{\noisesd}
\myinprod{u}{ \Nsum{i} \sqrt{\tuneParScaled_\numobs} \w_i \err_i}$
converges to a standard Gaussian.] to conclude that

\begin{align*}
\sum_{i = 1}^{\numobs} \sqrt{\tuneParScaled_\numobs} \w_i \err_i
\indistrb \Ncal(0, \noisesd^2 \Id_\Dim).
\end{align*}

Putting together the pieces, we conclude that

\begin{align*}
\sqrt{\tuneParScaled_\numobs} \cdot \SigmaMat_\numobs^{\frac{1}{2}}
(\thetaDecorr - \theta) = \Bias_\numobs + \sum_{i = 1}^{\numobs}
\sqrt{\tuneParScaled_\numobs} \w_i \err_i \indistrb \Ncal(0,
\noisesd^2 \Id_\Dim),
\end{align*}

which completes the proof of Theorem <ref>.

§.§.§ Proof of claim (<ref>):

For simplicity, let us assume $\noisesd$ is known. Recalling the decomposition (<ref>)
we have

\begin{align*}
\tuneParScaled_\numobs \cdot
\| \thetaDecorr - \thetastar \|_{\SigmaMat_\numobs}^2
&= \enorm{\Bias_\numobs}^2 +
2 \myinprod{\Bias_\numobs}{\ZeroMartin_\numobs} + \enorm{\ZeroMartin_\numobs}^2 \\
&\leq \enorm{\Bias_\numobs}^2
+ 2 \enorm{\Bias_\numobs} \cdot \enorm{\ZeroMartin_\numobs}
+ \enorm{\ZeroMartin_\numobs}^2
\end{align*}

Invoking the condition $\sqrt{\tuneParScaled_\numobs \log
\lambda_{\max}(\SigmaMat_\numobs) } \cdot \opnorm{\MyIdDim -
\W_\numobs \X_\numobs \SigmaMat_\numobs^{-\frac{1}{2}}}  \inEllOne
0$ we immediately have $\enorm{\Bias_\numobs}^2
\inEllOne 0$.
It suffices to show that $\enorm{\ZeroMartin_\numobs}^2 \leq \dims$
and $\Exs[\enorm{\ZeroMartin_\numobs}^2] \rightarrow \dims$. Observe that

\begin{align*}
\Exs[\enorm{\ZeroMartin_\numobs}^2]
&= \tuneParScaled_\numobs \cdot \sum_{i = 1}^{\numobs} \Exs [\error_i^2 \enorm{\w_i}^2]
+ \sum_{i \neq j} \tuneParScaled_\numobs \cdot \Exs [\error_i \error_j \w_i^\top \w_j] \\
& \stackrel{(i)}{=} \sum_{i = 1}^{\numobs} \tuneParScaled_\numobs \cdot \noisesd^2 \cdot \Exs [\enorm{\w_i}^2]
+  \sum_{i \neq j} \tuneParScaled_\numobs \cdot \Exs [\error_i \error_j \w_i^\top \w_j]
\end{align*}

The last line above follows from the assumption $\Exs[\error_i^2 \mid
\filtration_{i - 1}] = \noisesd^2$ and the fact (by construction)
that $\w_i \in \filtration_{i - 1}$. Taking trace on both sides of
equation (<ref>) we have that
, and
using the stronger $L_1$ version of assumption <ref>
along with the proof techniques of Lemma <ref> we
have $ \tuneParScaled_\numobs \cdot \sum_{i = 1}^\numobs
\Exs[\enorm{\w_i}^2] \rightarrow \noisesd^2 \dims$. It remains to show
that $\Exs [\error_i \error_j \w_i^\top \w_j] = 0$ for all $i \neq
j$. Without loss of generality, assume $i < j$. By construction of
$\w_i$ and the martingale assumption <ref> of the noise
$\error_i$, we have $\{ \w_i, \w_j, \error_i\} \in \filtration_{j -
1}$. As a result, we conclude that

\begin{align*}
\Exs [\error_i \error_j \w_i^\top \w_j] = \Exs \Big[ \error_i \cdot
\w_i^\top \w_j \Exs[ \error_j \mid \filtration_{j - 1}]\Big] = 0,
\end{align*}

thereby completing the proof of the claim (<ref>).

§.§ Proof of Theorem <ref>

We prove part (a) of Theorem <ref> in
Section <ref> and part (b) of
Theorem <ref> in
Section <ref>.

§.§.§ Proof of Theorem <ref>(a)

Throughout the proof, we use $\thetahat$ to denote a generic estimator
for $\thetastar$.  We assume that the estimator $\thetahat$ is a
function only of the $\numobs$ datapoints $\{(\x_i, \y_i) \}_{i =
1}^{\numobs}$ and family of selection algorithms
one for each $i \in [\numobs]$; of course, the estimator $\thetahat$
does not know the value of the true parameter $\thetastar$.  Consider
any positive semidefinite and potentially data-dependent matrix
$\mahal \in \reals^{d\times d}$, and define the nonnegative scalar
loss function

\begin{align}
\label{eqn:Loss}
\Loss(\thetahat, \thetastar) & \mydefn (\thetahat - \thetastar)^\top
\mahal (\thetahat - \thetastar).
\end{align}

§.§ From minimax to Bayes risk

In terms of the above notations, Theorem <ref>
(a) posits a lower bound on the minimax risk:

\begin{align}
\label{eqn:minimax-risk}
\inf_{\thetahat} \sup_{\thetastar \in \real^\dims} \;
\Exs[\Loss (\thetahat, \thetastar) \mid \thetastar],
\end{align}

where in the above expression, we have taken an expectation of the
loss $\Loss(\thetahat, \thetastar)$ over the randomness in the data
$(\X_\numobs, \Y_\numobs)$ conditioned on $\thetastar$. We establish
the lower bound Theorem <ref>(a) on the minimax
risk (<ref>) via the standard avenue of first lower
bounding the minimax risk by the Bayes risk, and then providing a
lower bound on the Bayes risk.  In order to do so, we make use of the

\begin{align}
\label{eqn:Minimax-to-Bayes}
\inf_{\thetahat} \sup_{\thetastar} \;
\Exs[
\Loss(\thetahat, \thetastar) \mid \thetastar]
\geq \inf_{\thetahat} \; \Exs_{\thetastar, \X_\numobs, \Y_\numobs}
\Loss(\thetahat, \thetastar),
\end{align}

where the expectation $\Exs_{\thetastar, \X_\numobs, \Y_\numobs}$
above is taken with respect the joint distribution on $(\thetastar,
\X_\numobs, \Y_\numobs)$. Note that this joint distribution is defined
by choosing a prior distribution over the parameter
$\thetastar$, and this choice of prior is a design parameter in our

§.§ Main argument

We claim that it suffices to prove that for any

\begin{align}
\label{eqn:bayes-risk-bound}
\Exs[\Loss(\thetahat, \thetastar) \mid \X_\numobs, \Y_\numobs] \geq
\sigma^2 \tr(\mahal \SigmaMat_\numobs^{-1}).
\end{align}

where the expectation $\Exs[\cdot \mid \X_n, \y_n]$ is taken with
respect to the conditional distribution of  Indeed, taking expectation over $\X_n, \y_n$ yields
the desired bound:

\begin{align}
\Exs_{\thetastar, \X_\numobs, \Y_\numobs} \Loss(\thetahat, \thetastar)
& = \Exs_{\X_\numobs, \Y_\numobs} \Exs[\Loss(\thetahat, \thetastar)
\mid \X_\numobs, \Y_\numobs] \geq \sigma^2 \Exs_{\X_\numobs,
\Y_\numobs} \tr(\mahal \SigmaMat_\numobs^{-1}).
\end{align}

Accordingly, it remains to prove the
bound (<ref>).

§.§.§ Proof of bound (<ref>)

We complete the proof of this bound by first computing the conditional
distribution of $\thetastar \mid \X_\numobs, \Y_\numobs $, and then
lower bounding the conditional expectation of the loss
$\Loss(\thetahat, \thetastar)$ given the data $(\X_n, \y_n)$.
Concretely, we show that under the prior distribution $\thetastar \sim
\Ncal(0, \rho^2 \Id_\Dim)$, we have

\begin{align}
\label{eqn:Conditional-dist-of-thetastar}
\thetastar \mid \X_\numobs, \Y_\numobs & \sim \Ncal( \mu_n,
\PopSigma_n), \quad \text{where} \notag \\ \mu_n = \PopSigma_n
\X_\numobs \Y_\numobs \quad &\text{and} \quad \PopSigma_n =
(\SigmaMat_\numobs/\sigma^2 + \Id_\Dim/\rho^2)^{-1}.
\end{align}

A simple calculation using these distributional properties yields that
for any  matrix $\mahal$—one that may
depend on the data $(\X_n, \y_n)$—the function  is
minimized[Here we have assumed that the prior distribution
$\thetastar \sim \Ncal(0, \rho^2)$ and the error variance $\sigma^2$
are known to the estimator $\thetahat$; this assumption is justified
since without the knowledge of the prior distribution on $\thetastar$
and error-variance $\sigma^2$, the minimum value of the expected loss
$\Exs[\Loss(\thetahat, \thetastar) \mid \X_\numobs, \Y_\numobs]$ can
only increase, which yields a (possibly) stronger lower bound.]  by
the choice .  Moreover, this choice of estimator yields the minimum
value $\sigma^2\tr(\mahal \PopSigma_n)$.  Finally, we are free to
choose the value of the prior error variance $\rho^2$; in particular,
taking the limit $\rho^2 \rightarrow \infty$ yields the
claim (<ref>).

It remains to prove the auxiliary
claim (<ref>).

§.§ Proof of
claim (<ref>)
We proceed via
induction on the number of datapoints $\numobs$.

*Base case For $n = 0$, we have

\begin{align}
\thetastar \mid \X_0, \Y_0 \equiv \thetastar \sim \Ncal(0, \rho^2
\Id_\Dim),
\end{align}

using the facts that $\thetastar \sim \Ncal(0, \rho^2 \Id_\Dim)$ by
our choice of prior, and the triple $(\X_0, \Y_0, \SigmaMat_0)$ are
defined as zeros of respective dimensions. This proves the
statement (<ref>) for $n = 0$, and
$\mu_0 = 0$, and $\PopSigma_0 = \rho^2 \Id_\Dim$.

*Induction step
Given some  assume that the
claim (<ref>) holds for $\numobs -
1$.  Here we show that the statement then holds for $n$.  Recall that
the query algorithm  is oblivious to the true
value $\thetastar$; thus, the conditional distribution $\x_n \mid
\filtration_n$ is independent of $\thetastar$ (see the discussion
before Theorem <ref>). Furthermore, from the
model (<ref>), it follows that the conditioned
random variable $y_n \mid \filtration_{n - 1}, \x_n, \thetastar$
follows a $\Ncal(\x_\numobs^\top \thetastar, \sigma^2)$ distribution,
and using the induction
hypothesis (<ref>), we conclude that

Now let $\frac{d \Prob (\thetastar \mid \X_n, \Y_n)}{d
\lebesgue^{d}(\thetastar)}$ denote the Radon-Nikodym derivative of
the conditional distribution defind by $\thetastar \mid \X_n, \y_n$
with respect to the Lebesgue measure $\lebesgue^{d}$ on $\real^\dims$.
With the last three observations in hand, an application of Bayes'
rule yields

\begin{align*}
\frac{d \Prob (\thetastar \mid \X_n, \Y_n)}{d
\lebesgue^{d}(\thetastar)} & \propto \exp \left\lbrace
-\frac{1}{2}(\thetastar - \mu_{n - 1})^\top \PopSigma_{n - 1}^{-1}
(\thetastar - \mu_{n - 1})^\top \right\rbrace \times \exp
\left\lbrace -\frac{1}{2\sigma^2}(\y_n - \x_n^\top \thetastar)^2
\right\rbrace \\ &\propto \exp \left\lbrace -\frac{1}{2}(\thetastar -
\mu_{n})^\top \PopSigma_{n} (\thetastar - \mu_{n})^\top
\right\rbrace,
\end{align*}

where the pair $(\mu_n, \PopSigma_n)$ are given by

\begin{align*}
\PopSigma_{n}^{-1} = \PopSigma_{n - 1}^{-1} + \frac{\x_n
\x_n^\top}{\sigma^2} \quad \text{and} \quad \mu_n =
\frac{1}{\sigma^2} \PopSigma_{n} \sum_{i = 1}^{n} \x_i \y_i.
\end{align*}

This completes the proof of the inductive step, and putting together
the pieces yields the claim of part (a)
of <Ref>.

§.§.§ Proof of Theorem <ref>(b)

This proof and construction follows by discretizing the construction in [19], which establishes a similar result in a kernelized version of the problem in continuous time. Regrettably, the error terms that arise as a consequence of the discretization lead to a rather unpleasant calculation. We may assume without loss of generality that $\sigma = 1$. The
general result with $\sigma > 0$ can be obtained by a rescaling
argument. For simplicity we also assume that $n$ is divisible by $d-1$, which can be relaxed by correctly rounding the indices of the many sums that appear in the calculations that follow. Let $d > 1$ and $\dir = e_d =
(\bm{0}, 1)^\top \in \reals^d$.
Our proof follows a standard Bayesian argument.
Consider a randomly generated vector $\thetastar$ that is equal to
$\zero$ with probability $1/2$, and otherwise sampled from a multivariate
Gaussian distribution with mean $\mu = (\bm{0}, 1)^\top \in \reals^d$ and degenerate

\begin{align*}
\Sigma = \begin{pmatrix} \MyIdDim & \zero \\ \zero^\top & 0 \end{pmatrix}\,.
\end{align*}

Lower bounding the supremum by an expectation over this prior, the
minimax risk can be lower bounded as

\begin{align}
\label{eq:lower-1}
\inf_{\cE} \sup_{\thetastar} \E_{\thetastar}\left[\frac{(\ip{\dir,
\thetastar} - \cE)^2}{\norm{\dir}^2_{\SigmaMat_n^{-1}}}\right] \geq
\inf_{\cE} \E\left[\frac{(\ip{\dir, \thetastar} -
\cE)^2}{\norm{\dir}^2_{\SigmaMat_n^{-1}}}\right] \,.
\end{align}

where the second expectation integrates over randomness in
$\thetastar$ as well as the observations $\{(\x_i, \y_i)\}_{i \in

We now provide a sequential definition of the selection algorithm that
yields the claimed lower bound. Each covariate $\x_i$ is supported on
the last coordinate as well as one of the first $d-1$ coordinates, chosen in round-robin fashion.
Let $u_i = 1 + (i \operatorname{mod} (d-1))$ and $v_i = \ceil{i/(d-1)}$, which are chosen so that

\begin{align*}
\{u_i\}_{i\geq 1} = \{1,2,\ldots,d-1,1,2,\ldots\} \text{ and } \{v_i\}_{i \geq 1} = \{\underbrace{1,1,\ldots,1}_{d-1 \text{ times}},2,2,\ldots\} \,.
\end{align*}

In other words, $u_i$ is the index of the first non-zero coordinate in $\x_i$
and $v_i$ is the number of times coordinate $u_i$ was non-zero in rounds $j \leq i$.
The first $d-1$ coordinates of the covariate process are deterministic and the last coordinate is chosen adaptively to maximize the difficulty of estimation. Precisely, $\{\x_i\}_{i
\geq 1}$ is given by

\begin{align*}
\x_i \mydefn b_{v_i} e_{u_i} + a_{u_i,v_i} e_d \in \reals^d \,,
\end{align*}

$b_v \mydefn v^{-1/4} / \sqrt{d}$, and the random sequence $\{a_{u,v}\}_{u \in \{1,\ldots,d-1\},v \geq 1}$ is to
defined momentarily.
Let $\y_{u,w} \mydefn \y_{u+(w-1)(d-1)}$, which is the observed response in the round $i$ where coordinate $u$ was non-zero for the $v$th time. Define

\begin{align*}
m_{u,v} \mydefn \sum_{w=1}^v b_w (\y_{u,w} - a_{u,w}), \quad \mbox{and} \quad d_v
\mydefn 1 + \sum_{w=1}^v b_w^2.
\end{align*}

Then $a_{u,v} \mydefn -b_v m_{u,v-1} / d_{v-1}$, noting that $a_{u,1} = 0$ for all $u$.

To provide some intuition, our construction is designed so as to make
estimation challenging. Let $\error_{u,v} = \error_{u+(v-1)(d-1)}$. On the event that $\thetastar \neq \zero$, we
have the equality $\y_{u,v} = a_{u,v} + \thetastarScalar_u b_v + \epsilon_{u,v}$, and
the ratio $m_{u,v} / d_v$ is the ridge regression estimate of
$\thetastarScalar_u$. For any vector $\thetastar \neq 0$, the choice $\x_i$ ensures that

\begin{align*}
\ip{\x_i, \thetastar} = a_{u_i,v_i} + b_{v_i} \thetastarScalar_{u_i} \approx 0 =
\ip{\x_i, \zero}\,,
\end{align*}

so that the observed responses are extremely similarly under either of
the events $\{\thetastar = \zero\}$ or $\{\thetastar \neq \zero\}$.
But $\ip{\dir,\thetastar} = \bm{1}(\thetastar \neq \zero)$, which means that any estimator of $\ip{\dir,\thetastar}$ must have large
error in expectation. What is missing is to formalize the above claims
and show that $\norm{\dir}^2_{\SigmaMat^{-1}_n}$ shrinks suitably fast.

Returning the proof, since $\norm{\dir}^2_{\SigmaMat_n^{-1}}$ is
$\filtration_n$-measurable, the infimum on the right-hand side of
equation (<ref>) is achieved by the estimator

\begin{align*}
\cE_\numobs \mydefn \E[\ip{\dir, \thetastar}|\filtration_n] =
\Prob\left(\thetastar \neq \zero|\filtration_n\right).
\end{align*}

Therefore, introducing the event $\Event = \{\thetastar \neq \zero\}$,
we have the lower bound

\begin{align*}
\textrm{Risk} & \geq \E \Biggr[ \frac{\big(\Prob \big[\thetastar \neq
\zero \mid \filtration_\numobs \big] - \ip{e_d,
\thetastar}\big)^2}{\norm{e_d}^2_{\SigmaMat_\numobs^{-1}}}
\Biggr] = \E \Biggr[
\frac{\cE_\numobs(1-\cE_\numobs)}{\norm{e_d}^2_{\SigmaMat_n^{-1}}}
\Biggr] \geq \frac{1}{2} \E \Biggr[
\frac{\cE_\numobs(1-\cE_\numobs)}{\norm{e_d}^2_{\SigmaMat_n^{-1}}}
\Bigg| \Event \Biggr] \,.
\end{align*}

In the remainder of the proof, we study the laws of $\cE_\numobs$ and
$\smash{\norm{e_d}^2_{\SigmaMat_n^{-1}}}$ under the measure
$\Prob(\cdot \mid \Event)$.  For a sequence of random variables $X_n$,
we use the notation $X_\numobs = \Omega_p(a_n)$ to mean that $\sup_n
\Prob(X_n / a_n < C_\epsilon \mid \Event) < \epsilon$. We will show below that

\begin{align}
\label{eq:lower-a}
\cE_\numobs(1-\cE_\numobs) &= \Omega_p(1), \qquad \mbox{and} \\
\label{eq:lower-b}
\frac{1}{\norm{e_d}_{\SigmaMat_n^{-1}}^2}
&= \Omega_p(d \log(n/d^3)) \,.
\end{align}

Therefore, there exists a universal constant $C > 0$ such that for $n \geq d^3/C$,

\begin{align*}
\Prob\left(1/\norm{e_d}_{\SigmaMat_n^{-1}}^2 \geq C d \log(n)\right) \geq 3/4 \text{ and }
\Prob\left(\cE_\numobs(1 - \cE_\numobs) \geq C\right) \geq 3/4\,.
\end{align*}

By a union bound and the positivity in the integrand of the risk,

\begin{align*}
\textrm{Risk} &\geq \frac{C^2 d \log(n/d^3)}{4} \,.
\end{align*}

It remains to show that the lower bounds (<ref>)
and (<ref>) hold, which we prove in Appendix <ref>.

§.§.§ Proof of <Ref>(c)

Let $\ExsData$ denote expectation over a data set drawn from the
distribution indexed by $\thetabf$, and define $\ProbData$ as the
analogous probability.  For a given dataset based on $\numobs$
samples, consider a confidence interval of the form $[\lowhat_\numobs,
\uphat_\numobs]$.  Introducing the shorthand $\smallv^2 \mydefn v^T
\SigmaMat_\numobs^{-1} v$, we then define the minimax risk

\begin{align*}
\CIrisk & \mydefn \inf_{\hackCI} \sup_{\thetastar} \Exs_\thetastar
\bigg\{ \frac{\uphat_\numobs - \lowhat_\numobs}{\smallv} \bigg\},
\end{align*}

where we take the infimum over all estimators $\hackCI$ such that
$\Prob_{\thetastar} \big( [\lowhat_\numobs, \uphat_\numobs] \ni v^\top \thetastar \big)
\geq 1 - \alpha$ for each value of $\thetastar$.

By the usual Bayesian argument, for any prior distribution $\pi$ on $\thetastar$, we have the lower

\begin{align*}
\CIrisk \geq \inf_{\hackCI} \E_{\thetastar \sim \pi} \Exs_\thetastar \bigg[
\frac{\uphat_\numobs - \lowhat_\numobs}{\smallv} \bigg].
\end{align*}

We now obtain a further lower bound by enlarging the space of possible
estimators $\hackCI$, in particular requiring only that $\hackCI$
belong to the set

\begin{align*}
\mathcal{A} & \mydefn \Big \{ \hackCI \; \mid \; \Exs_{\thetastar \sim
\pi} \Prob_\thetastar \big( [\lowhat_\numobs, \uphat_\numobs] \ni
v^\top \thetastar\big) \geq 1 - \alpha \Big \}.
\end{align*}

Since this allows for a larger collection of possible estimators, we
have the lower bound

\begin{align*}
\CIrisk & \ge \inf_{\hackCI \in \mathcal{A}} \E_{\thetastar \sim \pi}
\Exs_\thetastar \bigg\{ \frac{\uphat_\numobs - \lowhat_\numobs}{\smallv}
\bigg\}.
\end{align*}

We are now free to choose the prior.  In particular, we set $\pi$
equal to the density $\phi_\rho$ of the Gaussian random vector
$\Ncal(0, \rho^2 \Id_d)$.  From our previous calculations (<ref>) we have that conditional on the observed
data, the random vector $\thetastar$ is Gaussian with
covariance $(\SigmaMat_n/\sigma^2 +\Id_d/\rho^2)^{-1}$. With this
choice, the random variable $v^\top\thetastar$, conditioned on
the observed data, is a Gaussian random variable with variance $\smallvtil^2
\mydefn v^T (\SigmaMat_n/\sigma^2 + \Id_d/\rho^2)^{-1} v$.  Therefore, the width of
any confidence interval  $\hackCI \in \mathcal{A}$ is lower bounded by
$2 z_{1-\alpha/2} \smallvtil$, and we have

\begin{align*}
\CIrisk & \ge 2 z_{1-\alpha/2} \E_{\thetastar \sim \pi} \Exs_\thetastar \Big[
\frac{\smallvtil}{\smallv} \Big].
\end{align*}

It remains to show that $\lim_{\rho \to \infty} \E_{\thetastar\sim \pi}
\Exs_\thetastar \big[ \frac{\smallvtil}{\smallv} \big] = 1$.  Recalling that
$\phi_\rho$ denotes the Gaussian density with zero mean and covariance
$\rho^2 I_d$, we have

\begin{align*}
\E_{\thetastar \sim \pi} \Exs_\thetastar \Big[ \frac{\smallvtil}{\smallv} -1
\Big] & = \int \Exs_\thetastar \Big[\frac{\smallvtil}{\smallv} - 1 \Big]
\phi_\rho(\thetastar) d \thetastar,
\end{align*}

By the bounded convergence theorem, we have

\begin{align*}
\lim_{\rho \to \infty} \Exs_\thetastar \Big[\frac{\smallvtil}{\smallv} - 1
\Big] & = 0,
\end{align*}

pointwise for each $\thetabf$.  Consequently, the quantity in the
integral converges point-wise to zero, so that applying the bounded
convergence theorem again yields

\begin{align*}
\lim_{\rho \to \infty} \E_{\thetastar \sim \pi} \Exs_\thetastar \Big[
\frac{\smallvtil}{\smallv} -1 \Big] = \lim_{\rho \to \infty} \int
\Exs_\thetastar \Big[ \frac{\smallvtil}{\smallv} - 1 \Big] \phi_\rho(\thetastar)
d \thetastar & = 0,
\end{align*}

which completes the proof of part (c).

§ PROOFS OF COROLLARIES

We now turn to the proofs of our three corollaries, with
Sections <ref>, <ref>,
and <ref> devoted to the proofs of the
Corollaries <ref>, <ref>
and <ref>, respectively.

§.§ Proof of Corollary <ref>

In light of Theorem <ref>, it suffices to verify
Assumptions <ref>–<ref>.  The assumptions
stated in Corollary <ref> ensure that the error
sequence $\Nseq{\error}{i}$ satisfies Assumption <ref>.
The growth conditions in Assumption <ref> are satisfied
due to the minimum arm-pull
assumption (<ref>).  It remains to verify
the three conditions in Assumption <ref>.

Beginning with the asymptotic negligibility condition, we have

\begin{align*}
\max \limits_{i \in [\numobs]} \;\;
\frac{1}{\tuneParScaled_\numobs} \x_i^T \XscaleMat_i^{-1} \x_i \leq
\frac{1}{\tuneParScaled_\numobs} \frac{ \max_{i \in [\numobs]}
\|\x_i^2 \|}{(\log \numobs)^2} & = \frac{\log\log(\numobs)}{(\log
\numobs)} \rightarrow 0.
\end{align*}

The first inequality above uses the bound $\XscaleMat_i^{-1} \preceq
\frac{1}{\log(\numobs)^2} \cdot \Id_\Dim$ (see the
definition (<ref>)); the second equality uses $\|
\x_i \|_2^2 = 1$, and the final step follows by substituting
${\tuneParScaled_\numobs = 1/((\log\numobs) \cdot \log\log(\numobs))}$.

Turning to the vanishing bias condition in <ref>, we invoke
the operator norm bound (<ref>) on the matrix
$\Id_\Dim - \W_\numobs \X_\numobs \SigmaMat_\numobs^{-\frac{1}{2}}$ to
find that

\begin{align*}
\sqrt{\tuneParScaled_\numobs \log\lambda_{\max}(\SigmaMat_\numobs)}
\cdot \opnorm{ \Id_\Dim - \W_\numobs \X_\numobs
\SigmaMat_\numobs^{-\frac{1}{2}} } & \leq
\sqrt{\tuneParScaled_\numobs \log \numobs} \cdot \BigoP(1) \\
& =  \sqrt{\frac{\log\log(n)}{(\log \numobs)}} \cdot \BigoP(1)
\stackrel{\smallop{p}}{\longrightarrow} 0,
\end{align*}

where we have used the bound  in the above derivation.

Finally, we verify the variance stability condition in
Assumption <ref> with the help of the following lemma

[Commutative guarantee]

For any collection of matrices $\;\;$  that commute with
each other, we have

\begin{align*}
\opnorm{ \Id - \Nsum{i} \w_i \x_i^\top \XscaleMat_i^{- \frac{1}{2}} }
\leq \exp \left( - \frac{\lambda_{\min}(\Nsum{i}
\XscaleMat_i^{-\frac{1}{2}} \x_i \x_i^\top
\XscaleMat_i^{-\frac{1}{2}})}{\tuneParScaled_\numobs} \right).
\end{align*}

See the end of this subsection for the proof of this claim.

Let us complete the proof of <ref>
using <ref>.  In the multi-armed bandit setting
of Corollary <ref>, the matrices
$\NseqNosub{\XscaleMat_i^{-\frac{1}{2}} \x_i \x_i^\top
\XscaleMat_i^{-\frac{1}{2}}}{i}$ are all diagonal, and hence they
commute. Thus, invoking the operator norm bound from
Lemma <ref> yields

\begin{align*}
\opnorm{ \Id_\Dim - \Nsum{i} \w_i \x_i^\top
\XscaleMat_i^{-\frac{1}{2}} } & \leq \exp \left( - \frac{
\lambda_{\min}(\Nsum{i} \XscaleMat_i^{-\frac{1}{2}} \x_i
\x_i^\top \XscaleMat_i^{-\frac{1}{2}}) }{\tuneParScaled_\numobs}
\right).
\end{align*}

Recall that in the bandits model (<ref>), the matrices $\SigmaMat_\numobs$ and $\x_i \x_i^\top$ are diagonal.
By construction (<ref>) and the minimum arm-pull
condition (<ref>), the tuning matrix $\XscaleMat_i$ is also diagonal with diagonal entries upper bounded by the corresponding diagonal entries of the (diagonal) matrix $\SigmaMat_\numobs$. Combining these two observations we have that
$\SigmaMat_\numobs^{-\frac{1}{2}} \x_i \x_i^\top
\SigmaMat_\numobs^{-\frac{1}{2}} \preceq \XscaleMat_i^{-\frac{1}{2}}
\x_i \x_i^\top \XscaleMat_i^{-\frac{1}{2}}$.

Consequently, we find that

\begin{align*}
- \lambda_{\min}(\Nsum{i} \XscaleMat_i^{-\frac{1}{2}} \x_i \x_i^\top
\XscaleMat_i^{-\frac{1}{2}}) & \stackrel{(i)}{\leq} -
\lambda_{\min}(\Nsum{i} \SigmaMat_\numobs^{-\frac{1}{2}} \x_i
\x_i^\top \SigmaMat_\numobs^{-\frac{1}{2}})
= -1,
\end{align*}

where the final equality follows from the definition of
$\SigmaMat_\numobs$.  Substituting the value $\tuneParScaled_\numobs =
1 / (\log(\numobs) \cdot \log\log(\numobs))$ yields

\begin{align*}
\opnorm{\Id_\Dim - \Nsum{i} \w_i \x_i^\top
\XscaleMat_i^{-\frac{1}{2}} } \leq \exp(-1/\tuneParScaled_\numobs)
\; \leq \; \frac{1}{\numobs} \rightarrow 0.
\end{align*}

This verifies the variance stability condition from
Assumption <ref>, and applying
Theorem <ref> yields
Corollary <ref>.

The only remaining detail is to prove
Lemma <ref>.

§.§ Proof of Lemma <ref>

For notational convenience, we use the shorthands  and $\Z_i^\top
\mydefn \begin{bmatrix} \z_1 & \cdots & \z_i \end{bmatrix}$, as
previously introduced in Section <ref>. Substituting the
formula for the weight vector $\w_i$ from <ref>,
and performing some algebra yields

\begin{align*}
\left( \Id - \W_\numobs \Z_\numobs\right)^\top \left(\Id -
\W_\numobs \Z_\numobs\right) & = \prod \limits_{j = 1}^{\numobs}
\left( \Id_\Dim - \frac{\z_{\numobs + 1 -j} \z_{\numobs + 1 -
j}^\top}{\tuneParScaled_\numobs + \|\z_{\numobs +1 -j}\|^2}
\right) \prod \limits_{i = 1}^{\numobs} \left( \Id_\Dim - \frac{\z_i
\z_i^\top}{\tuneParScaled_\numobs + \|\z_i\|^2} \right) \\
\begin{multlined}[t]
\exp \left[ \sum_{j=1}^\numobs \log \left( \Id_\Dim -
\frac{\z_{\numobs + 1- j} \z_{\numobs + 1-
j}^\top}{\tuneParScaled_\numobs + \|\z_{\numobs + 1 - j}\|^2}
\right)\right. \\ + \left.\Nsum{i} \log \left( \Id_\Dim -
\frac{\z_{i} \z_{i}^\top}{\tuneParScaled_\numobs + \|\z_{i}\|^2}
\right) \right]
\end{multlined}\\
& \stackrel{(i)}{\preccurlyeq} \exp \left( - \Nsum{j}
\frac{\z_{\numobs + 1 -j} \z_{\numobs + 1
-j}^\top}{\tuneParScaled_\numobs + \|\z_{\numobs + 1 -j}\|^2} - \Nsum{i}
\frac{\z_i \z_i^\top}{\tuneParScaled_\numobs  + \|\z_{i}\|^2} \right) \\
& \preccurlyeq \exp \left( - 2\cdot\Nsum{i} \frac{\z_i
\z_i^\top}{\tuneParScaled_\numobs} \right),
\end{align*}

where step (i) above uses the fact that $\exp(\log(1 - a)) \leq
\exp(-a)$ for any scalar $a < 1$ and that the matrices $\{\z_i
\z_i^\top\}_{i \in [\numobs]}$ commute.  Via an inductive argument, it
can be verified that the entries of the matrix $\Id_\Dim -
\frac{\z_{i} \z_{i}^\top}{\tuneParScaled_\numobs + \|\z_{i}\|^2}$ are
all upper bounded by $1$.  Putting together the pieces, we conclude
that the operator norm satisfies the bound

\begin{align*}
\opnorm{ \Id - \W_\numobs \Z_\numobs } \leq \exp \left( -
\frac{\lambda_{\min}(\Nsum{i} \z_i
\z_i^\top)}{\tuneParScaled_\numobs} \right),
\end{align*}

as claimed.

§.§ Proof of Corollary <ref>

The proof of this claim is similar to that of
Corollary <ref>; in particular, we need to verify
Assumptions <ref>–<ref>.  Recall that the
time series model (<ref>) in
Corollary <ref> is a special case of the
stochastic linear regression model (<ref>) with
$(x_i, y_i) \equiv (y_{i - 1}, y_i)$; thus the covariance term based
on the data $\{(x_i, y_i) \}_{i \in [\numobs]}$ is given by $\Nsum{i}
x_{i}^2 \equiv \Nsum{i} y_{i - 1}^2$. Here we have used the convention
$y_0 = 0$.

The moment condition <ref> is satisfied since the
additive noise $\error_i$ in the autoregressive
model (<ref>) is assumed to have a standard Gaussian
distribution. Before we verify the remaining conditions, it is helpful
to deduce a few bounds regarding the sample covariance term $\Nsum{i}
y_{i - 1}^2$. In particular, we show that for any $\thetastarScalar
\in (-1, 1]$, the sample covariance term satisfies the following

\begin{align}
\label{eqn:AR1-lower-bound}
\Nsum{i} y_{i - 1}^2 &\rightarrow \infty \quad \text{almost surely},
\quad \\
\label{eqn:AR1_var_scaling}
\log ( \Nsum{i} y_{i - 1}^2 ) &= \BigoP(\log \numobs), \quad
\text{and} \quad (\log\numobs)^2  y_{\numobs}^2 =
\BigoP ( \Nsum{i} y_{i - 1}^2 ),
\end{align}

We prove
these bounds at the end of this sub-section, but let us complete the
proof of the Corollary using these bounds.

First, observe that the condition <ref> follows from the
growth condition (<ref>), and the asymptotic
negligibility condition in <ref> is satisfied by noting

\begin{align*}
\max_{i \in [\numobs]} \frac{1}{\tuneParScaled_\numobs} \cdot
\frac{y_{i - 1}^2}{ \max\{ (\log \numobs)^{2} y_{i -
1}^2 , \;\; \sum_{j =1}^{i-1} y_j^2 \} } & \leq \frac{(\log
\numobs) \cdot \log\log(n) }{(\log \numobs)^{2}}
\rightarrow 0.
\end{align*}

Next, in order to verify the vanishing bias condition in
Assumption <ref>, doing a calculation similar to
Proposition <ref> we find that (see the arguments
leading up to
bounds (<ref>)–(<ref>)
and their proofs)

\begin{align*}
\sqrt{\tuneParScaled_\numobs \log (\Nsum{i} y_{i - 1}^2 )} \cdot
\left| 1 - \Nsum{i} \frac{w_i y_{i -1}}{\XscaleMat_n^{\frac{1}{2}}}
\right | & \stackrel{(i)}{=} \BigoP \left( \frac{1}{\sqrt{\log\log(n)}} \right) \cdot \BigoP \left( 1 + \sqrt{
\frac{\XscaleMat_\numobs}{\Nsum{i} y_{i - 1}^2 } } \right) \\
& \stackrel{(ii)}{=} \BigoP \left( \frac{1}{\sqrt{\log\log(n)}} \right)
\cdot \BigoP(1) \\
& \stackrel{\smallop{p}}{\longrightarrow} 0,
\end{align*}

where step (i) follows by invoking the first part of the
bound (<ref>) and step (ii) uses the second part
of <ref>.

Finally, we verify the variance stability condition
in <ref> with the help of
Lemma <ref>, as previously stated and proved in
the proof of Corollary <ref>.  Note that in dimension
$\Dim = 1$, the commutativity condition in
Lemma <ref> holds trivially.  Consequently, we
may apply Lemma <ref> to the one-dimensional
autoregressive model (<ref>) so as to obtain the

\begin{align*}
\left | 1 - \Nsum{i} \frac{w_i y_{i - 1}}{ \XscaleMat_i^\frac{1}{2}}
\right| & \leq \frac{1}{\numobs}.
\end{align*}

See the calculations following the statement of
Lemma <ref> in the proof of
Corollary <ref> for details on this step.

This verifies the variance stability condition from
Assumption <ref>, and applying
Theorem <ref> yields
Corollary <ref>.

The only remaining detail is to prove the

§.§ Proofs of the
bounds (<ref>)–(<ref>)

The proof of the first part of the bound (<ref>)
follows by invoking Theorem 2 part (i) from the
paper [18].  Concretely, in the
paper [18], the authors showed that when
$|\thetastarScalar| \leq 1$, then there is some constant $a > 0$ such
that $y_\numobs = \BigoP(\numobs^a)$. Thus, we have the relation $\Nsum{i} y_{i -
1}^2 = \BigoP(\numobs^{2a + 1})$, and first part of the
bound (<ref>) follows.

We divide the proof of the remaining bounds into two parts, depending
on the value of $\thetastarScalar$.

§.§ Case 1
First, suppose that $\thetastarScalar = 1$.
Recall that in <ref> we argued that

\begin{align*}
\frac{1}{\numobs^2} \Nsum{i} \y_{i - 1}^2 &\indistrb \int_{0}^1
\wiener^2(t) dt.
\end{align*}

In light of the last relation, the growth
condition (<ref>) is immediate.  For the remaining
bounds, note that ; thus we have , and we conclude that

\begin{align*}
\log(\numobs)^{2}  \y_{\numobs}^2 = \BigoP
\left( \Nsum{i} \y_{i - 1}^2 \right),
\end{align*}

as claimed.

§.§ Case 2
Otherwise, we may assume that $|\thetastarScalar|
< 1$, in which case the term $\Nsum{i} \y_{i - 1}^2$ stabilizes [18];
concretely, we have

\begin{align*}
\frac{1}{\numobs} \Nsum{i} \y_{i - 1}^2 \almostSurely c, \quad
\text{where} \quad c>0 \;\; \text{is a non-random scalar}.
\end{align*}

The growth condition (<ref>) follows directly from
the above relation.  Moreover, we have $y_{\numobs - 1} = \Nsum{i}
{\thetastarScalar}^{i} \error_{\numobs - i} \sim \Ncal \left(0,
\frac{1}{1 - {\thetastarScalar}^2} \right)$.  Putting these two pieces
together yields ${(\log\numobs)^{2} \cdot y_{\numobs - 1}^2 = \BigoP
\left( \Nsum{i} \y_{i - 1}^2 \right)}$.

§.§ Proof of Corollary <ref>

We obtain the first claim of the Corollary <ref>
by applying Theorem <ref>, and the second part of
the Corollary <ref> follows
from <ref>.  We prove these two parts separately.

Proof of claim (<ref>):
In order to apply Theorem <ref> to the setup of
Corollary <ref> it suffices to verify the
Assumption <ref>. Recall that our choice of scaling
$\XscaleMat_i = \Nsum{j} \greedyEps_j \ExsXsqLB$ matrices does not
actually vary as a function of the round $i$.  For this reason, we simply
write $\XscaleMat$ from here onwards.

We begin by verifying the asymptotic negligibility condition
in <ref>.  Observe that

\begin{align}
\label{eqn:asymp-negligibility-suff-explore}
\Exs \left\lbrace \max_{i \in [\numobs]} \;\;
\frac{1}{\tuneParScaled_\numobs} \x_i^\top \XscaleMat^{-1} \x_i
\right\rbrace \leq \frac{1}{\tuneParScaled_\numobs} \cdot \frac{ \Exs
\left [ \max \limits_{i \in [\numobs]} \| \x_i \|_2^2 \right ]}{
\lambda_{\min}(\ExsXsqLB) \; (\Nsum{i} \greedyEps_i ) } \rightarrow
\end{align}

where the first inequality above follows by substituting the value of
the scaling matrix $\XscaleMat$, and the second step follows by
invoking the sufficient exploration
condition (<ref>).

Next, we verify the variance stability and vanishing bias conditions
in <ref>.  In doing, we make use of the following
auxiliary result:

Under the sufficient exploration
condition (<ref>), for any tuning parameter
$\tuneParScaled_\numobs \in (0, 1/ (\log(\numobs) \cdot \log\log(\numobs))]$ and
a sufficient large sample size $\numobs$, we have

\begin{align*}
\Exs \big[ \frobnorm{ \MyIdDim - \W_\numobs \X_\numobs
\XscaleMat^{-\frac{1}{2}} }^2 \big] & \leq \frac{\Dim}{K\numobs}.
\end{align*}

See Section <ref> for the proof of this

Taking <ref> as given, we now complete the
proof of Corollary <ref>.  Note that the
variance stability condition in <ref> follows directly
from the Frobenius norm bound in Lemma <ref> and
by letting the number of datapoints $n \rightarrow \infty$, keeping the
dimension $\Dim$ fixed.

In order to prove the vanishing bias condition in <ref>,
we first bound the operator norm of the matrix $\MyIdDim - \W_\numobs
\X_\numobs \SigmaMat_\numobs^{-\frac{1}{2}}$:

\begin{align}
\opnorm{ \MyIdDim - \W_\numobs \X_\numobs
\SigmaMat_\numobs^{-\frac{1}{2}}} & \leq 1 + \opnorm{ \W_\numobs
\X_\numobs \XscaleMat^{-\frac{1}{2}} } \cdot \opnorm{
\XscaleMat^{\frac{1}{2}} \SigmaMat_\numobs^{-\frac{1}{2}} } \notag \\
& = \BigoP(1)
\label{eqn:WNS-gen-bound},
\end{align}

where the derivation above uses the Frobenius norm upper bound from
Lemma <ref> and the fact that $\opnorm{
\XscaleMat^\frac{1}{2} \SigmaMat_\numobs^{-\frac{1}{2}}} = \BigoP(1)$
by the choice of the tuning parameter $\XscaleMat$; see the
bound (<ref>) for instance. Using the last
bound on $\opnorm{ \MyIdDim - \W_\numobs \X_\numobs
\SigmaMat_\numobs^{-\frac{1}{2}}}$, we then find that

\begin{align*}
\sqrt{\tuneParScaled_\numobs \log \lambda_{\max}(\SigmaMat_\numobs)}
\cdot \opnorm{ \MyIdDim - \W_\numobs \X_\numobs
\SigmaMat_\numobs^{-\frac{1}{2}} } \leq
\sqrt{\tuneParScaled_\numobs \log \lambda_{\max}(\SigmaMat_\numobs)}
\cdot \BigoP(1) \stackrel{\smallop{p}}{\longrightarrow} 0,
\end{align*}

where the last step above utilizes the choice $\tuneParScaled_\numobs
= \smalloP(\log K\numobs)$ and the bound $
\lambda_{\max}(\SigmaMat_\numobs) = \BigoP(\log K\numobs )$. (Recall
the uniform boundedness assumption (<ref>).)
All together, we have verified the assumptions of
Theorem <ref>, so that
Corollary <ref> follows.

It remains to prove Lemma <ref>.

§.§.§ Proof of Lemma <ref>

Throughout this proof, we use the shorthands , $\Z_i^\top = [\z_1, \z_2, \ldots \z_i]$,
$\W_i = [\w_1, \ldots, \w_i]$ and ${\DelMat_i \mydefn \MyIdDim - \W_i
\Z_i}$.  Substituting the expression (<ref>) for
the weight vector $\w_i$ we find that

\begin{align}
\frobnorm{ \DelMat_{i-1} }^{2} - \frobnorm{ \DelMat_{i} }^{2} & =
\frac{\tuneParScaled_\numobs + \| \z_{i}\|_{2}^{2}}{\left(
\tuneParScaled_\numobs/2 + \| \z_{i} \|_{2}^{2}\right)^{2}} \tr
\left\{ \DelMat_{i-1} \z_{i} \z_{i}^{\top}
\DelMat_{i-1}^{\top}\right\} \nonumber \\
\label{eqn:delta_i_decay}
& \geq \frac{1}{ \tfrac{\tuneParScaled_\numobs}{2} + \|\z_{i}
\|_{2}^{2}} \tr \left\{\DelMat_{i-1} \z_{i} \z_{i}^{\top}
\DelMat_{i-1}^{\top}\right\}.
\end{align}

In equation (<ref>), we proved
that the random variable  converges to zero in
probability; consequently, we may assume that

\begin{align*}
\Prob \big[
\max_{i \in \numobs} \| \z_i\|_2^2 \leq \tuneParScaled_\numobs/2
\big] \geq \frac{1}{2}
\end{align*}

for all sufficiently large values of the sample size
$\numobs$. Keeping this in mind, taking expectations conditional on
the sigma-field $\filtration_{i-1}$ on both sides in the
inequality (<ref>), and using the fact that
$\DelMat_{i} \in \filtration_{i - 1}$, we have

\begin{align*}
\Exs \left [ \frobnorm{\DelMat_{i-1} }^{2} \mid \filtration_{i-1}
\right ] - \Exs \left [ \frobnorm{\DelMat_{i}}^{2} \mid
\filtration_{i-1} \right ] \geq \frac{\greedyEps_i}{2
\tuneParScaled_\numobs \Nsum{i} \greedyEps_i} \Exs \left [
\frobnorm{\DelMat_{i - 1}}^{2} \mid \filtration_{i-1}\right].
\end{align*}

Rearranging the last inequality and using the upper bound $(1 - t)
\leq \exp(-t)$ for $t \geq 0$ we obtain

\begin{align*}
\Exs \left [ \frobnorm{\DelMat_{i}}^{2} \mid \filtration_{i-1}
\right] & \leq \exp \left(\frac{-\greedyEps_i}{2
\tuneParScaled_\numobs \Nsum{i} \greedyEps_i }\right)
\frobnorm{\DelMat_{i-1}}^{2}.
\end{align*}

Iterating the last bound $\numobs$ times and removing the conditioning
on the sigma filed $\filtration_{i -1}$, we find that

\begin{align*}
\Exs \left\{ \frobnorm{\DelMat_{\numobs}}^{2} \right\} \leq \Dim
\exp \left(-\frac{1}{2\tuneParScaled_\numobs} \right) \stackrel{(i)}{\leq}
\Dim \exp \left(-\frac{1}{2} \cdot \log(\numobs) \cdot \log\log(\numobs)
\right) \stackrel{(ii)}{\leq} \frac{\Dim}{K\numobs}.
\end{align*}

Here step (i) follows by using $\tuneParScaled_\numobs \leq
\tfrac{1}{\log(\numobs) \cdot \log\log(\numobs)}$, and step (ii) holds since
the sample size is assumed $\numobs$ to be sufficiently large. This
completes the proof of the claim (<ref>).

Proof of
claim (<ref>):  In order to
apply <ref>, it suffices to verify condition <ref> for $\{\Basis \x_i \}_{i = 1}^\numobs$
with the
choice of tuning parameters (<ref>). Note that in the proof
of (<ref>), we already verified that conditions (<ref>), (<ref>), and (<ref>)
ensures that assumption <ref> is satisfied. Fortunately, these three conditions
are not affected by the change of basis transformation, and are readily satisfied by the regressors $\{ \Basis \x_i \}_{i = 1}^\numobs$.

Indeed, for any orthonormal basis matrix $\Basis$, via linearity of expectation, we have

\begin{align*}
\Exs[\Basis \randDir_i \randDir_i^\top \Basis] = \Basis \Exs[ \randDir_i \randDir_i^\top] \Basis^\top
\succeq \Basis \ExsXsqLB \Basis^\top.
\end{align*}

Moreover, for any orthonormal basis matrix $\Basis$, we have

\begin{align*}
\lambda_{\min}(\Basis \ExsXsqLB \Basis^\top) = \lambda_{\min}
(\ExsXsqLB) \qquad \text{and} \qquad \enorm{\Basis \x_i} =
\enorm{\x_i} \leq K.
\end{align*}

Thus, following a proof similar
to (<ref>) we have that the
assumption <ref> parts (a), (c), and
condition <ref> part (b), modified for $\{ \Basis \x_i
\}_{i = 1}^{\numobs}$, are satisfied. Finally, from the bounded
covariates condition (<ref>) we have
$\lambda_{\max}(\SigmaMat_{\dir, \numobs}) \leq K \numobs$, and as a
result, . Combining
this observation with a calculation similar to
equation (<ref>), we deduce that
condition <ref>(b) holds, thereby completing the
proof of the claim (<ref>).

§ DISCUSSION

In this paper, we proposed a family of online debiasing estimators for
adaptive linear regression and analyze their asymptotic properties. We
introduced an online debiasing estimator, and proved that it admits a
Gaussian limit under considerably weaker conditions than the estimator. We highlighted its practical behavior using examples from
multi-armed bandits, time series modeling, and active learning in
which online debiasing yields asymptotic normality while does
not.  We also proved a minimax lower bound for the adaptive linear
regression model; in conjunction with our upper bounds, our results
reveal that the online debiasing estimator is minimax optimal.

This work opens up a number of directions for future research.  For
example, it would be interesting to characterize the non-asymptotic
behavior of estimators based on online debiasing.  Concretely, we
would like to investigate the rate of distributional convergence of
the online debiasing estimators to the appropriate Gaussian

§.§ Acknowledgements
This work was partially supported by
a BAIR-Microsoft research grant to MJW and LM, as well as DOD ONR
Office of Naval Research N00014-21-1-2842, National Science Foundation
DMS grant 2015454, and National Science Foundation CCF grant 1955450 to MJW.

§ PROOFS OF THE PROPOSITIONS

This section provides the proofs of the two propositions stated in
this paper.   <Ref> is devoted to the proof
of <Ref>, whereas <Ref> is
devoted to the proof of <Ref>.

§.§ Proof of <Ref>

Our proof is based on the following auxiliary result that
characterizes the asymptotic behavior of $\thetaBasis$.

Under , given any consistent estimator
$\widehat{\noisesd}^2$ of $\noisesd^2$, we have

\begin{align}
\label{eqn:confinv-asn}
\sqrt{\frac{\tuneParScaled_\numobs}{ \OpTuning_\numobs^2
\widehat{\noisesd}^2}} \cdot
\diagCov_\numobs^{\frac{1}{2}}(\thetaBasis - \Basis \thetastar)
\indistrb \Ncal(0, \MyIdDim).
\end{align}

We prove this claim shortly, but let us complete the proof of
Proposition <ref> using
Lemma <ref>.  Now, by construction we have $\Basis
\Basis^\top = \MyIdDim$ and $\Basis^\top\e_1 = \dir$.  Using these two
properties, we can write

\begin{align*}
\e_1^\top \Basis \thetastar = \dir^\top \thetastar  \;\;
\text{and} \;\; \y_i = \inprod{\Basis \x_i}{ \Basis \thetastar} +
\error_i \quad \text{for all} \;\; i = 1, \ldots, \numobs.
\end{align*}

Consequently, in this new basis, estimating the scalar $\dir^\top
\thetastar$ is same as estimating the first coordinate of transformed
vector $\Basis \thetastar$. Next, by construction of the matrix $\Basis$, we have

\begin{align}
\label{eqn:Cov-diag-entry}
\e_1^\top \diagCov_{\dir, \numobs}^{-1} \e_1 = \dir^\top \SigmaMat_{\numobs}^{-1} \dir
\qquad \text{and} \qquad \OpTuning_\numobs = \opnorm{\diagCov_{\dir, \numobs}^{-\frac{1}{2}} \SigmaMat_{\dir,\numobs}^{\frac{1}{2}}}.
\end{align}

Thus, we deduce

\begin{align*}
\e_1^\top \diagCov_{\dir,\numobs}^{\frac{1}{2}}(\thetaBasis - \Basis \thetastar)
& = \sqrt{(\diagCov_\numobs)_{11}} \cdot ( (\thetaBasis)_1 - \dir^\top \thetastar) \\
& = \sqrt{\frac{1}{(\diagCov_\numobs^{-1})_{11}}} \cdot (\e_1^\top \thetaBasis -  \dir^\top \thetastar) \\
& =  \sqrt{\frac{1}{\dir^\top \SigmaMat_\numobs^{-1} \dir} } \cdot
( \e_1^\top \thetaBasis -  \dir^\top \thetastar),
\end{align*}

The first equality above follows since the first row of the matrix $\diagCov_{\dir, \numobs}$ is proportional
to $\e_1$ by construction and the fact that $\e_1 \Basis \thetastar = \dir^\top \thetastar$.
The last line follows from the relation (<ref>). Thus, from  property (<ref>)
we deduce

\begin{align}
\label{eqn:one-dim-asn}
\sqrt{\frac{\tuneParScaled_\numobs}{\OpTuning_\numobs^2 \noisesd^2}} \cdot \sqrt{\frac{1}{\dir^\top \SigmaMat_\numobs^{-1} \dir} } \cdot
( \e_1^\top \thetaBasis -  \dir^\top \thetastar) \indistrb \Ncal(0, 1).
\end{align}

Define, the set $\CIset_{\dir, 1 - \alpha} \subseteq \real$ as

\begin{align*}
\CIset_{\dir, 1 - \alpha}  &\mydefn \left\{ \theta \in \real \mid  -\NormalQuantile{1 - \alpha/2} \leq \sqrt{\frac{\tuneParScaled_\numobs}{ \OpTuning_\numobs^2 \sigmahat^2}} \cdot \sqrt{\frac{1}{\dir^\top \SigmaMat_\numobs^{-1} \dir} } \cdot
( \e_1^\top \thetaBasis -  \theta) \leq \NormalQuantile{1 - \alpha/2} \right\} \\
&\equiv \Big[
\e_1^\top \thetaBasis -
\tfrac{\OpTuning_\numobs \cdot \sigmahat}{\sqrt{\tuneParScaled_\numobs}}
(\myinprod{\e_1}{\SigmaMat_\numobs^{-1} \e_1})^{\frac{1}{2}}
\NormalQuantile{1 - \alpha/2},
\quad
\e_1^\top \thetaBasis +
\tfrac{\OpTuning_\numobs \cdot \sigmahat}{\sqrt{\tuneParScaled_\numobs}}
(\myinprod{\e_1}{\SigmaMat_\numobs^{-1} \e_1})^{\frac{1}{2}}
\NormalQuantile{1 - \alpha/2} \Big]
\end{align*}

where, $\NormalQuantile{1 - \alpha/2}$ is the $1 - \alpha/2$ quantile of
the standard Gaussian random variable.
From the equation (<ref>) we have $ \lim_{\numobs
\rightarrow \infty} \Prob(\dir^\top \thetastar \in \CIset_{\dir, 1 -
\alpha}) = 1 - \alpha$, i.e., $\CIset_{\dir, 1 - \alpha}$ is an
asymptotically exact $1 - \alpha$ confidence intervals for $\dir^\top
\thetastar$. This completes the proof of the Proposition <ref>.
It remains to prove Lemma <ref>.

§.§ Proof of Lemma <ref>

Observe that

\begin{align*}
\y_i = \inprod{\Basis \x_i}{ \Basis \thetastar} +
\error_i \quad \text{for all} \;\; i = 1, \ldots, \numobs.
\end{align*}

The proof of the Lemma <ref>
is similar to the proof of Theorem <ref> but modified
for the data $\{ \Basis \x_i, \y_i \}_{i = 1}^\numobs$ and with $\thetastar$
replaced by $\Basis \thetastar$. Without loss of generality, we
assume that $\noisesd$ is known; thus, it suffices to

Recalling the expression for $\thetaBasis$ from the definition (<ref>) we have

\begin{align*}
\frac{\sqrt{\tuneParScaled_\numobs}}{\OpTuning_\numobs} \cdot
\diagCov_{\dir,\numobs}^{\frac{1}{2}}(\thetaBasis - \Basis\thetastar) &=
\sqrt{\tuneParScaled_\numobs} \cdot \sum_{i = 1}^\numobs \w_i \error_i\\
&  \quad \quad  +
\sqrt{\tuneParScaled_\numobs} \cdot
\left( \tfrac{1}{\OpTuning_\numobs} \cdot \diagCov_{\dir,\numobs}^{\frac{1}{2}} \SigmaMat_{\dir,\numobs}^{-\frac{1}{2}} -
\W_\numobs \X_{\dir,\numobs} \SigmaMat_{\dir,\numobs}^{-\frac{1}{2}}  \right)  \SigmaMat_{\dir,\numobs}^{\frac{1}{2}} (\thetaBasisLS - \Basis\thetastar) \\
&= \ZeroMartin_\numobs + \Bias_\numobs
\end{align*}

It remains to prove $\Bias_\numobs \stackrel{\smallop{p}}{\longrightarrow} 0$ and $\ZeroMartin_\numobs \indist \Ncal(0, \noisesd^2 \MyIdDim)$. Observe that

\begin{align*}
\enorm{\Bias_\numobs} &\leq \sqrt{\tuneParScaled_\numobs} \cdot \opnorm{ \tfrac{1}{\OpTuning_\numobs} \cdot  \diagCov_{\dir, \numobs}^{\frac{1}{2}} \SigmaMat_{\dir,\numobs}^{-\frac{1}{2}}-
\W_\numobs \X_{\dir,\numobs} \SigmaMat_{\dir,\numobs}^{-\frac{1}{2}}} \cdot
\enorm{\SigmaMat_{\dir, \numobs}^{\frac{1}{2}}(\thetaBasisLS - \Basis\thetastar)} \\
&\leq \sqrt{\tuneParScaled_\numobs \log
\lambda_{\max}(\SigmaMat_{\dir,\numobs})} \cdot \opnorm{ \tfrac{1}{\OpTuning_\numobs} \cdot \diagCov_\numobs^{\frac{1}{2}} \SigmaMat_{\dir,\numobs}^{-\frac{1}{2}}-
\W_\numobs \X_{\dir,\numobs} \SigmaMat_{\dir,\numobs}^{-\frac{1}{2}}}  \\
&\stackrel{\smallop{p}}{\longrightarrow} 0,
\end{align*}

where, the second inequality uses Theorem 1 from the paper [18], and the last step uses the
vanishing bias condition <ref>(b).

The analysis of the martingale term $\ZeroMartin_\numobs$ is exactly same as that of Theorem <ref> proof. This completes the proof of the Lemma <ref>.

§.§ Proof of Proposition <ref>

Recalling that $\maxnorm{\mymatrix{M}} \mydefn \max_{i,j}
|\mymatrix{M}_{ij}|$ denotes the maximum absolute entry of a matrix,
we claim that it suffices to show that $\maxnorm{\W_\numobs \X_\numobs
\XscaleMat_\numobs^{-\frac{1}{2}}} \leq 4$. Indeed, when this claim
holds, we have

\begin{align*}
\opnorm{\Id - \W_\numobs \X_\numobs \SigmaMat_\numobs^{-\frac{1}{2}}}
& \leq 1 + \opnorm{ \W_\numobs \X_\numobs
\XscaleMat_\numobs^{-\frac{1}{2}} } \opnorm{
\XscaleMat_\numobs^{\frac{1}{2}} \SigmaMat_\numobs^{-\frac{1}{2}}} \\
& \leq 1 + \BigoP(\sqrt{\Dim}) \opnorm{ \W_\numobs \X_\numobs
\XscaleMat_\numobs^{-\frac{1}{2}} } \\
& \leq 1 + \BigoP(\Dim^2) \maxnorm{\W_\numobs \X_\numobs
\XscaleMat_\numobs^{-\frac{1}{2}} }
\end{align*}

The second last inequality above follows by noting that the diagonal
entries of the matrix $\XscaleMat_\numobs^{\frac{1}{2}}
\SigmaMat_\numobs^{-1} \XscaleMat_\numobs^{\frac{1}{2}}$ is of the
order $\BigoP(1)$; this bound uses the expression of the scaling
matrix $\XscaleMat_\numobs$ from the
definition (<ref>), and the operator-norm bound
$\| \SigmaMatLB_\numobs^{\frac{1}{2}} \diag(\SigmaMat_\numobs^{-1})
\SigmaMatLB_\numobs^{\frac{1}{2}} \|_{op} = \BigoP(1)$ from
assumption (<ref>).  The last inequality above
follows from the fact that $\opnorm{ \mymatrix{A} } \leq
\Dim^{\frac{3}{2}} \maxnorm{ \mymatrix{A} }$, for any
$\Dim$-dimensional matrix $\mymatrix{A}$. This completes the proof of
Proposition <ref>. The remainder of the proof is
devoted to establishing an upper-bound on the max-norm of the matrix
$\W_\numobs \X_\numobs \XscaleMat_\numobs^{-\frac{1}{2}}$. We do so by
proving the following upper bounds

\begin{align}
\label{eqn:first-term-max-norm-bound}
\maxnorm{ \sum_{i =1 }^k \w_i \x_i^\top \XscaleMat_i^{-\frac{1}{2}}
} & \leq 2 \quad \text{for} \;\; k = 1, \ldots \numobs, \qquad
\text{and} \\
\label{eqn:second-term-max-norm-bound}
\maxnorm{ \Nsum{i} \w_i \x_i^\top \XscaleMat_i^{-\frac{1}{2}} (\Id
- \XscaleMat_i^{\frac{1}{2}} \XscaleMat_\numobs^{-\frac{1}{2}} ) }
& \leq 2.
\end{align}

Note that a combination of these two bounds implies that .

Accordingly, the remainder of our proof is devoted to
establishing the bounds (<ref>)
and (<ref>).

§.§ Proof of bound (<ref>)

Using the expression for the weight vector $\w_i$ from
equation (<ref>), we have

\begin{align*}
\Id - \sum_{i=1}^k \w_i \x_i \XscaleMat_i^{-\frac{1}{2}} = \prod_{i =
1}^k \left( \Id - \frac{\XscaleMat_i^{-\frac{1}{2}} \x_i \x_i^\top
\XscaleMat_i^{-\frac{1}{2}} }{\tfrac{\tuneParScaled_\numobs}{2} +
\|\XscaleMat_i^{-\frac{1}{2}} \x_i\|^2} \right).
\end{align*}

We claim

\begin{align*}
\left\| \Id - \frac{\XscaleMat_i^{-\frac{1}{2}} \x_i \x_i^\top
\XscaleMat_i^{-\frac{1}{2}} }{\tfrac{\tuneParScaled_\numobs}{2} +
\|\XscaleMat_i^{-\frac{1}{2}} \x_i\|^2} \right\|_{\mathrm{op}} \leq 1 \quad \text { for all } i \in[n] .
\end{align*}

It suffices to show that the maximum absolute eigenvalue of the symmetric matrix $\Id - \frac{\XscaleMat_i^{-\frac{1}{2}} \x_i \x_i^\top
\XscaleMat_i^{-\frac{1}{2}} }{\tfrac{\tuneParScaled_\numobs}{2} +
\|\XscaleMat_i^{-\frac{1}{2}} \x_i\|^2}$ is upper bounded by 1. Indeed, for any $u \in \real^d$ with $\|u\|_2 = 1$

\begin{align*}
0 \leq    u^\top \left( \Id - \frac{\XscaleMat_i^{-\frac{1}{2}} \x_i \x_i^\top
\XscaleMat_i^{-\frac{1}{2}} }{\tfrac{\tuneParScaled_\numobs}{2} +
\|\XscaleMat_i^{-\frac{1}{2}} \x_i\|^2}  \right) u = 1 - \frac{(\x_{i}^{\top} \XscaleMat_{i}^{-\frac{1}{2}} u)^2}{\gamma_{n} / 2+\left\|\XscaleMat_{i}^{-\frac{1}{2}} \x_{i}\right\|_{2}^{2}} \leq 1
\end{align*}

The above conclusions follow from the fact that for $\|u\|_2 \leq 1$

\begin{align*}
|\x_{i}^{\top} \XscaleMat_{i}^{-\frac{1}{2}} u| \leq
\left\|\XscaleMat_{i}^{-\frac{1}{2}} \x_{i}\right\|_{2}