# Grokking Modular Polynomials

Darshil Doshi a, b, †
&Tianyu He a, b
&Aritra Das b
&Andrey Gromov a, b

###### Abstract

Neural networks readily learn a subset of the modular arithmetic tasks, while
failing to generalize on the rest. This limitation remains unmoved by the
choice of architecture and training strategies. On the other hand, an
analytical solution for the weights of Multi-layer Perceptron (MLP) networks
that generalize on the modular addition task is known in the literature. In
this work, we (i) extend the class of analytical solutions to include modular
multiplication as well as modular addition with many terms. Additionally, we
show that real networks trained on these datasets learn similar solutions upon
generalization (grokking). (ii) We combine these “expert” solutions to
construct networks that generalize on arbitrary modular polynomials. (iii) We
hypothesize a classification of modular polynomials into learnable and non-
learnable via neural networks training; and provide experimental evidence
supporting our claims.

11footnotetext: Condensed Matter Theory Center, University of Maryland,
College Park22footnotetext: Department of Physics, University of Maryland,
College Park44footnotetext: Corresponding author

{ddoshi, tianyuh, aritrad<EMAIL_ADDRESS>
##  1 Introduction

Modular arithmetic finds applications in many important fields including
cryptography (Regev, 2009; Wenger et al., 2023), computer algebra (Fortin et
al., 2020), number theory, error detection in serial numbers and music (A.
AsKew & Klima, 2018). Recently, modular arithmetic has peaked the interest of
the machine learning community due to the phenomenon of _Grokking_ (Power et
al., 2022) – delayed and sudden generalization that occurs long after
memorization. Since it was first observed (Power et al., 2022), many works
have attempted to understand _grokking dynamics_ (Liu et al., 2022; 2023;
Thilak et al., 2022; Notsawo et al., 2023; Kumar et al., 2023; Lyu et al.,
2023; Davies et al., 2022; Minegishi et al., 2023). Another line of work
attempts to unravel the grokking behaviors using solvable models and
mechanistic interpretability (Nanda et al., 2023; Gromov, 2023; Zhong et al.,
2023; Barak et al., 2022; Merrill et al., 2023; Doshi et al., 2023; Varma et
al., 2023; Rubin et al., 2023). It was noted in Power et al. (2022) that
neural networks readily grok a subset of modular arithmetic polynomials while
failing to generalize on the rest. This limitation has remained unwavering to
variations in architectures and training strategies.

Notably, Gromov (2023) presented an “analytical solution” for the weights of a
2-layer MLP network that has a $100\%$ accuracy on modular addition dataset
($n_{1}+n_{2}\,\mathrm{mod}\,p$). They showed that real networks trained on
modular addition data finds similar solutions upon grokking. In this work, we
extend the class of analytical solutions to to include modular multiplication
(e.g. $n_{1}n_{2}\,\mathrm{mod}\,p$) and modular addition with many terms
(e.g. $n_{1}+n_{2}+\cdots n_{S}\,\mathrm{mod}\,p$).111The solution for modular
multiplication can be also readily extended to many terms. We present the
solution for two terms just for simplicity. Indeed, we show that training real
networks leads to similar network weights upon grokking for both of these
tasks. Using these analytical solutions as “experts”, we construct networks
that offers generalization on arbitrary modular polynomials, _including the
ones that are deemed un-learnable in literature_. This formulation opens up
the possibility of training Mixture-of-Experts (Jordan & Jacobs, 1993; Shazeer
et al., 2017; Lepikhin et al., 2021; Fedus et al., 2022) models that can be
trained to learn arbitrary modular arithmetic tasks. Based on our analysis, we
hypothesize a classification of modular polynomials into learnable and non-
learnable by 2-layer MLPs (and possibly more general architectures such as
Transformers and deeper MLPs).

##  2 Modular Addition with Many Terms

Consider the modular addition task with many terms with arbitrary
coefficients:
$\left(c_{1}n_{1}+c_{2}n_{2}+\cdots+c_{S}n_{S}\right)\,\mathrm{mod}\,p\,;$
where $c_{s}\in\mathbb{Z}_{p}\backslash\\{0\\}$ are the nonzero coefficients
of the modular variables $n_{s}\in\mathbb{Z}_{p}$. Note that this is a
generalization of the modular addition tasks generally considered in
literature: $(n_{1}+n_{2})\,\mathrm{mod}\,p$ (Power et al., 2022; Gromov,
2023). We consider a 2-layer MLP (of sufficient width) for this task.

$\displaystyle{\bm{f}}_{addS}({\bm{e}}_{n_{1}}\oplus\cdots\oplus{\bm{e}}_{n_{S}})={\bm{W}}\phi\left({\bm{U}}({\bm{e}}_{n_{1}}\oplus\cdots\oplus{\bm{e}}_{n_{S}})\right)={\bm{W}}\phi\left({\bm{U}}^{(1)}\,{\bm{e}}_{n_{1}}+\cdots+{\bm{U}}^{(S)}\,{\bm{e}}_{n_{S}}\right)\,,$
(1)

where ${\bm{e}}_{n_{1}},\dots,{\bm{e}}_{n_{S}}\in\mathbb{R}^{p}$ are one_hot
encoded numbers $n_{1},\dots,n_{S}$. “$\oplus$” denotes concatenation of
vectors
$({\bm{e}}_{n_{1}}\oplus\cdots\oplus{\bm{e}}_{n_{S}}\in\mathbb{R}^{Sp})$.
${\bm{U}}\in\mathbb{R}^{N\times Sp}$ and ${\bm{W}}\in\mathbb{R}^{p\times N}$
are the first and second layer weight matrices, respectively. $\phi(x)=x^{S}$
is the element-wise activation function. ${\bm{U}}$ is decomposed into $S$
blocks of $N\times P$:
${\bm{U}}={\bm{U}}^{(1)}\oplus\cdots\oplus{\bm{U}}^{(S)}$.
${\bm{U}}^{(1)},\dots,{\bm{U}}^{(S)}$ serve as embedding matrices for
$n_{1},\dots,n_{S}$.
${\bm{f}}({\bm{e}}_{n_{1}}\oplus\cdots\oplus{\bm{e}}_{n_{S}})\in\mathbb{R}^{p}$
is the network-output on one example datapoint $(n_{1},\dots,n_{S})$. The
targets are one_hot encoded answers
${\bm{e}}_{(c_{1}n_{1}+\cdots+c_{S}n_{S})\,\mathrm{mod}\,p}$.

###  2.1 Analytical solution

For a sufficiently wide network, it is possible to write a solution for the
weights of the network that generalizes with $100\%$ accuracy.

$\displaystyle\begin{aligned}
&U^{(s)}_{ki}=A\,\cos{\left[\frac{2\pi}{p}\sigma(k)c_{s}i+\psi^{(s)}_{k}\right]}&&\left(\forall\,s\in[0,S]\right)&\\\
&W_{qk}=A\,\cos{\left[-\frac{2\pi}{p}\sigma(k)q-\sum_{s=1}^{S}\psi^{(s)}_{k}\right]}\,,&\end{aligned}$
(2)

where $\sigma(k)$ denotes a random permutation of $k$ in $S_{N}$ – reflecting
the permutation symmetry of the hidden neurons. The phases $\psi_{k}^{(s)}$
uniformly i.i.d. sampled between $(-\pi,\pi]$. $A$ is the normalization factor
to ensure correct magnitude of the output logits, given by
$A={\left(2^{S}/(N\cdot S!)\right)}^{1/(S+1)}$. Note that this solution is a
generalization of that presented in Gromov (2023).

The rows (columns) of ${\bm{U}}^{(1)},\dots,{\bm{U}}^{(S)},{\bm{W}}$ are
periodic with frequencies
$\frac{2\pi}{p}\sigma(k)c_{1},\cdots,\frac{2\pi}{p}\sigma(k)c_{S},-\frac{2\pi}{p}\sigma(k)$.
We can use Inverse Participation Ratio (IPR) (Gromov, 2023; Doshi et al.,
2023) to quantify the similarity between the analytical solution and real
trained networks. We denote discrete Fourier transforms of the row- (column-)
vectors ${\bm{U}}^{(s)}_{k\cdot},{\bm{W}}_{\cdot k}$ by
$\mathcal{F}\left({\bm{U}}^{(s)_{k\cdot}}\right),\mathcal{F}({\bm{W}}_{\cdot
k})$.222In general IPR for a vector ${\bm{u}}$ is defined as
$(\lVert{\bm{u}}\rVert_{2r}/\lVert{\bm{u}}\rVert_{2})^{2r}$. We set $r=2$,
which is a common choice.

$\displaystyle\text{IPR}\left({\bm{U}}^{(s)}_{k\cdot}\right)\coloneqq\left(\frac{\left\lVert\mathcal{F}\left({\bm{U}}^{(s)}_{k\cdot}\right)\right\rVert_{4}}{\left\lVert\mathcal{F}\left({\bm{U}}^{(s)}_{k\cdot}\right)\right\rVert_{2}}\right)^{4}\,;$
$\displaystyle\text{IPR}\left({\bm{W}}_{\cdot
k}\right)\coloneqq\left(\frac{\left\lVert\mathcal{F}\left({\bm{W}}_{\cdot
k}\right)\right\rVert_{4}}{\left\lVert\mathcal{F}\left({\bm{W}}_{\cdot
k}\right)\right\rVert_{2}}\right)^{4}\,;$ (3)

where $\lVert\cdot\rVert_{P}$ denotes the $L^{P}$-norm of the vector.
Subscripts “$k\cdot$” and “$\cdot k$” denote $k^{th}$ row and column vectors
respectively. We also define the per-neuron IPR as

$\displaystyle\text{IPR}_{k}\coloneqq\frac{1}{S+1}\left(\text{IPR}\left({\bm{U}}^{(1)}_{k\cdot}\right)+\cdots+\text{IPR}\left({\bm{U}}^{(S)}_{k\cdot}\right)+\text{IPR}\left({\bm{W}}_{\cdot
k}\right)\right)\,,$ (4)

and the average IPR of the network by averaging it over all neurons as
$\overline{IPR}\coloneqq\mathbb{E}_{k}\left[\text{IPR}_{k}\right]$.

In Figure 1, we show that the analytical solution equation 2 indeed has
$100\%$ for sufficient network-width. We observe an exponential increase in
the required width $N$ due to an exponential increase in the number of cross-
terms upon expansion; which needs to be suppressed by the factor $1/N$. We
refer the reader to Appendix A for a detailed discussion.

Figure 1: Modular addition with many terms – analytical solution applied to
real 2-layer MLP networks ($p=23$). The solution equation 2 works for
sufficiently wide networks. Note that the x-axis is scaled logarithmically;
which suggests an exponential increase in the required width upon adding more
terms (as expected). The accuracies shown are calculated over a randomly
chosen subset of the entire dataset, consisting of 10k examples. The results
shown are the best out of 10 random seeds.

###  2.2 Comparison with trained networks

We compare the above solution to networks trained on modular addition data. We
see that 2-layer MLP with sufficient width can indeed grok on modular
arithmetic data. Moreover, the trained network learns qualitatively similar
weights to the analytical solution equation 2. The network generalizes by
learning the periodic weights, as evidenced by the initial and final
$\overline{\mathrm{IPR}}$ distributions in Figure 2.

Figure 2: Training on modular addition with many terms
($(n_{1}+n_{2}+n_{3}+n_{4})\,\mathrm{mod}\,p$). $p=11;N=5000$; Adam optimizer;
learning rate $=0.005$; weight decay $=5.0$; $50\%$ of the dataset used for
training. (a) MSE loss on train and test dataset. (b) Accuracy on train and
test dataset as well as average IPR of the network $\overline{\mathrm{IPR}}$.
The training curves show the well-known grokking phenomenon; while
$\overline{\mathrm{IPR}}$ monotonically increases. (c) Initial and final IPR
distributions, evidently showing periodic neurons in the grokked network,
confirming the similarity to equation 2. Note that the IPR for the analytical
solution (equation 2) is 1.

##  3 Modular Multiplication

Now consider the modular multiplication task in two variables333The analysis
can be readily extended to modular multiplication in more than two variables;
in similar vein as the modular addition analysis present in the previous
section.: $(n_{1}^{a}\,n_{2}^{b})\,\mathrm{mod}\,p$ (where
$a,b\in\mathbb{Z}_{p}\backslash\\{0\\}$; $p$ is a prime number). It can be
learned by a 2-layer MLP, in a similar fashion as modular addition (Gromov,
2023; Doshi et al., 2023).

$\displaystyle{\bm{f}}_{mul2}({\bm{e}}_{n_{1}}\oplus{\bm{e}}_{n_{2}})={\bm{Q}}\phi\left({\bm{P}}({\bm{e}}_{n_{1}}\oplus{\bm{e}}_{n_{2}})\right)={\bm{Q}}\phi\left({\bm{P}}^{(1)}\,{\bm{e}}_{n_{1}}+{\bm{P}}^{(2)}\,{\bm{e}}_{n_{2}}\right)\,,$
(5)

where we have used identical notation to equation 1.

###  3.1 Exponential and logarithmic maps over finite fields

Note that both modular addition and multiplication are defined over finite
field of prime order $p$, denoted by $\mathcal{GF}(p)$. However, there is a
subtle yet crucial difference in the two operations. Addition has a cycle of
length $p$; i.e.
$\left(p+p+\cdots(n\,\text{times})\right)\,\mathrm{mod}\,p=0\;\forall\;n\in\mathbb{Z}_{p}$).
On the other hand, multiplication has a cycle of length $p-1$; i.e.
$\left(g\cdot g\cdot\cdots(n-1\,\text{times})\right)\;\mathrm{mod}\;p=1$,
where $g$ is a primitive root of $\mathcal{GF}(p)$.444This is a well-known
result called Fermat’s little theorem. Note that multiplying any number by $0$
gives $0$; which warrants a separate treatment for $0$. The remaining elements
$r\in\mathbb{Z}_{p}\backslash\\{0\\}$ can be used to define the map
$g^{r}:\mathbb{Z}_{p}\backslash\\{0\\}\rightarrow\mathbb{Z}_{p}\backslash\\{0\\}$,
where $g$ is a primitive root. This exponential is bijective (invertible).
Inverting this map defines the logarithmic map with base $g$. We emphasise
that over the finite field $\mathcal{GF}(p)$, the exponential ($g^{r}$) and
logarithmic ($\log_{g}r$) are basically a reshuffling of the elements of
$\mathbb{Z}_{p}\backslash\\{0\\}$. Reshuffling the numbers according to the
logarithmic map turns the multiplication task into the addition task. This is
the equivalent of the result
$\log_{g}(n_{1}n_{2})=\log_{g}n_{1}+\log_{g}n_{2}$ on finite fields. Thus,
barring the application of this map, the solution for multiplication is
similar to equation 6. Note that the element $0$ needs to be treated
separately while constructing the analytical solution for the network weights.

###  3.2 Analytical solution

Using the above insights, we now construct the analytical solution for the
modular multiplication task. Note that in the following equation, $i\neq
0;j\neq 0;k\neq 0,q\neq 0$.

$\displaystyle\begin{aligned} &P^{(1)}_{00}=1\;,\qquad P^{(2)}_{00}=1\;,\qquad
Q_{00}=1\;,&\\\ &P^{(1)}_{0i}=P^{(1)}_{k0}=0\;,\qquad
P^{(2)}_{0j}=P^{(2)}_{k0}=0\;,\qquad Q_{q0}=Q_{k0}=0\;,&\\\
&{P}^{(1)}_{ki}=\left(\frac{2}{N-1}\right)^{-\frac{1}{3}}\cos{\left[\frac{2\pi}{p-1}\sigma(k)(a\log_{g}{i})+\psi^{(1)}_{k}\right]}\,,&\\\
&{P}^{(2)}_{kj}=\left(\frac{2}{N-1}\right)^{-\frac{1}{3}}\cos{\left[\frac{2\pi}{p-1}\sigma(k)(b\log_{g}{j})+\psi^{(2)}_{k}\right]}\,,&\\\
&{Q}_{qk}=\left(\frac{2}{N-1}\right)^{-\frac{1}{3}}\cos{\left[-\frac{2\pi}{p-1}\sigma(k)(\log_{g}{q})-\psi^{(1)}_{k}-\psi^{(2)}_{k}\right]}\;,&\end{aligned}$
(6)

where $\sigma(k)$ denotes a random permutation of $k$ in $S_{N}$ – reflecting
the permutation symmetry of the hidden neurons. The phases $\psi_{k}^{(1)}$
and $\psi_{k}^{(2)}$ are uniformly i.i.d. sampled between $(-\pi,\pi]$.

The rows (columns) of ${\bm{P}}^{(1)},{\bm{P}}^{(2)},{\bm{Q}}$ are periodic
with frequencies
$\frac{2\pi}{p-1}a\sigma(k),\frac{2\pi}{p-1}b\sigma(k),-\frac{2\pi}{p-1}\sigma(k)$,
upon performing exponential map on their column (row) indices $i,j,q$ (Note
that this excludes the $0^{th}$ entry in each row (column). Those entries deal
with the element $0$). This mapping is equivalent to reshuffling the columns
(rows) of ${\bm{P}}^{(1)},{\bm{P}}^{(2)},{\bm{Q}}$. Let us denote these new
re-shuffled, 0-excluded weight matrices as
$\overline{{\bm{P}}}^{(1)},\overline{{\bm{P}}}^{(2)},\overline{{\bm{Q}}}$.

$\displaystyle\overline{P}^{(1)}_{ki}=P^{(1)}_{k,g^{i}}\;,$
$\displaystyle\overline{P}^{(2)}_{kj}=P^{(2)}_{k,g^{j}}\;,$
$\displaystyle\overline{Q}_{qk}=Q_{g^{q},k}$ $\displaystyle(i\neq 0,j\neq
0,q\neq 0,k\neq 0)\,,$ (7)

where we denote the shuffled indices ($i,j,q$) by $g^{i},g^{j},g^{q}$ in the
subscript. Again, we can use IPR quantify the periodicity. We denote discrete
Fourier transforms of the row- (column-) vectors
$\overline{{\bm{P}}}^{(1)}_{k\cdot},\overline{{\bm{P}}}^{(2)}_{k\cdot},\overline{{\bm{Q}}}_{\cdot
k}$ by
$\mathcal{F}\left(\overline{{\bm{P}}}^{(1)}_{k\cdot}\right),\mathcal{F}\left(\overline{{\bm{P}}}^{(2)}_{k\cdot}\right),\mathcal{F}\left(\overline{{\bm{Q}}}_{\cdot
k}\right)$.

$\displaystyle\text{IPR}\left(\overline{{\bm{P}}}^{(t)}_{k\cdot}\right)\coloneqq\left(\frac{\left\lVert\mathcal{F}\left(\overline{{\bm{P}}}^{(t)}_{k\cdot}\right)\right\rVert_{4}}{\left\lVert\mathcal{F}\left(\overline{{\bm{P}}}^{(t)}_{k\cdot}\right)\right\rVert_{4}}\right)^{2}\,;$
$\displaystyle\text{IPR}\left(\overline{{\bm{Q}}}_{\cdot
k}\right)\coloneqq\left(\frac{\left\lVert\mathcal{F}\left(\overline{{\bm{Q}}}_{\cdot
k}\right)\right\rVert_{4}}{\left\lVert\mathcal{F}\left(\overline{{\bm{Q}}}_{\cdot
k}\right)\right\rVert_{2}}\right)^{4}\,,$ (8)

where $t\in\\{1,2\\}$. Per-neuron IPR as well as average IPR of the network
can be defined as before; with $\text{IPR}_{k}\coloneqq
1/3\left(\text{IPR}\left(\overline{{\bm{P}}}^{(1)}_{k\cdot}\right)+\text{IPR}\left(\overline{{\bm{P}}}^{(2)}_{k\cdot}\right)+\text{IPR}\left(\overline{{\bm{Q}}}_{k\cdot}\right)\right)$
and $\overline{\text{IPR}}\coloneqq\mathbb{E}_{k}\left[\text{IPR}_{k}\right]$.

###  3.3 Comparison with trained networks

Now, we show that the real networks trained on modular multiplication data
learn similar features as Equation 6.

Figure 3: Training on modular multiplication ($n_{1}n_{2}\,\mathrm{mod}\,p$).
$p=97;N=500$; Adam optimizer; learning rate $=0.005$; weight decay $=5.0$;
$50\%$ of the dataset used for training. (a) MSE loss on train and test
datset. (b) Accuracy on train and test dataset as well as average IPR of the
network $\overline{\mathrm{IPR}}$. The training curves show the well-known
grokking phenomenon; while $\overline{\mathrm{IPR}}$ monotonically increases.
(c) Initial and final IPR distributions, evidently showing periodic neurons in
the grokked network, confirming the similarity to equation 6. Note that the
IPR for the analytical solution (equation 6) is 1.

##  4 Arbitrary Modular Polynomials

Consider a general modular polynomial in two variables $(n_{1},n_{2})$
containing $S$ terms:

$\displaystyle\left(c_{1}n_{1}^{a_{1}}n_{2}^{b_{1}}+c_{2}n_{1}^{a_{2}}n_{2}^{b_{2}}+\cdots+c_{S}n_{1}^{a_{S}}n_{2}^{b_{S}}\right)\,\mathrm{mod}\,p\,.$
(9)

We utilize the solutions presented in previous sections to construct a simple
network that generalizes on this task. Each term (without the coefficients
$c_{s}$) can be solved by a 2-layer MLP expert (equations 5, 6) that
multiplies the appropriate powers of $n_{1}$ and $n_{2}$. The output of these
$S$ terms can be added, along with the coefficients, using another expert
network (equations 1. The expert networks are designed to perform best when
their inputs are close to one_hot vectors; so we apply (low temperature)
Softmax function to the outputs of the each term. As before, the input to the
network are stacked, one_hot represented numbers $n_{1},n_{2}$:
$({\bm{e}}_{n_{1}}\oplus{\bm{e}}_{n_{2}})$

$\displaystyle{\bm{t}}^{(s)}={\bm{f}}^{(s)}_{mul2}({\bm{e}}_{n_{1}}\oplus{\bm{e}}_{n_{2}})={\bm{Q}}^{(s)}\phi\left({\bm{P}}^{(s,1)}\,{\bm{e}}_{n_{1}}+{\bm{P}}^{(s,2)}\,{\bm{e}}_{n_{2}}\right)\,,$
$\displaystyle(s\in[1,S])$ (10)
$\displaystyle{\bm{u}}^{(s)}=\texttt{softmax}_{\beta}\left({\bm{t}}^{(s)}\right)\,,$
(11)
$\displaystyle{\bm{z}}={\bm{f}}_{addS}\left({\bm{u}}^{(1)}\oplus\cdots\oplus{\bm{u}}^{(S)}\right)={\bm{W}}\phi\left({\bm{U}}^{(1)}\,{\bm{u}}^{(1)}+\cdots+{\bm{U}}^{(S)}\,{\bm{u}}^{(S)}\right)$
(12)

The weights ${\bm{P}}^{(s,1)},{\bm{P}}^{(s,2)},{\bm{Q}}^{(s)}$ are given by
equation 6, with appropriate powers $a_{s},b_{s}$ for each term. Similarly,
the weights ${\bm{U}}^{(s)},{\bm{W}}$ are given by equation 2, with
coefficients $c_{s}$. We have used the temperature-scaled Softmax function

$\texttt{softmax}_{\beta}\left(t^{(s)}_{i}\right)\coloneqq\frac{e^{\beta
t^{(s)}_{i}}}{\sum_{j}^{p}e^{\beta t^{(s)}_{j}}}\,.$ (13)

We select a high value for the inverse temperature $\beta\sim 100$ so that the
intermediate outputs ${\bm{u}}$ get close to one_hot vectors. This results in
a higher accuracy in the summation of monomials to be performed in the
subsequent layers.

In Appendix A, we show the performance of real networks, with the solutions
equation 10, on general modular polynomials. We show that such a construction
is able to learn polynomials that are un-learnable via training standard
network architectures.

Instead of using the analytical weights, one can train the “experts”:
${\bm{f}}^{(s)}_{addS},{\bm{f}}_{addS}$ separately, on multiplication and
addition datasets, and then combine them into a Mixture-of-Experts model.

##  5 Discussion

We have presented the analytical solutions for the weights of 2-layer MLP
networks on modular multiplication using the bijective exponential and
logarithmic maps on $\mathcal{GF}(p)$. We have also presented solutions for
modular multiplication with many terms and arbitrary coefficients. In both
cases, We have shown that real networks trained on these tasks find similar
solutions.

Using these “expert” solutions, we have constructed a network that generalizes
on arbitrary modular polynomials, given sufficient width. This network is
reminiscent of Mixture-of-Expert models. Therefore, our construction opens the
possibility of building such network architectures, that can learn modular
arithmetic tasks that SOTA models have been unable to learn. We delegate the
exploration of this avenue for future work.

On the other hand, our general formulation of analytical solutions for 2-layer
MLP points to a potential classification of modular polynomials into learnable
and non-learnable ones.

######  Hypothesis 5.1 (Weak generalization).

Consider training a 2-layer MLP network trained on a dataset consisting of
modular polynomial on $\mathcal{GF}(p)$ in two variables $(n_{1},n_{2})$, with
commonly used optimizers (SGD, Adam etc.), loss functions (MSE, CrossEntropy
etc.) and regularization methods (weight decay, dropout, BatchNorm). The
network achieves 100% test accuracy if the modular polynomial is of the
following form:555The networks are naturally trained on $<100\%$ of the entire
dataset; and tested on the rest.

$\displaystyle h(g_{1}(n_{1})+g_{2}(n_{2}))\,\mathrm{mod}\,p\,,$ (14)

where $g_{1}$, $g_{2}$ and $h$ are functions on $\mathcal{GF}(p)$.

  * •

Note that equation 14 also includes modular multiplication: taking
$g_{1},g_{2}$ to be $\log_{g}(\cdot)$ and $h$ to be $exp_{g}(\cdot)$, where
$g$ is a primitive root of $\mathcal{GF}(p)$.

  * •

If the function $h$ is invertible, then it is also possible to construct an
analytical solution for this task in a similar fashion as equations 2 and
6.666This general construction of analytical solution was also hypothesized in
Gromov (2023).

In Appendix C, we provide experimental evidence supporting this claim. We show
the results of training experiments on various modular polynomials belonging
to both categories. We see a clear difference in test performance on tasks
that fall under the form equation 14 and those that do not.

It is possible to generalize this claim to other network architectures.
General architectures such as transformers and deeper MLPs have also been
unable to solve general modular arithmetic tasks that do not fall under
equation 14 (Power et al., 2022). Moreover, Doshi et al. (2023) hypothesized a
general framework for generalization on modular arithmetic tasks in these
architectures. Consequently, we conjecture that Hypothesis 5.1 can be extended
to general architectures.

######  Hypothesis 5.2 (Strong generalization).

Consider training a standard neural network (MLP, Transformer etc.; not pre-
trained) on a dataset consisting of modular polynomial on $\mathcal{GF}(p)$ in
two variables $(n_{1},n_{2})$, with commonly used optimizers (SGD, Adam etc.),
loss functions (MSE, CrossEntropy etc.) and regularization methods (weight
decay, dropout, BatchNorm). The network achieves 100% test accuracy if the
modular polynomial is of the from of equation 14.

While experiments validate these hypotheses, a comprehensive proof of the
claim remains an open challenge. Formulating such a proof would require a
robust classification of modular polynomials; as well as the understanding
biases of neural network training on these tasks. We defer such an analysis to
future works.

#### Acknowledgments

A.G.’s work at the University of Maryland was supported in part by NSF CAREER
Award DMR-2045181, Sloan Foundation and the Laboratory for Physical Sciences
through the Condensed Matter Theory Center.

## References

  * A. AsKew & Klima (2018) K. Kennedy A. AsKew and V. Klima.  Modular arithmetic and microtonal music theory.  _PRIMUS_ , 28(5):458–471, 2018.  doi: 10.1080/10511970.2017.1388314.
  * Barak et al. (2022) Boaz Barak, Benjamin L. Edelman, Surbhi Goel, Sham M. Kakade, eran malach, and Cyril Zhang.  Hidden progress in deep learning: SGD learns parities near the computational limit.  In Alice H. Oh, Alekh Agarwal, Danielle Belgrave, and Kyunghyun Cho (eds.), _Advances in Neural Information Processing Systems_ , 2022.  URL https://openreview.net/forum?id=8XWP2ewX-im.
  * Davies et al. (2022) Xander Davies, Lauro Langosco, and David Krueger.  Unifying grokking and double descent.  In _NeurIPS ML Safety Workshop_ , 2022.  URL https://openreview.net/forum?id=JqtHMZtqWm.
  * Doshi et al. (2023) Darshil Doshi, Aritra Das, Tianyu He, and Andrey Gromov.  To grok or not to grok: Disentangling generalization and memorization on corrupted algorithmic datasets, 2023.
  * Fedus et al. (2022) William Fedus, Barret Zoph, and Noam Shazeer.  Switch transformers: Scaling to trillion parameter models with simple and efficient sparsity, 2022.
  * Fortin et al. (2020) Pierre Fortin, Ambroise Fleury, François Lemaire, and Michael Monagan.  High performance simd modular arithmetic for polynomial evaluation, 2020\.
  * Gromov (2023) Andrey Gromov.  Grokking modular arithmetic, 2023.  URL https://arxiv.org/abs/2301.02679.
  * Jordan & Jacobs (1993) M.I. Jordan and R.A. Jacobs.  Hierarchical mixtures of experts and the em algorithm.  In _Proceedings of 1993 International Conference on Neural Networks (IJCNN-93-Nagoya, Japan)_ , volume 2, pp. 1339–1344 vol.2, 1993.  doi: 10.1109/IJCNN.1993.716791.
  * Kumar et al. (2023) Tanishq Kumar, Blake Bordelon, Samuel J. Gershman, and Cengiz Pehlevan.  Grokking as the transition from lazy to rich training dynamics, 2023.
  * Lepikhin et al. (2021) Dmitry Lepikhin, HyoukJoong Lee, Yuanzhong Xu, Dehao Chen, Orhan Firat, Yanping Huang, Maxim Krikun, Noam Shazeer, and Zhifeng Chen.  {GS}hard: Scaling giant models with conditional computation and automatic sharding.  In _International Conference on Learning Representations_ , 2021.  URL https://openreview.net/forum?id=qrwe7XHTmYb.
  * Liu et al. (2022) Ziming Liu, Ouail Kitouni, Niklas S Nolte, Eric Michaud, Max Tegmark, and Mike Williams.  Towards understanding grokking: An effective theory of representation learning.  In S. Koyejo, S. Mohamed, A. Agarwal, D. Belgrave, K. Cho, and A. Oh (eds.), _Advances in Neural Information Processing Systems_ , volume 35, pp. 34651–34663. Curran Associates, Inc., 2022.  URL https://proceedings.neurips.cc/paper_files/paper/2022/file/dfc310e81992d2e4cedc09ac47eff13e-Paper-Conference.pdf.
  * Liu et al. (2023) Ziming Liu, Eric J Michaud, and Max Tegmark.  Omnigrok: Grokking beyond algorithmic data.  In _The Eleventh International Conference on Learning Representations_ , 2023.  URL https://openreview.net/forum?id=zDiHoIWa0q1.
  * Lyu et al. (2023) Kaifeng Lyu, Jikai Jin, Zhiyuan Li, Simon S. Du, Jason D. Lee, and Wei Hu.  Dichotomy of early and late phase implicit biases can provably induce grokking, 2023.
  * Merrill et al. (2023) William Merrill, Nikolaos Tsilivis, and Aman Shukla.  A tale of two circuits: Grokking as competition of sparse and dense subnetworks, 2023.
  * Minegishi et al. (2023) Gouki Minegishi, Yusuke Iwasawa, and Yutaka Matsuo.  Grokking tickets: Lottery tickets accelerate grokking, 2023.
  * Nanda et al. (2023) Neel Nanda, Lawrence Chan, Tom Lieberum, Jess Smith, and Jacob Steinhardt.  Progress measures for grokking via mechanistic interpretability, 2023\.
  * Notsawo et al. (2023) Pascal Jr. Tikeng Notsawo, Hattie Zhou, Mohammad Pezeshki, Irina Rish, and Guillaume Dumas.  Predicting grokking long before it happens: A look into the loss landscape of models which grok, 2023.
  * Power et al. (2022) Alethea Power, Yuri Burda, Harri Edwards, Igor Babuschkin, and Vedant Misra.  Grokking: Generalization beyond overfitting on small algorithmic datasets, 2022.
  * Regev (2009) Oded Regev.  On lattices, learning with errors, random linear codes, and cryptography.  _J. ACM_ , 56(6), sep 2009.  ISSN 0004-5411.  doi: 10.1145/1568318.1568324.  URL https://doi.org/10.1145/1568318.1568324.
  * Rubin et al. (2023) Noa Rubin, Inbar Seroussi, and Zohar Ringel.  Droplets of good representations: Grokking as a first order phase transition in two layer networks, 2023.
  * Shazeer et al. (2017) Noam Shazeer, *Azalia Mirhoseini, *Krzysztof Maziarz, Andy Davis, Quoc Le, Geoffrey Hinton, and Jeff Dean.  Outrageously large neural networks: The sparsely-gated mixture-of-experts layer.  In _International Conference on Learning Representations_ , 2017.  URL https://openreview.net/forum?id=B1ckMDqlg.
  * Thilak et al. (2022) Vimal Thilak, Etai Littwin, Shuangfei Zhai, Omid Saremi, Roni Paiss, and Joshua Susskind.  The slingshot mechanism: An empirical study of adaptive optimizers and the grokking phenomenon, 2022.
  * Varma et al. (2023) Vikrant Varma, Rohin Shah, Zachary Kenton, János Kramár, and Ramana Kumar.  Explaining grokking through circuit efficiency, 2023.
  * Wenger et al. (2023) Emily Wenger, Mingjie Chen, François Charton, and Kristin Lauter.  Salsa: Attacking lattice cryptography with transformers, 2023.
  * Zhong et al. (2023) Ziqian Zhong, Ziming Liu, Max Tegmark, and Jacob Andreas.  The clock and the pizza: Two stories in mechanistic explanation of neural networks, 2023.

##  Appendix A Analytical solutions

###  A.1 Modular multiplication

Here we show that that the analytical solution presented in equation 6 solves
the modular multiplication task. Let us first calculate the network output for
$n_{1}\neq 0,n_{2}\neq 0$. We will treat the cases with $n_{1}=0$ and/or
$n_{2}=0$ separately.

$\displaystyle{\bm{f}}_{mul2}({\bm{e}}_{n_{1}}\oplus{\bm{e}}_{n_{2}})$
$\displaystyle={\bm{Q}}\phi\left({\bm{P}}^{(1)}\,{\bm{e}}_{n_{1}}+{\bm{P}}^{(2)}\,{\bm{e}}_{n_{2}}\right)$
$\displaystyle=\sum_{k=1}^{N}Q_{qk}(P^{(1)}_{kn_{1}}+P^{(2)}_{kn_{2}})^{2}$
$\displaystyle=\frac{2}{N-1}\sum_{k=1}^{N-1}\cos{\left(-\frac{2\pi}{p-1}\sigma(k)\log_{g}q-(\psi^{(1)}_{k}+\psi^{(2)}_{k})\right)}\cdot$
$\displaystyle\qquad\qquad\cdot\left(\cos{\left(\frac{2\pi}{p-1}\sigma(k)a\log_{g}n_{1}+\psi^{(1)}_{k}\right)}+\cos{\left(\frac{2\pi}{p-1}\sigma(k)b\log_{g}n_{2}+\psi^{(2)}_{k}\right)}\right)^{2}$
$\displaystyle=\frac{2}{N-1}\sum_{k=1}^{N-1}\left\\{\frac{1}{4}\cos{\left(\frac{2\pi}{p}\sigma(k)(2a\log_{g}n_{1}-\log_{g}q)+\psi^{(1)}_{k}-\psi^{(2)}_{k}\right)}\right.$
$\displaystyle\left.\qquad\qquad\quad+\frac{1}{4}\cos{\left(\frac{2\pi}{p-1}\sigma(k)(2a\log_{g}n_{1}+\log_{g}q)+3\psi^{(1)}_{k}+\phi^{(2)}_{k}\right)}\right.$
$\displaystyle\left.\qquad\qquad\quad+\frac{1}{4}\cos{\left(\frac{2\pi}{p-1}\sigma(k)(2b\log_{g}n_{2}-\log_{g}q)-\psi^{(1)}_{k}+\psi^{(2)}_{k}\right)}\right.$
$\displaystyle\left.\qquad\qquad\quad+\frac{1}{4}\cos{\left(\frac{2\pi}{p-1}\sigma(k)(2b\log_{g}n_{2}+\log_{g}q)+\psi^{(1)}_{k}+3\psi^{(2)}_{k}\right)}\right.$
$\displaystyle\left.\qquad\qquad\quad+\boxed{\frac{1}{2}\cos{\left(\frac{2\pi}{p-1}\sigma(k)(a\log_{g}n_{1}+b\log_{g}n_{2}-\log_{g}q)\right)}}\right.$
$\displaystyle\left.\qquad\qquad\quad+\frac{1}{2}\cos{\left(\frac{2\pi}{p-1}\sigma(k)(a\log_{g}n_{1}+b\log_{g}n_{2}+\log_{g}q)+2\psi^{(1)}_{k}+2\psi^{(2)}_{k}\right)}\right.$
$\displaystyle\left.\qquad\qquad\quad+\frac{1}{2}\cos{\left(\frac{2\pi}{p-1}\sigma(k)(a\log_{g}n_{1}-b\log_{g}n_{2}-\log_{g}q)-2\psi^{(1)}_{k}\right)}\right.$
$\displaystyle\left.\qquad\qquad\quad+\frac{1}{2}\cos{\left(\frac{2\pi}{p-1}\sigma(k)(a\log_{g}n_{1}-b\log_{g}n_{2}+\log_{g}q)+2\psi^{(1)}_{k}\right)}\right.$
$\displaystyle\left.\qquad\qquad\quad+\cos{\left(-\frac{2\pi}{p-1}\sigma(k)\log_{g}q-\psi^{(1)}_{k}-\psi^{(2)}_{k}\right)}\right\\}\,.$
(15)

We have highlighted the term that will give us the desired output with a box.
Note that the desired term is the only one that does not have additive phases
in the argument. Recall that the phases $\psi^{(1)}$ and $\psi^{(2)}$ are
randomly chosen – uniformly iid sampled between $(-\pi,\pi]$. Consequently, as
$N$ becomes large, all other terms will vanish due to random phase
approximation. The only surviving term will be the boxed term. We can write
the boxed term in a more suggestive form to make the analytical solution
apparent.

$\displaystyle{\bm{f}}_{mul2}({\bm{e}}_{n_{1}}\oplus{\bm{e}}_{n_{2}})$
$\displaystyle=\frac{1}{N-1}\sum_{k=1}^{N-1}\cos{\left(\frac{2\pi}{p-1}\sigma(k)(a\log_{g}n_{1}+b\log_{g}n_{2}-\log_{g}q)\right)}$
$\displaystyle=\frac{1}{N-1}\sum_{k=1}^{N-1}\cos{\left(\frac{2\pi}{p-1}\sigma(k)(\log_{g}n_{1}^{a}+\log_{g}n_{2}^{b}-\log_{g}q)\right)}$
$\displaystyle\sim\delta^{p}\left(n_{1}^{a}n_{2}^{b}-q\right)\,,$ (16)

were we have defined the _modular Kronecker Delta function_ $\delta(\cdot)$ as
Kronecker Delta function up to integer multiplies of the modular base $p$.

$\displaystyle\delta^{p}(x)=\begin{cases}1&\;x=rp\qquad(r\in\mathbb{Z})\\\
0&\;\text{otherwise}\end{cases}\,,$ (17)

where $\mathbb{Z}$ denotes the set of all integers. Note that
$\delta^{p}(n_{1}^{a}n_{2}^{b}-q)$ are the desired one_hot encoded labels for
the modular multiplication task, by definition. Thus our network output with
periodic weights is indeed a solution.

Note that if either $n_{1}=0$ or $n_{2}=0$ (or both), the $0^{th}$ output
logit will be the largest ($=1$). Consequently, the network will correctly
predict 0 output for such cases.

Hence, the solution presented in equation 6 gives $100\%$ accuracy .

###  A.2 Modular addition with many terms

Next, we show that that the analytical solution presented in equation 2 solves
the task of modular addition with many terms.

$\displaystyle{\bm{f}}_{addS}({\bm{e}}_{n_{1}}\oplus\cdots\oplus{\bm{e}}_{n_{S}})$
$\displaystyle={\bm{W}}\phi\left({\bm{U}}^{(1)}\,{\bm{e}}_{n_{1}}+\cdots+{\bm{U}}^{(S)}\,{\bm{e}}_{n_{S}}\right)$
$\displaystyle=\sum_{k=1}^{N}W_{qk}(U^{(1)}_{kn_{1}}+\cdots
P^{(S)}_{kn_{S}})^{S}$ $\displaystyle=\frac{2^{S}}{N\cdot
S!}\sum_{k=1}^{N}\cos{\left(-\frac{2\pi}{p}\sigma(k)q-\sum_{s=1}^{S}\psi^{(s)}_{k}\right)}\cdot$
$\displaystyle\qquad\qquad\cdot\left(\cos{\left(\frac{2\pi}{p}\sigma(k)c_{1}n_{1}+\psi^{(1)}_{k}\right)}+\cdots+\cos{\left(\frac{2\pi}{p}\sigma(k)c_{2}n_{2}+\psi^{(S)}_{k}\right)}\right)^{S}$
$\displaystyle=\frac{2^{S}}{N\cdot
S!}\sum_{k=1}^{N}\left\\{\cdots+\boxed{\frac{S!}{2^{S}}\cos{\left(\frac{2\pi}{p}\sigma(k)(c_{1}n_{1}+\cdots+c_{S}n_{S}-q)\right)}}+\cdots\right\\}$
(18)

Here, we have omitted all the additional terms that drop out due to random
phase approximation; showing only the desired surviving term. Note that the
number of these additional terms increases exponentially in $S$; and they are
suppressed by a factor of $1/N$. This provides an explanation for an
exponential behaviour in the required network width $N$ with increasing number
of terms in the addition task Figure 1.

$\displaystyle{\bm{f}}_{addS}({\bm{e}}_{n_{1}}\oplus\cdots\oplus{\bm{e}}_{n_{S}})$
$\displaystyle=\frac{1}{N}\sum_{k=1}^{N}\cos{\left(\frac{2\pi}{p}\sigma(k)(c_{1}n_{1}+\cdots+c_{S}n_{S}-q)\right)}$
$\displaystyle\sim\delta^{p}\left(c_{1}n_{1}+\cdots+c_{S}n_{S}-q\right)\,.$

Note that $\delta^{p}\left(c_{1}n_{1}+\cdots+c_{S}n_{S}-q\right)$ is indeed
the desired output for the modular addition task.

##  Appendix B Performance on arbitrary modular polynomials

In this appendix, we document the performance of the network constructed in
equation 10 on general modular polynomials. Note that none of the polynomials
listed in Appendix B are learnable by training 2-layer MLP or other network
architectures (such as Transformers). The MSE loss and accuracy are computed
on the entire dataset ($p^{2}$ examples). The expert-widths are taken to be
$N_{1}=500$ and $N_{2}=2000$. Since $S=3$ in all the examples, we use cubic
($\phi(x)=x^{3}$) activation function. We have truncated the MSE loss at 6
digits after the decimal point.

Modular polynomial | MSE loss | Accuracy
---|---|---
$\left(2n_{1}^{4}n_{2}+n_{1}^{2}n_{2}^{2}+3n_{1}n_{2}^{3}\right)\,mod\,97$ | 0.007674 | $100\%$
$\left(n_{1}^{5}n_{2}^{3}+4n_{1}^{2}n_{2}+5n_{1}^{2}n_{2}^{3}\right)\,mod\,97$ | 0.007660 | $100\%$
$\left(7n_{1}^{4}n_{2}^{4}+2n_{1}^{3}n_{2}^{2}+4n_{1}^{2}n_{2}^{5}\right)\,mod\,97$ | 0.007683 | $100\%$
$\left(2n_{1}^{4}n_{2}+n_{1}^{2}n_{2}^{2}+3n_{1}n_{2}^{3}\right)\,mod\,23$ | 0.009758 | $100\%$
$\left(n_{1}^{5}n_{2}^{3}+4n_{1}^{2}n_{2}+5n_{1}^{2}n_{2}^{3}\right)\,mod\,23$ | 0.009757 | $100\%$
$\left(7n_{1}^{4}n_{2}^{4}+2n_{1}^{3}n_{2}^{2}+4n_{1}^{2}n_{2}^{5}\right)\,mod\,23$ | 0.010201 | $100\%$

##  Appendix C Training on learnable and non-learnable modular polynomials

Here we present the training results on various modular polynomials. 2-layer
MLP networks with quadratic activations and width $N=5000$ are trained Adam
optimizer, MSE loss, learning rate $=0.005$, weight decay $=5.0$, on $50\%$ of
the total dataset.

We observe that training the network on polynomials of the form
$h(g_{1}(n_{1})+g_{2}(n_{2}))\,\mathrm{mod}\,p$ results in generalization.
Whereas, making a slight change to these polynomials results in inability of
the networks to generalize. This serves as evidence for Hypothesis 5.1.

Note that the following results remain qualitatively unchanged upon
changing/tuning hyperparameters.

Modular polynomial | Train loss | Test loss | Train acc | Test acc
---|---|---|---|---
$\left(4n_{1}+n_{2}^{2}\right)^{3}\,\mod\,97$ | 0.000344 | 0.000569 | $100\%$ | $100\%$
$\left(4n_{1}+n_{2}^{2}\right)^{3}+n_{1}n_{2}\,\mod\,97$ | 0.001963 | 0.011216 | $100\%$ | $2.27\%$
$\left(2n_{1}+3n_{2}\right)^{4}\,\mod\,97$ | 0.000139 | 0.000172 | $100\%$ | $100\%$
$\left(2n_{1}+3n_{2}\right)^{4}-n_{1}^{2}\,\mod\,97$ | 0.001916 | 0.011097 | $100\%$ | $3.93\%$
$\left(5n_{1}^{3}+2n_{2}^{4}\right)^{2}\,\mod\,97$ | 0.000146 | 0.000161 | $100\%$ | $100\%$
$\left(5n_{1}^{3}+2n_{2}^{4}\right)^{2}-n_{2}\,\mod\,97$ | 0.001523 | 0.006108 | $100\%$ | $72.32\%$
$\left(4n_{1}+n_{2}^{2}\right)^{3}\,\mod\,23$ | 0.000143 | 0.006274 | $100\%$ | $100\%$
$\left(4n_{1}+n_{2}^{2}\right)^{3}+n_{1}n_{2}\,\mod\,23$ | 0.000147 | 0.051866 | $100\%$ | $1.89\%$
$\left(2n_{1}+3n_{2}\right)^{4}\,\mod\,23$ | 0.000093 | 0.001010 | $100\%$ | $100\%$
$\left(2n_{1}+3n_{2}\right)^{4}-n_{1}^{2}\,\mod\,23$ | 0.000132 | 0.049618 | $100\%$ | $7.17\%$
$\left(5n_{1}^{3}+2n_{2}^{4}\right)^{2}\,\mod\,23$ | 0.000056 | 0.001004 | $100\%$ | $100\%$
$\left(5n_{1}^{3}+2n_{2}^{4}\right)^{2}-n_{2}\,\mod\,23$ | 0.000150 | 0.052047 | $100\%$ | $2.64\%$


k}^{\prime}+\frac{\lambda}{K}\boldsymbol{I}\right)^{-1}(\widehat{\boldsymbol{\alpha}}_{\sim
k}\circ\boldsymbol{Y}_{0}^{\prime}\boldsymbol{y}_{T}),$ (139)

where $\widehat{\boldsymbol{\alpha}}_{\sim
k}=\widehat{\boldsymbol{\alpha}}_{1}\circ\dots\circ\widehat{\boldsymbol{\alpha}}_{k-1}\circ\widehat{\boldsymbol{\alpha}}_{k+1}\circ\dots\circ\widehat{\boldsymbol{\alpha}}_{K}$.
Leveraging properties of the Hadamard product noted in Styan (1973), we
rewrite

$\displaystyle\boldsymbol{Y}_{0}^{\prime}\boldsymbol{Y}_{0}\circ\widehat{\boldsymbol{\alpha}}_{\sim
k}\widehat{\boldsymbol{\alpha}}_{\sim k}^{\prime}$
$\displaystyle=\boldsymbol{D}(\widehat{\boldsymbol{\alpha}}_{\sim
k})\boldsymbol{Y}_{0}^{\prime}\boldsymbol{Y}_{0}\boldsymbol{D}(\widehat{\boldsymbol{\alpha}}_{\sim
k})$ (140)
$\displaystyle\boldsymbol{Y}_{0}^{\prime}\boldsymbol{y}_{T}\circ\widehat{\boldsymbol{\alpha}}_{\sim
k}$ $\displaystyle=\boldsymbol{D}(\widehat{\boldsymbol{\alpha}}_{\sim
k})\boldsymbol{Y}_{0}^{\prime}\boldsymbol{y}_{T},$ (141)

where $\boldsymbol{D}(\widehat{\boldsymbol{\alpha}}_{\sim k})$ is the diagonal
matrix formed from $\widehat{\boldsymbol{\alpha}}_{\sim k}$. Leveraging these
equalities, we simplify (139) as

$\displaystyle\widehat{\boldsymbol{\alpha}}_{k}$
$\displaystyle=\Big{(}\boldsymbol{D}(\widehat{\boldsymbol{\alpha}}_{\sim
k})\boldsymbol{Y}_{0}^{\prime}\boldsymbol{Y}_{0}\boldsymbol{D}(\widehat{\boldsymbol{\alpha}}_{\sim
k})+\frac{\lambda}{K}\boldsymbol{I}\Big{)}^{-1}\boldsymbol{D}(\widehat{\boldsymbol{\alpha}}_{\sim
k})\boldsymbol{Y}_{0}^{\prime}\boldsymbol{y}_{T}.$ (142)

We now turn to VT regression. Following the arguments above, for every
$k\in[K]$,

$\displaystyle\widehat{\boldsymbol{\beta}}_{k}$
$\displaystyle=\Big{(}\boldsymbol{D}(\widehat{\boldsymbol{\beta}}_{\sim
k})\boldsymbol{Y}_{0}\boldsymbol{Y}_{0}^{\prime}\boldsymbol{D}(\widehat{\boldsymbol{\beta}}_{\sim
k})+\frac{\lambda}{K}\boldsymbol{I}\Big{)}^{-1}\boldsymbol{D}(\widehat{\boldsymbol{\beta}}_{\sim
k})\boldsymbol{Y}_{0}\boldsymbol{y}_{N},$ (143)

where $\widehat{\boldsymbol{\beta}}_{\sim
k}=\widehat{\boldsymbol{\beta}}_{1}\circ\dots\circ\widehat{\boldsymbol{\beta}}_{k-1}\circ\widehat{\boldsymbol{\beta}}_{k+1}\circ\dots\circ\widehat{\boldsymbol{\beta}}_{K}$
and $\boldsymbol{D}(\widehat{\boldsymbol{\beta}}_{\sim k})$ is the diagonal
matrix formed from $\widehat{\boldsymbol{\beta}}_{\sim k}$. This completes the
proof. ∎

##  Appendix D Proofs for Inference

###  D.1 Proof of Theorem 3

To establish Theorem 3, we first state a few useful results.

######  Lemma 12 (Theorem 2.7.1 of Lehmann (2000)).

Let $X_{i}$ for $i=1,\dots,n$ be independently distributed with means
$\mathbb{E}[X_{i}]=\zeta_{i}$ and variances $\sigma^{2}_{i}$, and with finite
third moments. Let $\bar{X}=(1/n)\sum_{i=1}^{n}X_{i}$. Then

$\displaystyle\frac{\bar{X}-\mathbb{E}[\bar{X}]}{\operatorname{Var}(\bar{X})^{1/2}}\xrightarrow{d}\mathcal{N}(0,1),$
(144)

provided

$\displaystyle\left(\sum_{i=1}^{n}\mathbb{E}\left[|X_{i}-\zeta_{i}|^{3}\right]\right)^{2}=o\left(\Big{(}\sum_{i=1}^{n}\sigma_{i}^{2}\Big{)}^{3}\right).$
(145)

######  Lemma 13.

Consider a random vector $\boldsymbol{x}$ and random matrix $\boldsymbol{A}$.
Let $\mathbb{E}[\boldsymbol{x}|\boldsymbol{A}]=\boldsymbol{0}$ and
$\operatorname{Cov}(\boldsymbol{x}|\boldsymbol{A})=\boldsymbol{\Sigma}$. Then
$\mathbb{E}[\boldsymbol{x}^{\prime}\boldsymbol{A}\boldsymbol{x}|\boldsymbol{A}]=\operatorname{tr}(\boldsymbol{A}\boldsymbol{\Sigma})$.

###### Proof.

(i) [HZ model] Let Assumptions 1–2 hold. By (53), Lemma 12 yields

$\displaystyle\frac{\widehat{Y}_{NT}(0)-\mathbb{E}[\widehat{Y}_{NT}(0)|\boldsymbol{y}_{N},\boldsymbol{Y}_{0}]}{\operatorname{Var}(\widehat{Y}_{NT}(0)|\boldsymbol{y}_{N},\boldsymbol{Y}_{0})^{1/2}}\xrightarrow{d}\mathcal{N}(0,1).$
(146)

To evaluate
$\mathbb{E}[\widehat{Y}_{NT}(0)|\boldsymbol{y}_{N},\boldsymbol{Y}_{0}]$, we
first observe that

$\displaystyle\mathbb{E}[\widehat{Y}_{NT}(0)|\boldsymbol{y}_{N},\boldsymbol{Y}_{0}]$
$\displaystyle=\mathbb{E}[\langle\boldsymbol{y}_{N},\widehat{\boldsymbol{\alpha}}\rangle|\boldsymbol{y}_{N},\boldsymbol{Y}_{0}]$
(147)
$\displaystyle=\mathbb{E}[\langle\boldsymbol{y}_{N},\boldsymbol{Y}_{0}^{\dagger}\boldsymbol{y}_{T}\rangle|\boldsymbol{y}_{N},\boldsymbol{Y}_{0}]$
(148)
$\displaystyle=\boldsymbol{y}^{\prime}_{N}\boldsymbol{Y}_{0}^{\dagger}\mathbb{E}[\boldsymbol{y}_{T}|\boldsymbol{y}_{N},\boldsymbol{Y}_{0}]$
(149)
$\displaystyle=\boldsymbol{y}^{\prime}_{N}\boldsymbol{Y}_{0}^{\dagger}\boldsymbol{Y}_{0}\boldsymbol{\alpha}^{*}$
(150)
$\displaystyle=\boldsymbol{y}^{\prime}_{N}\boldsymbol{H}^{v}\boldsymbol{\alpha}^{*}.$
(151)

Moving to the variance term, we note that

$\displaystyle\operatorname{Var}(\widehat{Y}_{NT}(0)|\boldsymbol{y}_{N},\boldsymbol{Y}_{0})$
$\displaystyle=\boldsymbol{y}^{\prime}_{N}\operatorname{Cov}(\widehat{\boldsymbol{\alpha}}|\boldsymbol{y}_{N},\boldsymbol{Y}_{0})\boldsymbol{y}_{N}.$
(152)

Towards evaluating the above, we note that

$\displaystyle\operatorname{Cov}(\widehat{\boldsymbol{\alpha}}|\boldsymbol{y}_{N},\boldsymbol{Y}_{0})$
$\displaystyle=\operatorname{Cov}(\boldsymbol{Y}_{0}^{\dagger}\boldsymbol{y}_{T}|\boldsymbol{y}_{N},\boldsymbol{Y}_{0})$
(153)
$\displaystyle=\boldsymbol{Y}_{0}^{\dagger}\operatorname{Cov}(\boldsymbol{y}_{T}|\boldsymbol{y}_{N},\boldsymbol{Y}_{0})(\boldsymbol{Y}^{\prime}_{0})^{\dagger}$
(154)
$\displaystyle=\boldsymbol{Y}_{0}^{\dagger}\operatorname{Cov}(\boldsymbol{\varepsilon}_{T}|\boldsymbol{y}_{N},\boldsymbol{Y}_{0})(\boldsymbol{Y}^{\prime}_{0})^{\dagger}$
(155)
$\displaystyle=\boldsymbol{Y}_{0}^{\dagger}\boldsymbol{\Sigma}^{\text{hz}}_{T}(\boldsymbol{Y}^{\prime}_{0})^{\dagger}.$
(156)

Plugging (156) into (152), we obtain

$\displaystyle\operatorname{Var}(\widehat{Y}_{NT}(0)|\boldsymbol{y}_{N},\boldsymbol{Y}_{0})=\boldsymbol{y}_{N}^{\prime}\boldsymbol{Y}_{0}^{\dagger}\boldsymbol{\Sigma}^{\text{hz}}_{T}(\boldsymbol{Y}^{\prime}_{0})^{\dagger}\boldsymbol{y}_{N}=\widehat{\boldsymbol{\beta}}^{\prime}\boldsymbol{\Sigma}^{\text{hz}}_{T}\widehat{\boldsymbol{\beta}},$
(157)

where we recall that
$\widehat{\boldsymbol{\beta}}=(\boldsymbol{Y}_{0}^{\prime})^{\dagger}\boldsymbol{y}_{N}$.
Putting it all together, we conclude

$\displaystyle\frac{\widehat{Y}_{NT}(0)-\langle\boldsymbol{y}_{N},\boldsymbol{H}^{v}\boldsymbol{\alpha}^{*}\rangle}{(\widehat{\boldsymbol{\beta}}^{\prime}\boldsymbol{\Sigma}^{\text{hz}}_{T}\widehat{\boldsymbol{\beta}})^{1/2}}\xrightarrow{d}\mathcal{N}(0,1).$
(158)

(ii) [VT model] Let Assumptions 3–4 hold. Following the arguments above, we
have

$\displaystyle\frac{\widehat{Y}_{NT}(0)-\langle\boldsymbol{y}_{T},\boldsymbol{H}^{u}\boldsymbol{\beta}^{*}\rangle}{(\widehat{\boldsymbol{\alpha}}^{\prime}\boldsymbol{\Sigma}^{\text{vt}}_{N}\widehat{\boldsymbol{\alpha}})^{1/2}}\xrightarrow{d}\mathcal{N}(0,1).$
(159)

(iii) [Mixed model] Let Assumptions 5–6 hold. We will find it useful to write

$\displaystyle\widehat{Y}_{NT}(0)$
$\displaystyle=\boldsymbol{y}^{\prime}_{N}\boldsymbol{Y}_{0}^{\dagger}\boldsymbol{y}_{T}=\sum_{i\leq
N_{0}}\sum_{t\leq T_{0}}(\boldsymbol{Y}_{0}^{\dagger})_{it}Y_{iT}Y_{Nt}.$
(160)

By Assumption 6, (160) is a sum of independent random variables with

$\displaystyle\mathbb{E}[Y_{iT}Y_{Nt}|\boldsymbol{Y}_{0}]$
$\displaystyle=\mathbb{E}[Y_{iT}|\boldsymbol{Y}_{0}]\mathbb{E}[Y_{Nt}|\boldsymbol{Y}_{0}]$
(161) $\displaystyle\operatorname{Var}(Y_{iT}Y_{Nt}|\boldsymbol{Y}_{0})$
$\displaystyle=\mathbb{E}[Y_{iT}|\boldsymbol{Y}_{0}]^{2}\sigma^{2}_{Nt}+\mathbb{E}[Y_{Nt}|\boldsymbol{Y}_{0}]^{2}\sigma^{2}_{iT}+\sigma^{2}_{iT}\sigma^{2}_{Nt}.$
(162)

Lemma 12 then establishes that

$\displaystyle\frac{\widehat{Y}_{NT}(0)-\mathbb{E}[\widehat{Y}_{NT}(0)|\boldsymbol{Y}_{0}]}{\operatorname{Var}(\widehat{Y}_{NT}(0)|\boldsymbol{Y}_{0})^{1/2}}\xrightarrow{d}\mathcal{N}(0,1).$
(163)

Our aim is to evaluate $\mathbb{E}[\widehat{Y}_{NT}(0)|\boldsymbol{Y}_{0}]$
and $\operatorname{Var}(\widehat{Y}_{NT}(0)|\boldsymbol{Y}_{0})$. Towards the
former, we use Assumptions 5–6 with the law of total expectation to obtain

$\displaystyle\mathbb{E}[\widehat{Y}_{NT}(0)|\boldsymbol{Y}_{0}]$
$\displaystyle=\mathbb{E}\left[\mathbb{E}[\langle\boldsymbol{y}_{N},\boldsymbol{Y}_{0}^{\dagger}\boldsymbol{y}_{T}\rangle|\boldsymbol{\varepsilon}_{N},\boldsymbol{Y}_{0}]|\boldsymbol{Y}_{0}\right]$
(164)
$\displaystyle=\mathbb{E}\left[\mathbb{E}[\boldsymbol{y}_{N}^{\prime}\boldsymbol{Y}_{0}^{\dagger}(\boldsymbol{Y}_{0}\boldsymbol{\alpha}^{*}+\boldsymbol{\varepsilon}_{T})|\boldsymbol{\varepsilon}_{N},\boldsymbol{Y}_{0}]|\boldsymbol{Y}_{0}\right]$
(165)
$\displaystyle=\mathbb{E}\left[(\boldsymbol{Y}^{\prime}_{0}\boldsymbol{\beta}^{*}+\boldsymbol{\varepsilon}_{N})^{\prime}\boldsymbol{Y}_{0}^{\dagger}\boldsymbol{Y}_{0}\boldsymbol{\alpha}^{*}|\boldsymbol{Y}_{0}\right]$
(166)
$\displaystyle=\langle\boldsymbol{\beta}^{*},\boldsymbol{Y}_{0}\boldsymbol{\alpha}^{*}\rangle.$
(167)

Note that we have used the fact that $\boldsymbol{y}_{N}$ is deterministic
given $(\boldsymbol{\varepsilon}_{N},\boldsymbol{Y}_{0})$. Similarly, by the
law of total variance,

$\displaystyle\operatorname{Var}(\widehat{Y}_{NT}(0)|\boldsymbol{Y}_{0})$
$\displaystyle=\mathbb{E}[\operatorname{Var}(\widehat{Y}_{NT}(0)|\boldsymbol{\varepsilon}_{N},\boldsymbol{Y}_{0})|\boldsymbol{Y}_{0}]+\operatorname{Var}(\mathbb{E}[\widehat{Y}_{NT}(0)|\boldsymbol{\varepsilon}_{N},\boldsymbol{Y}_{0}]|\boldsymbol{Y}_{0}).$
(168)

Following the derivation of (157), we have

$\displaystyle\mathbb{E}[\operatorname{Var}(\widehat{Y}_{NT}(0)|\boldsymbol{\varepsilon}_{N},\boldsymbol{Y}_{0})|\boldsymbol{Y}_{0}]$
(169)
$\displaystyle\qquad\qquad=\mathbb{E}[\boldsymbol{y}^{\prime}_{N}\boldsymbol{Y}_{0}^{\dagger}\boldsymbol{\Sigma}^{\text{mix}}_{T}(\boldsymbol{Y}^{\prime}_{0})^{\dagger}\boldsymbol{y}_{N}|\boldsymbol{Y}_{0}]$
(170)
$\displaystyle\qquad\qquad=(\boldsymbol{Y}^{\prime}_{0}\boldsymbol{\beta}^{*})^{\prime}\boldsymbol{A}(\boldsymbol{Y}^{\prime}_{0}\boldsymbol{\beta}^{*})+\mathbb{E}[\boldsymbol{\varepsilon}^{\prime}_{N}\boldsymbol{A}\boldsymbol{\varepsilon}_{N}|\boldsymbol{Y}_{0}]+2\mathbb{E}[\boldsymbol{\varepsilon}^{\prime}_{N}\boldsymbol{Y}^{\prime}_{0}\boldsymbol{\beta}^{*}|\boldsymbol{Y}_{0}],$
(171)

where
$\boldsymbol{A}=\boldsymbol{Y}^{\dagger}_{0}\boldsymbol{\Sigma}^{\text{mix}}_{T}(\boldsymbol{Y}^{\prime}_{0})^{\dagger}$.
Notice that Assumption 6 gives
$\mathbb{E}[\boldsymbol{\varepsilon}^{\prime}_{N}\boldsymbol{Y}^{\prime}_{0}\boldsymbol{\beta}^{*}|\boldsymbol{Y}_{0}]=0$.
Since $\boldsymbol{A}$ is deterministic given $\boldsymbol{Y}_{0}$, Lemma 13
yields

$\displaystyle\mathbb{E}[\boldsymbol{\varepsilon}^{\prime}_{N}\boldsymbol{A}\boldsymbol{\varepsilon}_{N}|\boldsymbol{Y}_{0}]$
$\displaystyle=\operatorname{tr}(\boldsymbol{A}\boldsymbol{\Sigma}^{\text{mix}}_{N}).$
(172)

Following the arguments that led to the derivation of (151), we have

$\displaystyle\operatorname{Var}(\mathbb{E}[\widehat{Y}_{NT}(0)|\boldsymbol{\varepsilon}_{N},\boldsymbol{Y}_{0}]|\boldsymbol{Y}_{0})$
$\displaystyle=\operatorname{Var}(\boldsymbol{y}^{\prime}_{N}\boldsymbol{H}^{v}\boldsymbol{\alpha}^{*}|\boldsymbol{Y}_{0})=(\boldsymbol{H}^{v}\boldsymbol{\alpha}^{*})^{\prime}\boldsymbol{\Sigma}^{\text{mix}}_{N}(\boldsymbol{H}^{v}\boldsymbol{\alpha}^{*}).$
(173)

Plugging (171), (172), and (173) into (168), we arrive at

$\displaystyle\operatorname{Var}(\widehat{Y}_{NT}(0)|\boldsymbol{Y}_{0})$
(174)
$\displaystyle\qquad=(\boldsymbol{H}^{v}\boldsymbol{\alpha}^{*})^{\prime}\boldsymbol{\Sigma}^{\text{mix}}_{N}(\boldsymbol{H}^{v}\boldsymbol{\alpha}^{*})+(\boldsymbol{H}^{u}\boldsymbol{\beta}^{*})^{\prime}\boldsymbol{\Sigma}^{\text{mix}}_{T}(\boldsymbol{H}^{u}\boldsymbol{\beta}^{*})+\operatorname{tr}(\boldsymbol{Y}^{\dagger}_{0}\boldsymbol{\Sigma}^{\text{mix}}_{T}(\boldsymbol{Y}^{\prime}_{0})^{\dagger}\boldsymbol{\Sigma}^{\text{mix}}_{N}).$
(175)

This completes the proof. ∎

###  D.2 Proofs for Model-Based Confidence Intervals

We first state a useful lemma to prove Lemmas 1–3.

######  Lemma 14.

_[Mixed model]_ Let Assumptions 5–6 hold. Then,

$\displaystyle\mathbb{E}[\widehat{v}^{\emph{mix}}_{0}|\boldsymbol{Y}_{0}]$
$\displaystyle=(\boldsymbol{H}^{u}\boldsymbol{\beta}^{*})^{\prime}\mathbb{E}[\widehat{\boldsymbol{\Sigma}}_{T}|\boldsymbol{Y}_{0}](\boldsymbol{H}^{u}\boldsymbol{\beta}^{*})+\operatorname{tr}(\boldsymbol{Y}_{0}^{\dagger}\mathbb{E}[\widehat{\boldsymbol{\Sigma}}_{T}|\boldsymbol{Y}_{0}](\boldsymbol{Y}^{\prime}_{0})^{\dagger}\boldsymbol{\Sigma}^{\emph{mix}}_{N})$
(176)
$\displaystyle\quad+(\boldsymbol{H}^{v}\boldsymbol{\alpha}^{*})^{\prime}\mathbb{E}[\widehat{\boldsymbol{\Sigma}}_{N}|\boldsymbol{Y}_{0}](\boldsymbol{H}^{v}\boldsymbol{\alpha}^{*})+\operatorname{tr}(\boldsymbol{Y}_{0}^{\dagger}\boldsymbol{\Sigma}^{\emph{mix}}_{T}(\boldsymbol{Y}^{\prime}_{0})^{\dagger}\mathbb{E}[\widehat{\boldsymbol{\Sigma}}_{N}|\boldsymbol{Y}_{0}])$
(177)
$\displaystyle\quad-\operatorname{tr}(\boldsymbol{Y}_{0}^{\dagger}\mathbb{E}[\widehat{\boldsymbol{\Sigma}}_{T}|\boldsymbol{Y}_{0}](\boldsymbol{Y}^{\prime}_{0})^{\dagger}\mathbb{E}[\widehat{\boldsymbol{\Sigma}}_{N}|\boldsymbol{Y}_{0}]).$
(178)

####  D.2.1 Proof of Lemma 1

###### Proof.

(i) [HZ model] Let Assumptions 1–2 hold. Taking note that
$\boldsymbol{H}^{u}_{\perp}\boldsymbol{Y}_{0}=\boldsymbol{0}$,

$\displaystyle\|\boldsymbol{H}^{u}_{\perp}\boldsymbol{y}_{T}\|_{2}^{2}$
$\displaystyle=\boldsymbol{y}^{\prime}_{T}\boldsymbol{H}^{u}_{\perp}\boldsymbol{y}_{T}$
(179)
$\displaystyle=(\boldsymbol{Y}_{0}\boldsymbol{\alpha}+\boldsymbol{\varepsilon}_{T})^{\prime}\boldsymbol{H}^{u}_{\perp}(\boldsymbol{Y}_{0}\boldsymbol{\alpha}+\boldsymbol{\varepsilon}_{T})$
(180)
$\displaystyle=\boldsymbol{\varepsilon}^{\prime}_{T}\boldsymbol{H}^{u}_{\perp}\boldsymbol{\varepsilon}_{T}.$
(181)

Applying Lemma 13 then gives

$\displaystyle\mathbb{E}[\boldsymbol{\varepsilon}^{\prime}_{T}\boldsymbol{H}^{u}_{\perp}\boldsymbol{\varepsilon}_{T}|\boldsymbol{y}_{N},\boldsymbol{Y}_{0}]=\operatorname{tr}(\boldsymbol{H}^{u}_{\perp})(\sigma^{\text{hz}}_{T})^{2}=(N_{0}-R)(\sigma^{\text{hz}}_{T})^{2},$
(182)

where the final equality follows because the trace of a projection matrix
equals its rank. Taken altogether, we have
$\mathbb{E}[\widehat{\boldsymbol{\Sigma}}_{T}^{\text{homo}}|\boldsymbol{y}_{N},\boldsymbol{Y}_{0}]=\boldsymbol{\Sigma}^{\text{hz}}_{T}$.
Therefore,

$\displaystyle\mathbb{E}[\widehat{v}_{0}^{\text{hz},\text{homo}}|\boldsymbol{y}_{N},\boldsymbol{Y}_{0}]$
$\displaystyle=\widehat{\boldsymbol{\beta}}^{\prime}\mathbb{E}[\widehat{\boldsymbol{\Sigma}}^{\text{homo}}_{T}|\boldsymbol{y}_{N},\boldsymbol{Y}_{0}]\widehat{\boldsymbol{\beta}}=v_{0}^{\text{hz}}.$
(183)

(ii) [VT model] Let Assumptions 3–4 hold. Following the arguments above, we
conclude that
$\mathbb{E}[\widehat{\boldsymbol{\Sigma}}_{N}^{\text{homo}}|\boldsymbol{y}_{T},\boldsymbol{Y}_{0}]=\boldsymbol{\Sigma}^{\text{vt}}_{N}$
and
$\mathbb{E}[\widehat{v}_{0}^{\text{vt},\text{homo}}|\boldsymbol{y}_{T},\boldsymbol{Y}_{0}]=v_{0}^{\text{vt}}$.

(iii) [Mixed model] Let Assumptions 5–6 hold. Following the arguments that led
to (182), we obtain
$\mathbb{E}[\widehat{\boldsymbol{\Sigma}}^{\text{homo}}_{T}|\boldsymbol{Y}_{0}]=\boldsymbol{\Sigma}^{\text{mix}}_{T}$
and
$\mathbb{E}[\widehat{\boldsymbol{\Sigma}}^{\text{homo}}_{N}|\boldsymbol{Y}_{0}]=\boldsymbol{\Sigma}^{\text{mix}}_{N}$.
Applying Lemma 14 then gives
$\mathbb{E}[\widehat{v}^{\text{mix},\text{homo}}_{0}|\boldsymbol{Y}_{0}]=v^{\text{mix}}_{0}$.
The proof is complete. ∎

###  D.3 Proof of Lemma 2

###### Proof.

Before we establish the biases of
$(\widehat{\boldsymbol{\Sigma}}^{\text{jack}}_{T},\widehat{\boldsymbol{\Sigma}}^{\text{jack}}_{N})$,
we first justify their forms. As noted in Section 4.1.3, jackknife is a
popular approach to estimate the covariances of
$(\widehat{\boldsymbol{\alpha}},\widehat{\boldsymbol{\beta}})$. Below, we
follow the standard techniques to derive the jackknife estimate of these
objects, which will then be used to derive
$(\widehat{\boldsymbol{\Sigma}}^{\text{jack}}_{T},\widehat{\boldsymbol{\Sigma}}^{\text{jack}}_{N})$.
Without loss of generality, we begin with $\widehat{\boldsymbol{\alpha}}$.
Notably, while standard derivations consider $\boldsymbol{Y}_{0}$ with full
column rank, we consider a general matrix $\boldsymbol{Y}_{0}$ that may be
rank deficient. This difference is subtle so the following proof is by no
means novel. We provide it simply for completeness.

To describe the jackknife, we define $\widehat{\boldsymbol{\alpha}}_{\sim i}$
as the minimum $\ell_{2}$-norm solution to (1), where
$\lambda_{1}=\lambda_{2}=0$, without the $i$th observation, i.e.,

$\displaystyle\widehat{\boldsymbol{\alpha}}_{\sim i}$
$\displaystyle=(\boldsymbol{Y}^{\prime}_{0,\sim i}\boldsymbol{Y}_{0,\sim
i})^{\dagger}\boldsymbol{Y}^{\prime}_{0,\sim i}\boldsymbol{y}_{T,\sim i},$
(184)

where $\boldsymbol{Y}_{0,\sim i}$ and $\boldsymbol{y}_{T,\sim i}$ correspond
to $\boldsymbol{Y}_{0}$ and $\boldsymbol{y}_{T}$ without the $i$th
observation. We define the pseudo-estimator as
$\tilde{\boldsymbol{\alpha}}_{i}=T_{0}\widehat{\boldsymbol{\alpha}}-(T_{0}-1)\widehat{\boldsymbol{\alpha}}_{\sim
i}$. With these quantities defined, we write the jackknife variance estimator
as

$\displaystyle\widehat{\boldsymbol{V}}^{\text{jack}}=\frac{1}{(T_{0}-1)^{2}}\sum_{i\leq
N_{0}}(\tilde{\boldsymbol{\alpha}}_{i}-\widehat{\boldsymbol{\alpha}})(\tilde{\boldsymbol{\alpha}}_{i}-\widehat{\boldsymbol{\alpha}})^{\prime}.$
(185)

To evaluate this quantity, we will rewrite
$\widehat{\boldsymbol{\alpha}}_{\sim i}$ in a more convenient form. In
particular,

$\displaystyle\boldsymbol{Y}^{\prime}_{0,\sim i}\boldsymbol{Y}_{0,\sim
i}=\boldsymbol{Y}^{\prime}_{0}\boldsymbol{Y}^{\prime}_{0}-\boldsymbol{y}_{i}\boldsymbol{y}^{\prime}_{i}$
(186) $\displaystyle\boldsymbol{Y}^{\prime}_{0,\sim i}\boldsymbol{y}_{T,\sim
i}=\boldsymbol{Y}^{\prime}_{0}\boldsymbol{y}_{T}-\boldsymbol{y}_{i}Y_{iT},$
(187)

where $\boldsymbol{y}_{i}=[Y_{it}:t\leq T_{0}]$ is the $i$th row of
$\boldsymbol{Y}_{0}$. We do not assume that
$\boldsymbol{Y}^{\prime}_{0}\boldsymbol{Y}_{0}$ is nonsingular. As such, we
use a generalized form of the Sherman-Morrison formula (Cline, 1965; Meyer,
1973) to obtain

$\displaystyle(\boldsymbol{Y}^{\prime}_{0,\sim i}\boldsymbol{Y}_{0,\sim
i})^{\dagger}=(\boldsymbol{Y}^{\prime}_{0}\boldsymbol{Y}_{0})^{\dagger}+(1-H^{u}_{ii})^{-1}(\boldsymbol{Y}^{\prime}_{0}\boldsymbol{Y}_{0})^{\dagger}\boldsymbol{y}_{i}\boldsymbol{y}^{\prime}_{i}(\boldsymbol{Y}^{\prime}_{0}\boldsymbol{Y}_{0})^{\dagger}.$
(188)

Recall
$\widehat{\boldsymbol{\alpha}}=(\boldsymbol{Y}^{\prime}_{0}\boldsymbol{Y}_{0})^{\dagger}\boldsymbol{Y}^{\prime}_{0}\boldsymbol{y}_{T}$
and note $Y_{iT}-\boldsymbol{y}^{\prime}_{i}\widehat{\boldsymbol{\alpha}}$ is
the $i$th element of
$\widehat{\boldsymbol{\varepsilon}}_{T}=\boldsymbol{H}^{u}_{\perp}\boldsymbol{y}_{T}$.
Using these facts, we plug (188) into (184) to yield

$\displaystyle\widehat{\boldsymbol{\alpha}}_{\sim i}$
$\displaystyle=\left[(\boldsymbol{Y}^{\prime}_{0}\boldsymbol{Y}_{0})^{\dagger}+(1-H^{u}_{ii})^{-1}(\boldsymbol{Y}^{\prime}_{0}\boldsymbol{Y}_{0})^{\dagger}\boldsymbol{y}_{i}\boldsymbol{y}^{\prime}_{i}(\boldsymbol{Y}^{\prime}_{0}\boldsymbol{Y}_{0})^{\dagger}\right](\boldsymbol{Y}^{\prime}_{0}\boldsymbol{y}_{T}-\boldsymbol{y}_{i}Y_{iT})$
(189)
$\displaystyle=\widehat{\boldsymbol{\alpha}}-(\boldsymbol{Y}^{\prime}_{0}\boldsymbol{Y}_{0})^{\dagger}\boldsymbol{y}_{i}Y_{iT}+(1-H^{u}_{ii})^{-1}(\boldsymbol{Y}^{\prime}_{0}\boldsymbol{Y}_{0})^{\dagger}\boldsymbol{y}_{i}\boldsymbol{y}^{\prime}_{i}\widehat{\boldsymbol{\alpha}}-H^{u}_{ii}(1-H^{u}_{ii})^{-1}(\boldsymbol{Y}^{\prime}_{0}\boldsymbol{Y}_{0})^{\dagger}\boldsymbol{y}_{i}Y_{iT}$
(190)
$\displaystyle=\widehat{\boldsymbol{\alpha}}-(1-H^{u}_{ii})^{-1}(\boldsymbol{Y}^{\prime}_{0}\boldsymbol{Y}_{0})^{\dagger}\boldsymbol{y}_{i}\widehat{\varepsilon}_{iT}.$
(191)

Inserting (191) into our pseudo-estimate, we have

$\displaystyle\tilde{\boldsymbol{\alpha}}_{i}$
$\displaystyle=T_{0}\widehat{\boldsymbol{\alpha}}-(T_{0}-1)\left(\widehat{\boldsymbol{\alpha}}-(1-H^{u}_{ii})^{-1}(\boldsymbol{Y}^{\prime}_{0}\boldsymbol{Y}_{0})^{\dagger}\boldsymbol{y}_{i}\widehat{\varepsilon}_{iT}\right)$
(192)
$\displaystyle=\widehat{\boldsymbol{\alpha}}+(T_{0}-1)(1-H^{u}_{ii})^{-1}(\boldsymbol{Y}^{\prime}_{0}\boldsymbol{Y}_{0})^{\dagger}\boldsymbol{y}_{i}\widehat{\varepsilon}_{iT}.$
(193)

Inserting (193) into (185), we have

$\displaystyle\widehat{\boldsymbol{V}}^{\text{jack}}$
$\displaystyle=(\boldsymbol{Y}^{\prime}_{0}\boldsymbol{Y}_{0})^{\dagger}\left(\sum_{i\leq
N_{0}}\frac{\widehat{\varepsilon}_{iT}^{2}}{(1-H^{u}_{ii})^{2}}\boldsymbol{y}_{i}\boldsymbol{y}^{\prime}_{i}\right)(\boldsymbol{Y}^{\prime}_{0}\boldsymbol{Y}_{0})^{\dagger}$
(194)
$\displaystyle=(\boldsymbol{Y}^{\prime}_{0}\boldsymbol{Y}_{0})^{\dagger}\boldsymbol{Y}^{\prime}_{0}\boldsymbol{\Omega}\boldsymbol{Y}_{0}(\boldsymbol{Y}^{\prime}_{0}\boldsymbol{Y}_{0})^{\dagger},$
(195)

where $\boldsymbol{\Omega}$ is a diagonal matrix with
$\Omega_{ii}=\widehat{\varepsilon}_{iT}^{2}(1-H^{u}_{ii})^{-2}$. Equivalently,
$\boldsymbol{\Omega}=\text{diag}([\boldsymbol{H}^{u}_{\perp}\circ\boldsymbol{H}^{u}_{\perp}\circ\boldsymbol{I}]^{\dagger}[\widehat{\boldsymbol{\varepsilon}}_{T}\circ\widehat{\boldsymbol{\varepsilon}}_{T}])$.
It then follows that

$\displaystyle\boldsymbol{y}^{\prime}_{N}\widehat{\boldsymbol{V}}^{\text{jack}}\boldsymbol{y}_{N}$
$\displaystyle=\widehat{\boldsymbol{\beta}}^{\prime}\boldsymbol{\Omega}\widehat{\boldsymbol{\beta}}.$
(196)

To arrive at (38), we define
$\widehat{\boldsymbol{\Sigma}}^{\text{jack}}_{T}=\boldsymbol{\Omega}$. This
corresponds to the EHW estimator with the jackknife correction. We derive (39)
for $\widehat{\boldsymbol{\beta}}$ by applying the same arguments above. Now,
we will evaluate the biases of
$(\widehat{\boldsymbol{\Sigma}}^{\text{jack}}_{T},\widehat{\boldsymbol{\Sigma}}^{\text{jack}}_{N})$.

(i) [HZ model] Let Assumptions 1–2 hold. We define
$(\sigma^{\text{hz}}_{iT})^{2}=\operatorname{Var}(\varepsilon_{iT}|\boldsymbol{y}_{N},\boldsymbol{Y}_{0})$
for $i=1,\dots,N_{0}$. Observe that

$\displaystyle\mathbb{E}[(\boldsymbol{H}^{u}_{\perp}\circ\boldsymbol{H}^{u}_{\perp}\circ\boldsymbol{I})^{\dagger}(\widehat{\boldsymbol{\varepsilon}}_{T}\circ\widehat{\boldsymbol{\varepsilon}}_{T})|\boldsymbol{y}_{N},\boldsymbol{Y}_{0}]$
$\displaystyle=(\boldsymbol{H}^{u}_{\perp}\circ\boldsymbol{H}^{u}_{\perp}\circ\boldsymbol{I})^{\dagger}\mathbb{E}[\widehat{\boldsymbol{\varepsilon}}_{T}\circ\widehat{\boldsymbol{\varepsilon}}_{T}|\boldsymbol{y}_{N},\boldsymbol{Y}_{0}].$
(197)

To evaluate (197), we follow the derivations of (151) and (156) to obtain

$\displaystyle\mathbb{E}[\widehat{\boldsymbol{\varepsilon}}_{T}|\boldsymbol{y}_{N},\boldsymbol{Y}_{0}]$
$\displaystyle=\boldsymbol{H}^{u}_{\perp}\boldsymbol{Y}_{0}\boldsymbol{\alpha}^{*}=\boldsymbol{0}$
(198)
$\displaystyle\operatorname{Cov}(\widehat{\boldsymbol{\varepsilon}}_{T}|\boldsymbol{y}_{N},\boldsymbol{Y}_{0})$
$\displaystyle=\boldsymbol{H}^{u}_{\perp}\boldsymbol{\Sigma}^{\text{hz}}_{T}\boldsymbol{H}^{u}_{\perp}.$
(199)

Recall that $\mathbb{E}[X^{2}]=\operatorname{Var}(X)+\mathbb{E}[X]^{2}$ for
any random variable $X$. Thus, combining (198) with (199) gives

$\displaystyle\mathbb{E}[\widehat{\boldsymbol{\varepsilon}}_{T}\circ\widehat{\boldsymbol{\varepsilon}}_{T}|\boldsymbol{y}_{N},\boldsymbol{Y}_{0}]$
$\displaystyle=(\boldsymbol{H}^{u}_{\perp}\boldsymbol{\Sigma}^{\text{hz}}_{T}\boldsymbol{H}^{u}_{\perp}\circ\boldsymbol{I})\boldsymbol{1}.$
(200)

Let
$\widehat{\boldsymbol{\gamma}}=\mathbb{E}[\widehat{\boldsymbol{\varepsilon}}_{T}\circ\widehat{\boldsymbol{\varepsilon}}_{T}|\boldsymbol{y}_{N},\boldsymbol{Y}_{0}]$.
By (200), the $\ell$th entry of $\widehat{\boldsymbol{\gamma}}$ can be written
as

$\displaystyle\widehat{\gamma}_{\ell}=\sum_{j\neq\ell}(H^{u}_{j\ell})^{2}(\sigma^{\text{hz}}_{jT})^{2}+(1-H^{u}_{\ell\ell})^{2}(\sigma^{\text{hz}}_{\ell
T})^{2},$ (201)

where $H^{u}_{j\ell}$ is the $(j,\ell)$th entry of $\boldsymbol{H}^{u}$. In
turn, this allows us to rewrite (200) as

$\displaystyle\widehat{\boldsymbol{\gamma}}$
$\displaystyle=(\boldsymbol{H}^{u}_{\perp}\circ\boldsymbol{H}^{u}_{\perp})\boldsymbol{\Sigma}^{\text{hz}}_{T}\boldsymbol{1}.$
(202)

Next, let
$\widehat{\boldsymbol{\zeta}}=(\boldsymbol{H}^{u}_{\perp}\circ\boldsymbol{H}^{u}_{\perp}\circ\boldsymbol{I})^{-1}\widehat{\boldsymbol{\gamma}}$.
Notice that the $\ell$th entry of $\widehat{\boldsymbol{\zeta}}$ is given by

$\displaystyle\widehat{\zeta}_{\ell}$ $\displaystyle=(\sigma^{\text{hz}}_{\ell
T})^{2}+\sum_{j\neq\ell}\frac{(H^{u}_{\ell
j})^{2}}{(1-H^{u}_{\ell\ell})^{2}}(\sigma^{\text{hz}}_{jT})^{2}.$ (203)

Therefore,
$\text{diag}(\widehat{\boldsymbol{\zeta}})=\boldsymbol{\Sigma}^{\text{hz}}_{T}+\boldsymbol{\Delta}^{\text{hz}}$,
where
$\Delta^{\text{hz}}_{\ell\ell}=\sum_{j\neq\ell}(\sigma^{\text{hz}}_{jT})^{2}(H^{u}_{\ell
j})^{2}(1-H^{u}_{\ell\ell})^{-2}$ for $\ell=1,\dots,N_{0}$. Notice if
$\max_{\ell}H^{u}_{\ell\ell}<1$, then
$(\boldsymbol{H}^{u}_{\perp}\circ\boldsymbol{H}^{u}_{\perp}\circ\boldsymbol{I})$
is nonsingular, i.e., the pseudo-inverse is precisely the inverse. In this
situation, plugging the above into (197) gives

$\displaystyle\mathbb{E}[\widehat{\boldsymbol{\Sigma}}^{\text{jack}}_{T}|\boldsymbol{y}_{N},\boldsymbol{Y}_{0}]$
$\displaystyle=\text{diag}\left((\boldsymbol{H}^{u}_{\perp}\circ\boldsymbol{H}^{u}_{\perp}\circ\boldsymbol{I})^{-1}\mathbb{E}[\widehat{\boldsymbol{\varepsilon}}_{T}\circ\widehat{\boldsymbol{\varepsilon}}_{T}|\boldsymbol{y}_{N},\boldsymbol{Y}_{0}]\right)$
(204)
$\displaystyle=\text{diag}\left((\boldsymbol{H}^{u}_{\perp}\circ\boldsymbol{H}^{u}_{\perp}\circ\boldsymbol{I})^{-1}\widehat{\boldsymbol{\gamma}}\right)$
(205) $\displaystyle=\text{diag}(\widehat{\boldsymbol{\zeta}})$ (206)
$\displaystyle=\boldsymbol{\Sigma}^{\text{hz}}_{T}+\boldsymbol{\Delta}^{\text{hz}}.$
(207)

From this, we conclude that

$\displaystyle\mathbb{E}[\widehat{v}_{0}^{\text{hz},\text{jack}}|\boldsymbol{y}_{N},\boldsymbol{Y}_{0}]$
$\displaystyle=\widehat{\boldsymbol{\beta}}^{\prime}\mathbb{E}[\widehat{\boldsymbol{\Sigma}}_{T}^{\text{jack}}|\boldsymbol{y}_{N},\boldsymbol{Y}_{0}]\widehat{\boldsymbol{\beta}}$
(208)
$\displaystyle=\widehat{\boldsymbol{\beta}}^{\prime}(\boldsymbol{\Sigma}^{\text{hz}}_{T}+\boldsymbol{\Delta}^{\text{hz}})\widehat{\boldsymbol{\beta}}$
(209)
$\displaystyle=v_{0}^{\text{hz}}+\widehat{\boldsymbol{\beta}}^{\prime}\boldsymbol{\Delta}^{\text{hz}}\widehat{\boldsymbol{\beta}},$
(210)

where we note that
$\widehat{\boldsymbol{\beta}}^{\prime}\boldsymbol{\Delta}^{\text{hz}}\widehat{\boldsymbol{\beta}}\geq
0$.

(ii) [VT model] Let Assumptions 3–4 hold. Following the arguments above, we
conclude
$\mathbb{E}[\widehat{\boldsymbol{\Sigma}}^{\text{jack}}_{N}|\boldsymbol{y}_{T},\boldsymbol{Y}_{0}]=\boldsymbol{\Sigma}^{\text{vt}}_{N}+\boldsymbol{\Gamma}^{\text{vt}}$,
where
$\Gamma^{\text{vt}}_{\ell\ell}=\sum_{j\neq\ell}(\sigma^{\text{vt}}_{Nj})^{2}(H^{v}_{\ell
j})^{2}(1-H^{v}_{\ell\ell})^{-2}$ for $\ell=1,\dots,T_{0}$. Thus,
$\mathbb{E}[\widehat{v}_{0}^{\text{vt},\text{jack}}|\boldsymbol{y}_{T},\boldsymbol{Y}_{0}]=v^{\text{vt}}_{0}+\widehat{\boldsymbol{\alpha}}^{\prime}\boldsymbol{\Gamma}^{\text{vt}}\widehat{\boldsymbol{\alpha}}$,
where we note that
$\widehat{\boldsymbol{\alpha}}^{\prime}\boldsymbol{\Gamma}^{\text{vt}}\widehat{\boldsymbol{\alpha}}\geq
0$.

(ii) [Mixed model] Let Assumptions 5–6 hold. We define
$(\sigma^{\text{mix}}_{iT})^{2}=\operatorname{Var}(\varepsilon_{iT}|\boldsymbol{Y}_{0})$
for $i=1,\dots,N_{0}$ and
$(\sigma^{\text{mix}}_{Nt})^{2}=\operatorname{Var}(\varepsilon_{Nt}|\boldsymbol{Y}_{0})$
for $t=1,\dots,T_{0}$. Following the arguments that led to (207), we obtain
$\mathbb{E}[\widehat{\boldsymbol{\Sigma}}^{\text{jack}}_{T}|\boldsymbol{Y}_{0}]=\boldsymbol{\Sigma}^{\text{mix}}_{T}+\boldsymbol{\Delta}^{\text{mix}}$,
where
$\Delta^{\text{mix}}_{\ell\ell}=\sum_{j\neq\ell}(\sigma^{\text{mix}}_{jT})^{2}(H^{u}_{\ell
j})^{2}(1-H^{u}_{\ell\ell})^{-2}$ for $\ell=1,\dots,N_{0}$. Similarly, we
obtain
$\mathbb{E}[\widehat{\boldsymbol{\Sigma}}^{\text{jack}}_{N}|\boldsymbol{Y}_{0}]=\boldsymbol{\Sigma}^{\text{mix}}_{N}+\boldsymbol{\Gamma}^{\text{mix}}$,
where
$\Gamma^{\text{mix}}_{\ell\ell}=\sum_{j\neq\ell}(\sigma^{\text{mix}}_{Nj})^{2}(H^{v}_{\ell
j})^{2}(1-H^{v}_{\ell\ell})^{-2}$ for $\ell=1,\dots,T_{0}$. Applying Lemma 14
then gives

$\displaystyle\mathbb{E}[\widehat{v}^{\text{mix},\text{jack}}_{0}|\boldsymbol{Y}_{0}]$
(211)
$\displaystyle=v^{\text{mix}}_{0}+(\boldsymbol{H}^{u}\boldsymbol{\beta}^{*})^{\prime}\boldsymbol{\Delta}^{\text{mix}}(\boldsymbol{H}^{u}\boldsymbol{\beta}^{*})+(\boldsymbol{H}^{v}\boldsymbol{\alpha}^{*})^{\prime}\boldsymbol{\Gamma}^{\text{mix}}(\boldsymbol{H}^{v}\boldsymbol{\alpha}^{*})+\operatorname{tr}(\boldsymbol{Y}_{0}^{\dagger}\boldsymbol{\Delta}^{\text{mix}}(\boldsymbol{Y}^{\prime}_{0})^{\dagger}\boldsymbol{\Gamma}^{\text{mix}}).$
(212)

The proof is complete. ∎

###  D.4 Proof of Lemma 3

###### Proof.

We adopt the strategy of Hartley et al. (1969) to prove our desired result.

(ii) [HZ model] Let Assumptions 1–2 hold. As in the proof of Lemma 2, we
define
$\widehat{\boldsymbol{\varepsilon}}_{T}=\boldsymbol{H}^{u}_{\perp}\boldsymbol{y}_{T}$.
Observe

$\displaystyle\mathbb{E}[(\boldsymbol{H}^{u}_{\perp}\circ\boldsymbol{H}^{u}_{\perp})^{-1}(\widehat{\boldsymbol{\varepsilon}}_{T}\circ\widehat{\boldsymbol{\varepsilon}}_{T})|\boldsymbol{y}_{N},\boldsymbol{Y}_{0}]$
$\displaystyle=(\boldsymbol{H}^{u}_{\perp}\circ\boldsymbol{H}^{u}_{\perp})^{-1}\mathbb{E}[\widehat{\boldsymbol{\varepsilon}}_{T}\circ\widehat{\boldsymbol{\varepsilon}}_{T}|\boldsymbol{y}_{N},\boldsymbol{Y}_{0}].$
(213)

To evaluate (213), we plug in (202) to obtain

$\displaystyle\mathbb{E}[(\boldsymbol{H}^{u}_{\perp}\circ\boldsymbol{H}^{u}_{\perp})^{-1}(\widehat{\boldsymbol{\varepsilon}}_{T}\circ\widehat{\boldsymbol{\varepsilon}}_{T})|\boldsymbol{y}_{N},\boldsymbol{Y}_{0}]$
$\displaystyle=(\boldsymbol{H}^{u}_{\perp}\circ\boldsymbol{H}^{u}_{\perp})^{-1}(\boldsymbol{H}^{u}_{\perp}\circ\boldsymbol{H}^{u}_{\perp})\boldsymbol{\Sigma}^{\text{hz}}_{T}\boldsymbol{1}=\boldsymbol{\Sigma}^{\text{hz}}_{T}\boldsymbol{1}.$
(214)

Plugging (214) into (213) yields

$\displaystyle\mathbb{E}[\widehat{\boldsymbol{\Sigma}}^{\text{HRK}}_{T}|\boldsymbol{y}_{N},\boldsymbol{Y}_{0}]$
$\displaystyle=\text{diag}\left((\boldsymbol{H}^{u}_{\perp}\circ\boldsymbol{H}^{u}_{\perp})^{-1}\mathbb{E}[\widehat{\boldsymbol{\varepsilon}}_{T}\circ\widehat{\boldsymbol{\varepsilon}}_{T}|\boldsymbol{y}_{N},\boldsymbol{Y}_{0}]\right)=\boldsymbol{\Sigma}^{\text{hz}}_{T}.$
(215)

It then follows that
$\mathbb{E}[\widehat{v}_{0}^{\text{hz},\text{HRK}}|\boldsymbol{y}_{N},\boldsymbol{Y}_{0}]=v_{0}^{\text{hz}}$.

(ii) [VT model] Let Assumptions 3–4 hold. Following the same arguments as
above, we conclude
$\mathbb{E}[\widehat{\boldsymbol{\Sigma}}^{\text{HRK}}_{N}|\boldsymbol{y}_{T},\boldsymbol{Y}_{0}]=\boldsymbol{\Sigma}^{\text{vt}}_{N}$
and
$\mathbb{E}[\widehat{v}_{0}^{\text{vt},\text{HRK}}|\boldsymbol{y}_{T},\boldsymbol{Y}_{0}]=v^{\text{vt}}_{0}$.

(ii) [Mixed model] Let Assumptions 5–6 hold. Following the arguments that led
to (215), we obtain
$\mathbb{E}[\widehat{\boldsymbol{\Sigma}}^{\text{HRK}}_{T}|\boldsymbol{Y}_{0}]=\boldsymbol{\Sigma}^{\text{mix}}_{T}$
and
$\mathbb{E}[\widehat{\boldsymbol{\Sigma}}^{\text{HRK}}_{N}|\boldsymbol{Y}_{0}]=\boldsymbol{\Sigma}^{\text{mix}}_{N}$.
Applying Lemma 14 then gives
$\mathbb{E}[\widehat{v}^{\text{mix},\text{HRK}}_{0}|\boldsymbol{Y}_{0}]=v^{\text{mix}}_{0}$.
The proof is complete. ∎

###  D.5 Proof of Lemma 14

###### Proof.

By linearity of expectations,

$\displaystyle\mathbb{E}[\widehat{v}^{\text{mix}}_{0}|\boldsymbol{Y}_{0}]$
$\displaystyle=\mathbb{E}[\widehat{v}_{0}^{\text{hz}}|\boldsymbol{Y}_{0}]+\mathbb{E}[\widehat{v}_{0}^{\text{vt}}|\boldsymbol{Y}_{0}]-\mathbb{E}[\operatorname{tr}(\boldsymbol{Y}_{0}^{\dagger}\widehat{\boldsymbol{\Sigma}}_{T}(\boldsymbol{Y}^{\prime}_{0})^{\dagger}\widehat{\boldsymbol{\Sigma}}_{N})|\boldsymbol{Y}_{0}].$
(216)

We evaluate each term in (216).

Beginning with the first term, note that the randomness in
$\widehat{\boldsymbol{\Sigma}}_{T}$ stems from $\boldsymbol{\varepsilon}_{T}$
and $\widehat{\boldsymbol{\beta}}$ is deterministic given
$(\boldsymbol{\varepsilon}_{N},\boldsymbol{Y}_{0})$. As such, Assumptions 5–6
with Lemma 13 gives

$\displaystyle\mathbb{E}[\widehat{v}_{0}^{\text{hz}}|\boldsymbol{Y}_{0}]$
$\displaystyle=\mathbb{E}[\widehat{\boldsymbol{\beta}}^{\prime}\widehat{\boldsymbol{\Sigma}}_{T}\widehat{\boldsymbol{\beta}}|\boldsymbol{Y}_{0}]$
(217)
$\displaystyle=\mathbb{E}\left[\mathbb{E}[\widehat{\boldsymbol{\beta}}^{\prime}\widehat{\boldsymbol{\Sigma}}_{T}\widehat{\boldsymbol{\beta}}|\boldsymbol{\varepsilon}_{N},\boldsymbol{Y}_{0}]|\boldsymbol{Y}_{0}\right]$
(218)
$\displaystyle=\mathbb{E}\left[\boldsymbol{y}^{\prime}_{N}\boldsymbol{Y}^{\dagger}_{0}\mathbb{E}[\widehat{\boldsymbol{\Sigma}}_{T}|\boldsymbol{Y}_{0}](\boldsymbol{Y}^{\prime}_{0})^{\dagger}\boldsymbol{y}_{N}|\boldsymbol{Y}_{0}\right]$
(219)
$\displaystyle=\mathbb{E}\left[(\boldsymbol{Y}^{\prime}_{0}\boldsymbol{\beta}^{*}+\boldsymbol{\varepsilon}_{N})\boldsymbol{Y}^{\dagger}_{0}\mathbb{E}[\widehat{\boldsymbol{\Sigma}}_{T}|\boldsymbol{Y}_{0}](\boldsymbol{Y}^{\prime}_{0})^{\dagger}(\boldsymbol{Y}^{\prime}_{0}\boldsymbol{\beta}^{*}+\boldsymbol{\varepsilon}_{N})|\boldsymbol{Y}_{0}\right]$
(220)
$\displaystyle=(\boldsymbol{H}^{u}\boldsymbol{\beta}^{*})^{\prime}\mathbb{E}[\widehat{\boldsymbol{\Sigma}}_{T}|\boldsymbol{Y}_{0}](\boldsymbol{H}^{u}\boldsymbol{\beta}^{*})+\operatorname{tr}(\boldsymbol{Y}_{0}^{\dagger}\mathbb{E}[\widehat{\boldsymbol{\Sigma}}_{T}|\boldsymbol{Y}_{0}](\boldsymbol{Y}^{\prime}_{0})^{\dagger}\boldsymbol{\Sigma}^{\text{mix}}_{N}).$
(221)

By an analogous argument, we derive

$\displaystyle\mathbb{E}[\widehat{v}_{0}^{\text{vt}}|\boldsymbol{Y}_{0}]=(\boldsymbol{H}^{v}\boldsymbol{\alpha}^{*})^{\prime}\mathbb{E}[\widehat{\boldsymbol{\Sigma}}_{N}|\boldsymbol{Y}_{0}](\boldsymbol{H}^{v}\boldsymbol{\alpha}^{*})+\operatorname{tr}(\boldsymbol{Y}_{0}^{\dagger}\boldsymbol{\Sigma}^{\text{mix}}_{T}(\boldsymbol{Y}^{\prime}_{0})^{\dagger}\mathbb{E}[\widehat{\boldsymbol{\Sigma}}_{N}|\boldsymbol{Y}_{0}]).$
(222)

Finally, we use the linearity of the trace operator with Assumption 6 to
obtain

$\displaystyle\mathbb{E}[\operatorname{tr}(\boldsymbol{Y}_{0}^{\dagger}\widehat{\boldsymbol{\Sigma}}_{T}(\boldsymbol{Y}^{\prime}_{0})^{\dagger}\widehat{\boldsymbol{\Sigma}}_{N})|\boldsymbol{Y}_{0}]$
$\displaystyle=\mathbb{E}\left[\mathbb{E}[\operatorname{tr}(\boldsymbol{Y}_{0}^{\dagger}\widehat{\boldsymbol{\Sigma}}_{T}(\boldsymbol{Y}^{\prime}_{0})^{\dagger}\widehat{\boldsymbol{\Sigma}}_{N})|\boldsymbol{\varepsilon}_{N},\boldsymbol{Y}_{0}]|\boldsymbol{Y}_{0}\right]$
(223)
$\displaystyle=\mathbb{E}\left[\operatorname{tr}(\boldsymbol{Y}_{0}^{\dagger}\mathbb{E}[\widehat{\boldsymbol{\Sigma}}_{T}|\boldsymbol{Y}_{0}](\boldsymbol{Y}^{\prime}_{0})^{\dagger}\widehat{\boldsymbol{\Sigma}}_{N})|\boldsymbol{Y}_{0}\right]$
(224)
$\displaystyle=\operatorname{tr}(\boldsymbol{Y}_{0}^{\dagger}\mathbb{E}[\widehat{\boldsymbol{\Sigma}}_{T}|\boldsymbol{Y}_{0}](\boldsymbol{Y}^{\prime}_{0})^{\dagger}\mathbb{E}[\widehat{\boldsymbol{\Sigma}}_{N}|\boldsymbol{Y}_{0}]).$
(225)

Putting everything together completes the proof. ∎


# Sign changes of the partial sums of a random multiplicative function

Marco Aymone  Marco Aymone, Departamento de Matemática, Universidade Federal
de Minas Gerais (UFMG), Brazil<EMAIL_ADDRESS>,  Winston Heap  Winston
Heap<EMAIL_ADDRESS>and  Jing Zhao  Jing Zhao<EMAIL_ADDRESS>
###### Abstract.

We provide a simple proof that the partial sums $\sum_{n\leq x}f(n)$ of a
Rademacher random multiplicative function $f$ change sign infinitely often as
$x\to\infty$, almost surely.

##  1\. Introduction

The behaviour of partial sums of multiplicative functions has a long history
exhibiting deep connections with arithmetic. The partial sums of the Möbius
function for instance, $M(x)=\sum_{n\leqslant x}\mu(n)$, are intimately
related to analytic properties of the Riemann zeta function. Indeed, it is
well known that the bound $M(x)\ll x^{1/2+\epsilon}$ for all $\epsilon>0$ is
equivalent to the Riemann Hypothesis. By standard methods one easily obtains
the complementary bound $M(x)=\Omega(\sqrt{x})$ unconditionally. The current
gap between these two bounds stands at a factor much smaller than
$x^{\epsilon}$ due to Soundararajan [10] who showed that $M(x)\ll
x^{1/2}\exp((\log x)^{1/2}(\log\log x)^{14})$. The true size of these
fluctuations may be yet smaller. Gonek (unpublished) and Ng [9] have
conjectured that the correct order should be around $x^{1/2}(\log\log\log
x)^{5/4}$, although this seems rather deep and out of reach at the moment. A
related, but more tractable, problem in this area is that of sign changes of
$M(x)$. By relatively simple arguments involving the Dirichlet series of
$1/\zeta(s)$, one can show that $M(x)$ changes sign infinitely often as
$x\to\infty$111We say that any function $M(x)$ changes its sign infinitely
often as $x\to\infty$ if neither of the two following inequalities holds for
all sufficiently large $x$: $M(x)\geq 0$ or $M(x)\leq 0$..

In 1944, Wintner [11] introduced a random model for $M(x)$. Let
$(f(p))_{p\text{ prime}}$ be a set of i.i.d. random variables taking values
$\pm 1$ with probability $1/2$ (Rademacher variables) and extend this by
multiplicativity to $\mathbb{N}$ with support on the squarefree integers. So
for example $f(6)=f(2)f(3)$ and $f(4)=0$. Wintner showed that the bound

$\mathcal{M}(x):=\sum_{n\leqslant x}f(n)\ll x^{1/2+\epsilon}$

holds, for all $\epsilon>0$, almost surely. This bound has since seen several
improvements, the current best being due independently to Basquin [4] and
Lau–Tenenbaum–Wu [8] who showed that

$\mathcal{M}(x)\ll x^{1/2}(\log\log x)^{2+\epsilon}$

almost surely. A very recent result due to Harper [7] shows there almost
surely exist arbitrarily large values of $x$ for which

$|\mathcal{M}(x)|\geqslant\sqrt{x}(\log\log x)^{1/4+o(1)}.$

Harper has indicated222Private communication. that his methods may be modified
to produce both large positive and negative values of $\mathcal{M}(x)$ for
arbitrarily large $x$, thus proving that $\mathcal{M}(x)$ has an infinite
number of sign changes. In this note we aim to give a simple proof of this
fact akin to the proof in the deterministic case.

######  Theorem 1.1.

Let $f$ be a Rademacher random multiplicative function and $\mathcal{M}(x)$
its partial sums up to $x$. Then $\mathcal{M}(x)$ changes its sign infinitely
often as $x\to\infty$, almost surely.

A classical tool in analytic number theory to analyse questions of sign
changes is Landau’s oscillation Theorem, which requires that the Dirichlet
series $F(s):=\sum_{n=1}^{\infty}\frac{f(n)}{n^{s}}$ is analytic in some open
set containing its abscissa of convergence. The problem here is that for a
Rademacher random multilpicative function $f$, the Dirichlet series $F(s)$ is
quite irregular near its abcissa of convergence, which is equal to $1/2$. Our
proof of Theorem 1.2 leverages two facts: that almost surely, $F(\sigma)\to 0$
as $\sigma\to 1/2^{+}$, and for any real $t\neq 0$, $\limsup_{\sigma\to
1/2^{+}}|F(\sigma+it)|=\infty$. After partial summation this roughly
translates to the almost sure limits

$\int_{1}^{\infty}\frac{\mathcal{M}(x)}{x^{1+\sigma}}dx\to
0,\qquad\int_{1}^{\infty}\frac{|\mathcal{M}(x)|}{x^{1+\sigma}}dx\to\infty$

as $\sigma\to 1/2^{+}$ which visibly captures the sign changes.

This method allow us to deduce a result relatively stronger than Theorem 1.1.

######  Theorem 1.2.

Let $f$ be a Rademacher random multiplicative function and
$\lambda:[1,\infty)\to\mathbb{R}$ be any function such that
$\int_{1}^{\infty}\frac{|\lambda(u)|}{u^{3/2}}du<\infty$. Then,
$\mathcal{M}(x)+\lambda(x)$ changes sign infinitely often as $x\to\infty$,
almost surely.

We can take, for instance, in the nearly extremal case,
$\lambda(x)=\pm\frac{\sqrt{x}}{(\log x)(\log\log x)^{2}}$ and conclude that
for these choices $\mathcal{M}(x)+\lambda(x)$ changes sign infinitely often,
almost surely. As a consequence, we acquire a relatively simple proof of the
fact that

$\mathcal{M}(x)=\Omega_{\pm}\left(\frac{\sqrt{x}}{(\log x)(\log\log
x)^{2}}\right),$

almost surely. Since $\mathcal{M}(x)$ is always an integer and
$f(n)\in\\{0,\pm 1\\}$ for all positive integers $n$, we deduce the following
corollary.

######  Corollary 1.1.

For any integer $z$, $\mathcal{M}(x)=z$ for an infinite number of integers
$x$, almost surely.

In probabilistic language, we can view $\mathcal{M}(x)$ as a (multiplicative)
random walk. We say that a Markov chain is recurrent if it visits at least one
of its site infinitely often, almost surely. Otherwise, we say that the Markov
chain is transient. An interesting phenomenon is that the simple random walk
is recurrent in dimensions $1$ and $2$, but transient in dimension $3$ and
above. Corollary 1.1 can therefore be interpreted as saying that
$\mathcal{M}(x)$ is recurrent (although strictly speaking $\mathcal{M}(x)$ is
not a Markov chain due to the dependence structure of $f(n)$).

##  2\. Proof of the main result

Clearly Theorem 1.1 is a consequence of Theorem 1.2 by choosing $\lambda(x)=0$
for all $x$. Thus, we will focus in the proof of Theorem 1.2.

In what follows, $p$ denotes a generic prime number, and $\sigma$ a parameter
bigger than $1/2$. By the Kolmogorov one-series Theorem, the Dirichlet series
$\sum_{p}\frac{f(p)}{p^{s}}$ converges almost surely, provided that
$Re(s)=\sigma>1/2$. This result will be implicit in all of the following
Lemmas.

###  2.1. Some results for random Dirichlet series with independent summands

######  Lemma 2.1.

Let $(f(p))_{p\text{ prime }}$ be Rademacher random variables. Then, for any
$\epsilon>0$ there exists a sequence $\sigma_{k}\to 1/2^{+}$ such that

$\sum_{p}\frac{f(p)}{p^{\sigma_{k}}}\ll\left(\log\left(\frac{1}{2\sigma_{k}-1}\right)\right)^{1/2+\epsilon},$

almost surely.

###### Proof.

The series $\sum_{p}\frac{f(p)}{p^{\sigma}}$ has variance

$\sum_{p}\frac{1}{p^{2\sigma}}=\log\zeta(2\sigma)+O(1)=\log\left(\frac{1}{2\sigma-1}\right)+O(1).$

Hoeffding’s inequality for infinite summands (for a proof see e.g.[2], Lemmas
2.2 and 2.3) states that for any $\lambda>0$

$\mathbb{P}\left(\sum_{p}\frac{f(p)}{p^{\sigma}}\geq\lambda\right)\leq\exp\left(-\frac{\lambda^{2}}{2\sum_{p}\frac{1}{p^{2\sigma}}}\right).$

Since the random variables $(f(p))_{p}$ are symmetric, we can replace in the
probability above, the infinite sum by its absolute value at a cost of twice
the upper bound. Therefore, by choosing
$\lambda=\sqrt{2}(\log(1/(2\sigma-1)))^{1/2+\epsilon}$, we obtain that for
some constant $c>0$

$\mathbb{P}\left(\left|\sum_{p}\frac{f(p)}{p^{\sigma}}\right|\geq\sqrt{2}(\log(1/(2\sigma-1)))^{1/2+\epsilon}\right)\leq
c\exp\left(-(\log(1/(2\sigma-1)))^{2\epsilon}\right).$

Thus, by the Borel-Cantelli lemma, for any sequence $\sigma_{k}\to 1/2^{+}$
such that

$\sum_{k=1}^{\infty}\exp\left(-(\log(1/(2\sigma_{k}-1)))^{2\epsilon}\right)<\infty,$

we have that the target bound holds almost surely along this sequence. ∎

######  Remark 2.1.

It is interesting to observe that the upper bound in Lemma 2.1 holds for all
$\sigma$ sufficiently close to $1/2^{+}$, not just for a sequence. A proof of
this can be done by following the steps of [2] from which we can deduce a
sharp upper bound at the level of the law of the iterated logarithm. For our
purposes, we only need this weaker bound for a sequence $\sigma_{k}$, we thank
Adam Harper for pointing this out.

######  Lemma 2.2.

Let $(f(p))_{p\text{ prime}}$ be Rademacher random variables. Then, for each
fixed $t\in\mathbb{R}\setminus\\{0\\}$, the infinite sum

$\sqrt{\frac{2}{\log\left(\frac{1}{2\sigma-1}\right)}}\sum_{p}\frac{f(p)\cos(t\log
p)}{p^{\sigma}}$

converges in probability distribution, as $\sigma\to 1/2^{+}$, to a standard
Gaussian random variable.

###### Proof.

Let

$\displaystyle S(\sigma)$ $\displaystyle:=\sum_{p}\frac{f(p)\cos(t\log
p)}{p^{\sigma}},$ $\displaystyle V(\sigma)$
$\displaystyle:=\mathbb{E}S(\sigma)^{2}=\sum_{p}\frac{\cos^{2}(t\log
p)}{p^{2\sigma}}.$

Now, $\cos^{2}(t\log p)=\frac{1+\cos(2t\log p)}{2}$ and

$\sum_{p}\frac{\cos(2t\log
p)}{2p^{2\sigma}}=\frac{1}{2}Re(\log\zeta(2\sigma+2it))+O(1),$

which, by the classical estimates for the Riemann zeta function, is $O(1)$
provided that $t\neq 0$ is fixed. Thus, we have that the variance

$V(\sigma)=\frac{1}{2}\log\left(\frac{1}{2\sigma-1}\right)+O(1).$

Now we argue as in the proof of Lemma 3.3 of [1] to show the convergence to
the Gaussian. We do this by the method of characteristic functions.

By the independence of $(f(p))_{p\text{ prime}}$ and the dominated convergence
theorem:

$\displaystyle\varphi_{\sigma}(s)$
$\displaystyle:=\mathbb{E}\exp\left(\frac{isS(\sigma)}{\sqrt{V(\sigma)}}\right)=\prod_{p}\cos\left(\frac{s\cos(t\log
p)}{\sqrt{V(\sigma)}p^{\sigma}}\right).$

Our aim is to show that for each fixed $s\in\mathbb{R}$,
$\varphi_{\sigma}(s)\to\exp(-s^{2}/2)$ as $\sigma\to 1/2^{+}$. Observe that
$\varphi_{\sigma}(s)$ is an even function of $s$, so we may assume $s\geq 0$.
Also, note that for each $s\geq 0$ we may choose $\sigma>1/2$ such that

$\left|\frac{s\cos(t\log
p)}{\sqrt{V(\sigma)}p^{\sigma}}\right|\leq\frac{1}{100},\,\,\,\,\,\,\,\,\,\,0\leq
1-\cos\left(\frac{s\cos(t\log
p)}{\sqrt{V(\sigma)}p^{\sigma}}\right)\leq\frac{1}{100}$

for all primes $p$ since $V(\sigma)$ gets large as $\sigma\to 1/2^{+}$.

For $|x|\leq 1/100$, we have that $\log(1-x)=-x+O(x^{2})$ and
$\cos(x)=1-\frac{x^{2}}{2}+O(x^{4})$. Thus, we have:

$\displaystyle\log\varphi_{\sigma}(s)$
$\displaystyle=\sum_{p}\log\cos\left(\frac{s\cos(t\log
p)}{\sqrt{V(\sigma)}p^{\sigma}}\right)$
$\displaystyle=\sum_{p}\log\left(1-\left(1-\cos\left(\frac{s\cos(t\log
p)}{\sqrt{V(\sigma)}p^{\sigma}}\right)\right)\right)$
$\displaystyle=-\sum_{p}\left(1-\cos\left(\frac{s\cos(t\log
p)}{\sqrt{V(\sigma)}p^{\sigma}}\right)\right)+O\left(\sum_{p}\left(1-\cos\left(\frac{s\cos(t\log
p)}{\sqrt{V(\sigma)}p^{\sigma}}\right)\right)^{2}\right)$
$\displaystyle=-\sum_{p}\frac{s^{2}\cos^{2}(t\log
p)}{2V(\sigma)p^{2\sigma}}+O\left(\sum_{p}\frac{s^{4}}{V^{2}(\sigma)p^{4\sigma}}\right)$
$\displaystyle=-\frac{s^{2}}{2}+O\left(\frac{s^{4}}{V^{2}(\sigma)}\right).$

We conclude that $\varphi_{\sigma}(s)\to\exp(-s^{2}/2)$ as $\sigma\to
1/2^{+}$. ∎

######  Lemma 2.3.

Let $(f(p))_{p\text{ prime}}$ be Rademacher random variables. Then, for each
fixed $t\in\mathbb{R}\setminus\\{0\\}$, and any sequence $\sigma_{k}\to
1/2^{+}$, we have that

$\limsup_{k\to\infty}\sum_{p}\frac{f(p)\cos(t\log p)}{p^{\sigma_{k}}}=\infty,$

almost surely.

###### Proof.

Let $S(\sigma)$ and $V(\sigma)$ be as in the proof of Lemma 2.2. Observe that
the event in which

$\limsup_{k\to\infty}\frac{S(\sigma_{k})}{\sqrt{V(\sigma_{k})}}\geq 1$

is a tail event, in the sense that it does not depend on the random variables
$(f(p))_{p\leq y}$, for any fixed and large $y$. Now,

$\displaystyle\mathbb{P}\left(\limsup_{k\to\infty}\frac{S(\sigma_{k})}{\sqrt{V(\sigma_{k})}}\geq
1\right)$
$\displaystyle\geq\mathbb{P}\left(\frac{S(\sigma_{k})}{\sqrt{V(\sigma_{k})}}\geq
1,\text{ for infinitely many }k\right)$
$\displaystyle=\lim_{n\to\infty}\mathbb{P}\left(\bigcup_{k=n}^{\infty}\left[\frac{S(\sigma_{k})}{\sqrt{V(\sigma_{k})}}\geq
1\right]\right)$
$\displaystyle\geq\frac{1}{\sqrt{2\pi}}\int_{1}^{\infty}e^{-s^{2}/2}ds,$

where in the last inequality above we used the convergence to the Gaussian
from Lemma 2.2. Thus, our tail event has positive probability and hence by the
Kolmogorov $0-1$ law this probability must be $1$. ∎

###  2.2. Some results for Rademacher random multiplicative functions

The next result is essentially due to Wintner [11], and for a proof we refer
the reader to Lemma 2.1 of [3].

######  Lemma 2.4.

Let $f$ be a Rademacher random multiplicative function and
$F(s):=\sum_{n=1}^{\infty}\frac{f(n)}{n^{s}}$. Then $F(s)$ converges for all
$s$ in the half plane $Re(s)>1/2$ and

$\log
F(s)=\sum_{p}\frac{f(p)}{p^{s}}-\frac{1}{2}\sum_{p}\frac{1}{p^{2s}}+\sum_{p}\sum_{m=3}^{\infty}\frac{(-1)^{m+1}f(p)^{m}}{mp^{ms}},$

almost surely.

######  Lemma 2.5.

Let $f$ be a Rademacher random multiplicative function and
$\lambda:[1,\infty)\to\mathbb{R}$ be such that
$\int_{1}^{\infty}\frac{|\lambda(u)|}{u^{3/2}}du<\infty$. Then, almost surely

$\int_{1}^{\infty}\frac{1}{u^{3/2}}\bigg{|}\lambda(u)+\sum_{n\leq
u}f(n)\bigg{|}du=\infty.$

###### Proof.

Assume, by contradiction, that with positive probability the integral above
converges. Then, with positive probability, for each $t\in\mathbb{R}$, the
integral

$\int_{1}^{\infty}\frac{1}{u^{3/2+it}}\left(\lambda(u)+\sum_{n\leq
u}f(n)\right)du$

converges absolutely. Thus, by the dominated convergence theorem

$\lim_{\sigma\to
1/2^{+}}\int_{1}^{\infty}\frac{1}{u^{1+\sigma+it}}\left(\lambda(u)+\sum_{n\leq
u}f(n)\right)du=\int_{1}^{\infty}\frac{1}{u^{3/2+it}}\left(\lambda(u)+\sum_{n\leq
u}f(n)\right)du.$

By the dominated convergence theorem again, the limit

$\lim_{\sigma\to
1/2^{+}}\int_{1}^{\infty}\frac{\lambda(u)}{u^{1+\sigma+it}}du$

exists and is finite. Thus, the following limit exists and is finite:

$\lim_{\sigma\to
1/2^{+}}\int_{1}^{\infty}\frac{1}{u^{1+\sigma+it}}\left(\sum_{n\leq
u}f(n)\right)du.$

Now we recall that by partial summation the Dirichlet series
$F(\sigma+it):=\sum_{n=1}^{\infty}\frac{f(n)}{n^{\sigma+it}}$ satisfies

$F(\sigma+it)=(\sigma+it)\int_{1}^{\infty}\frac{1}{u^{1+\sigma+it}}\sum_{n\leq
u}f(n)du,$

where this integral converges absolutely for each $\sigma>1/2$ almost surely.
Thus we have that for each real $t\neq 0$, the limit $\lim_{\sigma\to
1/2^{+}}F(\sigma+it)$ exists and is finite with positive probability. Now, by
Lemma 2.4, as $\sigma\to 1/2^{+}$, almost surely we have that

$\displaystyle\log F(\sigma+it)$
$\displaystyle=\sum_{p}\frac{f(p)}{p^{\sigma+it}}-\frac{1}{2}\sum_{p}\frac{1}{p^{2(\sigma+it)}}+\sum_{p}\sum_{m=3}^{\infty}\frac{(-1)^{m+1}f(p)^{m}}{mp^{m(\sigma+it)}}$
$\displaystyle=\sum_{p}\frac{f(p)}{p^{\sigma+it}}-\frac{\log\zeta(2\sigma+2it)}{2}+O(1)$
$\displaystyle=\sum_{p}\frac{f(p)}{p^{\sigma+it}}+O(1),$

provided that $t$ is real, fixed and $\neq 0$. Since
$Re\sum_{p}\frac{f(p)}{p^{\sigma+it}}=\sum_{p}\frac{f(p)\cos(t\log
p)}{p^{\sigma}}$ and $|F(\sigma+it)|=\exp(Re\log F(\sigma+it))$, by Lemma 2.3,
almost surely $|F(\sigma+it)|$ tends to $\infty$ along a subsequence
$\sigma=\sigma_{k}\to 1/2^{+}$. This gives the desired contradiction. ∎

######  Remark 2.2.

At this point it is interesting to compare our method with that of Halász [5].
Concerning omega bounds for $|\sum_{n\leq x}f(n)|$, Halász exploited the fact
that we can obtain a much larger bound if we consider

$\sup_{1\leq t\leq 2}\sum_{p}\frac{f(p)\cos(t\log p)}{p^{\sigma}}.$

Indeed, he obtained a lower bound for this supremum that is slightly smaller
than twice the variance of this random sum evaluated at any fixed $t\neq 0$,
with high probability. Later, this idea was exploited by Harper in [6] where
he improved Halász’s omega bound for the partial sums. In the case of sign
changes, we stress that we only need to show that this random sum blows up
along a sequence $\sigma_{k}\to 1/2^{+}$ for a fixed value of $t\neq 0$,
almost surely.

######  Lemma 2.6.

Let $f$ be a Rademacher random multiplicative function and $F(\sigma)$ its
Dirichlet series. Then, there exists a sequence $\sigma_{k}\to 1/2^{+}$ such
that $F(\sigma_{k})=(2\sigma_{k}-1)^{1/2+o(1)}$, almost surely.

###### Proof.

Indeed, by Lemmas 2.1 and 2.4, for any $0<\epsilon<1/2$, there exists a
sequence $\sigma_{k}\to 1/2^{+}$ such that almost surely:

$\displaystyle\log F(\sigma_{k})$
$\displaystyle=\sum_{p}\frac{f(p)}{p^{\sigma_{k}}}-\frac{1}{2}\sum_{p}\frac{1}{p^{2\sigma_{k}}}+O(1)$
$\displaystyle=\sum_{p}\frac{f(p)}{p^{\sigma_{k}}}-\frac{\log\zeta(2\sigma_{k})}{2}+O(1)$
$\displaystyle=O\left(\left(\log\left(\frac{1}{2\sigma_{k}-1}\right)\right)^{1/2+\epsilon}\right)-\frac{1}{2}\log\left(\frac{1}{2\sigma_{k}-1}\right)$
$\displaystyle=(1/2+o(1))\log(2\sigma_{k}-1).$

∎

###  2.3. Proof of Theorem 1.2

###### Proof.

Assume by contradiction that with positive probability,
$\lambda(x)+\sum_{n\leq x}f(n)$ does not changes sign infinitely often as
$x\to\infty$, and moreover, say that for some (random) $x_{0}$,
$\lambda(x)+\sum_{n\leq x}f(n)\geq 0$, for all $x\geq x_{0}$. Thus, for
sufficiently large $x$ we have

$\int_{1}^{x}\frac{1}{u^{3/2}}\left(\lambda(u)+\sum_{n\leq
u}f(n)\right)du=O(1)+\int_{1}^{x}\frac{1}{u^{3/2}}\left|\lambda(u)+\sum_{n\leq
u}f(n)\right|du.$

Now we recall that $F(\sigma):=\sum_{n=1}^{\infty}\frac{f(n)}{n^{\sigma}}$ is
equal to $\sigma\int_{1}^{\infty}\frac{1}{u^{1+\sigma}}\sum_{n\leq u}f(n)du$.
Thus, by Lemma 2.6, along a sequence $\sigma_{k}\to 1/2^{+}$, we have that
with positive probability:

$\displaystyle(2\sigma_{k}-1)^{1/2+o(1)}+\sigma_{k}\int_{1}^{\infty}\frac{\lambda(u)}{u^{1+\sigma_{k}}}du$
$\displaystyle=F(\sigma_{k})+\sigma_{k}\int_{1}^{\infty}\frac{\lambda(u)}{u^{1+\sigma_{k}}}du$
$\displaystyle\geq\sigma_{k}\int_{1}^{x}\frac{1}{u^{1+\sigma_{k}}}\left(\lambda(u)+\sum_{n\leq
u}f(n)\right)du.$

Making the limit $\sigma_{k}\to 1/2^{+}$ in the inequality above, we obtain
that

$O(1)\geq\int_{1}^{x}\frac{1}{u^{3/2}}\left|\lambda(u)+\sum_{n\leq
u}f(n)\right|du,$

which is by Lemma 2.5, a contradiction for sufficiently large $x$.

Now assume that with positive probability $\lambda(x)+\sum_{n\leq x}f(n)\leq
0$, for all $x\geq x_{0}$. Arguing as above, we have that

$\int_{1}^{x}\frac{1}{u^{3/2}}\left(\lambda(u)+\sum_{n\leq
u}f(n)\right)du=O(1)-\int_{1}^{x}\frac{1}{u^{3/2}}\left|\lambda(u)+\sum_{n\leq
u}f(n)\right|du.$

Thus, with positive probability:

$\displaystyle(2\sigma_{k}-1)^{1/2+o(1)}+\sigma_{k}\int_{1}^{\infty}\frac{\lambda(u)}{u^{1+\sigma}}du$
$\displaystyle=F(\sigma_{k})+\sigma_{k}\int_{1}^{\infty}\frac{\lambda(u)}{u^{1+\sigma_{k}}}du$
$\displaystyle\leq\sigma_{k}\int_{1}^{x}\frac{1}{u^{1+\sigma_{k}}}\left(\lambda(u)+\sum_{n\leq
u}f(n)\right)du.$

Making the limit $\sigma_{k}\to 1/2^{+}$ in the inequality above, we obtain
that

$O(1)\leq-\int_{1}^{x}\frac{1}{u^{3/2}}\left|\lambda(u)+\sum_{n\leq
u}f(n)\right|du,$

which is again by Lemma 2.5, a contradiction for sufficiently large $x$. ∎

##  3\. The multiplicative random walk on $\mathbb{Z}^{2}$

As said before, Corollary 1.1 can be interpreted as saying that the
multiplicative random walk is recurrent. A natural question could be what
happens in dimension $2$.

By defining $(f(p))_{p\text{ prime}}$ to be i.i.d. with $f(p)$ uniformly
distributed over the set $\\{\pm 1,\pm i\\}$, and extending $f$ to the
positive integers multiplicatively with support on the squarefree integers,
then for any positive integer $n$, $f(n)\in\\{0,\pm 1,\pm i\\}$. Therefore,
$\mathcal{M}(x)=\sum_{n\leq x}f(n)$ is always a number of the form $a+bi$,
where $a$ and $b$ are integers, and hence can be seen as the multiplicative
random walk in $2$ dimensions.

Notice that in two dimensions we lose the notion of sign changes, but we can
still ask about recurrence. We conclude this section with the following
question.

Question. In two dimensions, does the equation $\mathcal{M}(x)=0$ have an
infinite number of integer solutions?

##  4\. Some simulations

Below we plot a sample from the multiplicative random walk in dimensions $1$
and $2$.

The multiplicative random walk in dimension $1$

Figure 1. The dashed curve is given by $x\mapsto\pm x^{1/2}$, and the
continuous curve is given by the multiplicative random walk in dimension $1$:
$x\mapsto\sum_{n\leq x}f(n)$, where $1\leq x\leq 10^{7}$.

The multiplicative random walk in dimension $2$

Figure 2. Let $M(x):=\sum_{n\leq x}f(n)$. The blue part is given by set of
points visited by the walk, more precisely $\\{M(x)\in\mathbb{C}:1\leq x\leq
10^{6}\\}$.

Acknowledgements. We would like to thank the referee for a careful reading of
the paper and for important and useful suggestions and corrections. MA is
supported by CNPq, grant Universal no. 403037/2021-2.

## References

  * [1] M. Aymone, Real zeros of random Dirichlet series, Electron. Commun. Probab., 24 (2019), pp. Paper No. 54, 8.
  * [2] M. Aymone, S. Frómeta, and R. Misturini, Law of the iterated logarithm for a random Dirichlet series, Electron. Commun. Probab., 25 (2020), pp. Paper No. 56, 14.
  * [3] M. Aymone and V. Sidoravicius, Partial sums of biased random multiplicative functions, J. Number Theory, 172 (2017), pp. 343–382.
  * [4] J. Basquin, Sommes friables de fonctions multiplicatives aléatoires, Acta Arith., 152 (2012), pp. 243–266.
  * [5] G. Halász, On random multiplicative functions, in Hubert Delange colloquium (Orsay, 1982), vol. 83 of Publ. Math. Orsay, Univ. Paris XI, Orsay, 1983, pp. 74–96.
  * [6] A. J. Harper, Bounds on the suprema of Gaussian processes, and omega results for the sum of a random multiplicative function, Ann. Appl. Probab., 23 (2013), pp. 584–616.
  * [7] A. J. Harper, Almost Sure Large Fluctuations of Random Multiplicative Functions, International Mathematics Research Notices, (2021).  rnab299.
  * [8] Y.-K. Lau, G. Tenenbaum, and J. Wu, On mean values of random multiplicative functions, Proc. Amer. Math. Soc., 141 (2013), pp. 409–420.
  * [9] N. Ng, The distribution of the summatory function of the Möbius function, Proc. London Math. Soc. (3), 89 (2004), pp. 361–389.
  * [10] K. Soundararajan, Partial sums of the Möbius function, J. Reine Angew. Math., 631 (2009), pp. 141–152.
  * [11] A. Wintner, Random factorizations and Riemann’s hypothesis, Duke Math. J., 11 (1944), pp. 267–275.