# Hole radii for the Kac polynomials and derivatives

Hoi H. Nguyen and Oanh Nguyen  Department of Mathematics
The Ohio State University
231 W 18th Ave
Columbus, OH 43210 USA<EMAIL_ADDRESS>Division of Applied Mathematics
Brown University
Providence, RI 02906, USA<EMAIL_ADDRESS>
###### Abstract.

The Kac polynomial

$f_{n}(x)=\sum_{i=0}^{n}\xi_{i}x^{i}$

with independent coefficients of variance 1 is one of the most studied models
of random polynomials.

It is well-known that the empirical measure of the roots converges to the
uniform measure on the unit disk. On the other hand, at any point on the unit
disk, there is a hole in which there are no roots, with high probability. In a
beautiful work [11], Michelen showed that the holes at $\pm 1$ are of order
$1/n$. We show that in fact, all the hole radii are of the same order. The
same phenomenon is established for the derivatives of the Kac polynomial as
well.

H.N. is supported by NSF CAREER grant DMS 1752345, O.N. is supported by NSF
grants DMS–1954174 and DMS–2246575. This work was initiated under the SQuaREs
2021 program of AIM, we thank the Institute for the generous support.

##  1\. Introduction

Approximation by roots of polynomials of coefficients $\\{-1,0,1\\}$ is a
classical and interesting topic in analysis, with fascinating pictures and
conjectures. For instance, it follows from a result of Borwein and Pinner [2,
Theorem 1] that for any given $\zeta$ of $d$-th root of unity, the distance
from it to any root $z$ of $\\{-1,0,1\\}$ polynomials of degree $n$ which
vanish at $\zeta$ of order at most $k$ (i.e. $f^{(k+1)}(\zeta)\neq 0$) can be
bounded by

$|z-\zeta|\geq
e^{-1}\frac{(k!)^{\lceil\phi(d)/2\rceil}}{(n+1)^{(k+1)\lceil\phi(d)/2\rceil+1}},$

where $\phi(d)$ is the usual Euler phi-function. This result is asymptotically
optimal. On the other hand, the smallest distance can be sub-exponentially
small if $\zeta$ is on the unit circle and not a root of unity (such as when
$\zeta$ is an algebraic number of small Mahler measure), see for instance [2,
Corollary 1, Theorem 3]. The situation at 1 is also interesting, it was shown
from the same paper [2, Corollary 4, Theorem 6] that
$|1-z|\geq\frac{1}{n^{k+2}}$ (which is again near optimal) for any real roots
$z$ of $\\{-1,0,1\\}$ polynomials of degree $n$ which vanish at $1$ of order
exactly $k$. Note that the distance is significantly larger if $z$ is purely
complex. See Figure 1. We also refer the reader to [2, 4, 15] and the
references therein for further interesting discussions and problems.

Figure 1. Zeros of all polynomials with $\pm 1$ coefficients and degree at
most eight ([2])

Our goal in this note is to study the distances from some probabilistic
viewpoint. More generally, consider the Kac polynomial

$f_{0,n}(x)=\sum_{i=0}^{n}\xi_{i}x^{i}$

where $\xi_{i}$ are independent (not necessarily identically distributed),
real-valued random variables with mean 0 and variance 1.

For this random polynomial, it is well-known that the empirical distribution
of the roots converges to the uniform distribution on the unit circle ([9]).
So, the roots concentrate near the unit circle. And in particular, the real
roots concentrate near $\pm 1$. However, precisely at $\pm 1$, there are holes
that do not contain any roots. It was conjectured by Shepp and Vabderbei [16],
and confirmed recently by Michelen [11, Theorem 1.2] that the typical distance
of real roots to 1 for random Kac polynomial is of order $O(1/n)$.

######  Theorem 1.1.

Let $\xi_{i}$ be iid with mean zero and variance one. For any constant
$\delta>0$, there exists a constant $C>0$ so that

${\mathbf{P}}(\mbox{there exists a real root in $[1-C/n,1+C/n]$})\geq
1-\delta$

for all $n$ sufficiently large.

It is not hard to establish the lower bound and conclude that the hole at $1$
(and $-1$) has radius of order $\Theta(1/n)$.

How about the hole radius at other points on the unit circle? A recent result
by Cook, Yakir, Zeitouni and the first author [5] (see also Michelen and
Sahasrabudhe [12] for the Gaussian case) shows that the distance between the
zero set of $f_{n,0}$ and the unit circle is of order $\frac{1}{n^{2}}$. So,
this is a lower bound for all hole radii.

From the result in [2] and Figure 1, it is natural to predict that the hole
radii exhibit different orders at different points. For instance, in Figure 2
(source [1]) where all roots of polynomials with coefficients $\pm 1$ and
degree at most 24, one can observe that there are largest holes at $\pm 1$,
smaller holes possibly at the roots of unity, and barely visible holes at
other points. In Figure 3, we draw sampled roots of random Kac polynomials
with $\pm 1$ coefficients and degree $n=1000$. Note that as $n=1000$ is large
compared to $24$, the holes are no longer visible in the figure unless being
zoomed in properly. However, the striking similarities between the two figures
would suggest that the same observation would remain true for large $n$. In
Figure 4, we display sampled roots of the first derivative polynomial with
$\pm 1$ coefficients and degree $n=1000$ which has the same pattern as in
Figure 3.

Figure 2. All roots of polynomials with coefficients $\pm 1$ and degree at
most 24 (source [1]) Figure 3. Sampled roots of the Kac polynomials with
Rademacher coefficients Figure 4. Sampled roots of the first derivative of Kac
polynomials with Rademacher coefficients

Disproving this prediction, in this paper, we show that for every point
$\zeta$ on the unit circle, the precise order of the hole at $\zeta$ is
$\Theta(1/n)$. Moreover, we show that this holds also for the derivatives of
the Kac polynomial. For a positive integer constant $\rho$, let us define

$f_{\rho,n}=\sum_{i=0}^{n}a_{i,\rho,n}\xi_{i}z^{i}$

where $a_{i,\rho,n}=(1+o_{n}(1))i(i-1)\dots(i-\rho+1)\in{\mathbb{R}}$ with the
$o_{n}(1)$ converges to $0$ uniformly in $i$, as $n\to\infty$ and $\rho$
fixed. For any set $S\subset{\mathbb{C}}$, let $N_{f_{\rho,n}}(S)$ be the
number of roots of $f_{\rho,n}$ in $S$. If $\rho=0$, we get the Kac
polynomial. The $\rho$-th derivative of the Kac polynomial corresponds to
$x^{-\rho}f_{\rho,n}(x)$. Here is our main result.

######  Theorem 1.2.

Let $\rho$ be any positive integer constant. Assume that the random variables
$\xi_{i}$ are independent with mean 0, variance 1 and bounded
$(2+\varepsilon_{0})$-moment for some $\varepsilon_{0}>0$. For every $\zeta\in
S^{1}$, the radius of the hole at $\zeta$ is $\Omega(1/n)$. In particular, for
every $\varepsilon>0$, there exist positive constants $c_{\rho}$ and
$C_{\rho}$ such that for all $\zeta\in S^{1}$

$\text{[Upper bound]}\hskip
56.9055pt{\mathbf{P}}(N_{f_{\rho,n}}(B(\zeta,C_{\rho}/n))>0)\geq
1-\varepsilon$ (1)

and

$\text{[Lower bound]}\hskip
56.9055pt{\mathbf{P}}(N_{f_{\rho,n}}(B(\zeta,c_{\rho}/n))=0)\geq
1-\varepsilon$ (2)

for sufficiently large $n$.

As far as we understand, the proof in [11] is restricted to real roots and
cannot be applied to complex roots.

For a related discussion concerning the hole probability, we refer to [3], [7,
Section 5.1.3] and the references therein.

##  2\. Proof sketch and ingredients

For the upper bound, we need to split into two cases: $\zeta=\pm 1$ and
$\zeta\neq\pm 1$. For the former, we show that the proof in [11] can be
adapted to cover $f_{\rho,n}$ for general $\rho\in{\mathbb{N}}$. This method
relies on the simple observation that if a polynomial $f$ changes sign in an
interval on the real line then it has at least one root there. Since this
observation only holds for real roots, it merely works for $\zeta=\pm 1$ where
one can reduce the upper-bound problem to showing that there exists a real
root in the interval centered at $\zeta$ and radius $C/n$. For $\zeta\neq\pm
1$, we need a different approach.

To this end, we note that the expected value of $N_{f_{\rho,n}}(B(\zeta,C/n))$
being large does not imply that the number of real roots is non-zero with high
probability. However, it can be achieved via Chebyshev inequality if we can
show that its variance is of smaller order than its mean squared. To do so,
the high-level idea is to show that $f_{\rho,n}$, after rescaled properly,
converges to a Gaussian process, say $f_{\infty}$ which we can estimate the
growth of $\mathbf{Var}\left(N_{f_{\infty}}(B(\zeta,C/n))\right)$ in terms of
$C$ and then pass the result back to $f_{\rho,n}$. So, we consider a rescaled
version of $f_{\rho,n}$ by zooming in at the local neighborhood of $\zeta$ as
follows

$g_{n}(z)=\frac{1}{n^{\rho+1/2}}f_{\rho,n}\left(\zeta+\frac{1}{n}\zeta
z\right).$

When the random variables are Gaussian, we know that $g_{n}$ is a Gaussian
process with covariance

$\displaystyle{\mathbf{E}}g_{n}(z)\overline{g_{n}(w)}$ $\displaystyle=$
$\displaystyle\frac{1}{n^{2\rho+1}}\sum_{k=0}^{n}|\zeta|^{2k}a_{k,\rho,n}^{2}\left(1+\frac{1}{n}z\right)^{k}\left(1+\frac{1}{n}\bar{w}\right)^{k}$
$\displaystyle=$
$\displaystyle\frac{1}{n^{2\rho+1}}\sum_{k=0}^{n}a_{k,\rho,n}^{2}\left(1+\frac{1}{n}z\right)^{k}\left(1+\frac{1}{n}\bar{w}\right)^{k}.$

We note that $\zeta$ disappears on the right-most side and could potentially
account for why the hole radii are of the same order.

To show the cancellation
$\mathbf{Var}\left(N_{f_{\infty}}(B(\zeta,C/n))\right)=o\left({\mathbf{E}}\left(N_{f_{\infty}}(B(\zeta,C/n))\right)\right)^{2}$,
we note that the variance is an integral of
$\rho_{2}(z_{1},z_{2})-\rho_{1}(z_{1})\rho_{1}(z_{2})$ which is in turn of
order $|z_{1}-z_{2}|^{-1}$. So, when $z_{1}$ and $z_{2}$ are far away, the
integrand is small, accounting for the cancellation. This suggests that the
numbers of roots in far-away regions are weakly correlated which is consistent
with previously established results for real roots ([13]) and radius of
complex roots ([5]). To handle the diagonal region when $z_{1}$ is near
$z_{2}$, we come up with a simple argument, though via rather long and tedious
algebraic manipulations, showing that the function
$\rho_{2}(z_{1},z_{2})-\rho_{1}(z_{1})\rho_{1}(z_{2})$ is indeed continuous
everywhere and hence the diagonal has negligible contribution, see Lemma 4.6.
To carry out this strategy, we actually replace the whole ball $B(\zeta,C/n)$
by a subset, denoted by $U_{C}$, which is a thin strip along the unit circle.
This is a major device that allows us to reduce from $|z_{1}-z_{2}|$ to
$|{\operatorname{Im}}(z_{1})-{\operatorname{Im}}(z_{2})|$ which reduces the
dimension and facilitates the rather elegant proof that follows.

To pass from $f_{\infty}$ to $f_{\rho,n}$, we need to show some sort of
uniform integrability of $N_{f_{\infty}}(B(\zeta,C/n)$. For that, we adapt a
double Taylor expansion argument used in [10] (see Lemma 4.7).

Finally, to establish the lower bound (2), we will show that the expected
number of roots in the ball $B(\zeta,c/n)$ is small for sufficiently small $c$
and then apply Markov’s inequality. The derivation of the expected number of
roots is first reduced to the Gaussian setting when all coefficients $\xi_{i}$
are iid standard Gaussian, via the universality properties of the random
polynomials. To do the calculation for the Gaussian case, there are two
possible ways. The first way is to directly apply the classical Kac-Rice
formula to $f_{\rho,n}$. The second way, which is what we perform here, is to
derive it through $f_{\infty}$ using the limits that we already establish for
the upper bound.

Notations. For the rest of the paper, to simplify the notation, we will often
drop the subscript $\rho$. For instance, we write $f_{n}$ in place of
$f_{\rho,n}$. For a function $f$, let
$\mathcal{Z}(f)=\\{z\in{\mathbb{C}}:f(z)=0\\}$ be the zero set of $f$.

We use standard asymptotic notations under the assumption that $n$ tends to
infinity. For two positive sequences $(a_{n})$ and $(b_{n})$, we say that
$a_{n}\gg b_{n}$ or $b_{n}\ll a_{n}$ if there exists a constant $C$ such that
$b_{n}\leq Ca_{n}$. If $|c_{n}|\ll a_{n}$ for some sequence $(c_{n})$, we also
write $c_{n}\ll a_{n}$.

If $a_{n}\ll b_{n}\ll a_{n}$, we say that $b_{n}=\Theta(a_{n})$. If
$\lim_{n\to\infty}\frac{a_{n}}{b_{n}}=0$, we say that $a_{n}=o(b_{n})$. We
also write that $a_{n}=O_{C}(b_{n})$ if the implied constant depends on a
given parameter $C$.

##  3\. Proof of Theorem 1.2: upper bound for $\zeta=\pm 1$

When $\zeta=\pm 1$, Michelen [11] already showed the stated upper bound for
the Kac polynomial $f_{0,n}$. We will show that this proof can be easily
adapted to cover the general case $f_{\rho,n}$. We assume that $\zeta=1$ as
the case $\zeta=-1$ is completely similar. It suffices to show that with
probability $\geq 1-\varepsilon$, there is at least one root of $f_{\rho,n}$
in the interval $J:=[1-C/n,1+C/n]$ for some large constant $C$. Let $f_{0,n}$
be the $\rho$-th anti-derivative of $z^{-\rho}f_{\rho,n}$, then

$f_{0,n}=\sum_{i=0}^{n}(1+o(1))\xi_{i}z^{i}$

is basically the Kac polynomial (if disregarding the $1+o(1)$ terms). By
interlacing, this can be deduced from showing that there are at least $\rho+1$
roots of $f_{0,n}$ in the same interval. To this end, we show that we can find
$(\rho+1)$ sub-intervals of $J$ each of which observes a sign change of
$f_{0,n}$ and hence contains at least one root.

Consider the rescaled polynomial

$h_{n}(x)=\frac{1}{\sqrt{n}}f_{0,n}(1+x/n),\quad x\in{\mathbb{R}}.$

Let $M$ be a large constant and $x_{1},\dots,x_{M}$ be deterministic points in
$J$. By [11, Lemma 5] (which is a rather direct application of the Lindeberg-
Feller Central Limit Theorem), the random vector
$(h_{n}(x_{1}),\dots,h_{n}(x_{M}))$ converges to the Gaussian vector
$(h(x_{1}),\dots,h(x_{M}))$ where $h$ is a centered, real Gaussian process
with covariance

${\mathbf{E}}h(x)h(y)=\int_{0}^{1}e^{(x+y)t}dt.$

(In fact, in [11], this result is established for the Kac polynomial without
the $1+o(1)$ terms as above but the proof can easily go through without
changes when these terms are present.)

Let $M=(\rho+1)K$ where $K$ is a large constant to be chosen. By [11, Lemma
6], there exists a constant $\gamma>0$ such that any centered Gaussian vector
$(Z_{1},\dots,Z_{K})$ with variances ${\mathbf{E}}Z_{i}^{2}=1$ for all $i$ and
covariances $|{\mathbf{E}}Z_{i}Z_{j}|\leq\gamma$ for all $i\neq j$ satisfies

${\mathbf{P}}(\text{all $Z_{1},\dots,Z_{K}$ have the same sign})\leq
2^{-K+2}.$

Direct calculation shows that if $y=\alpha x$ with $x>1$ and $\alpha>1$ then

$\left|{\mathbf{E}}\frac{h(x)h(y)}{\sqrt{\mathbf{Var}h(x)\mathbf{Var}h(y)}}\right|=o_{\alpha\to\infty}(1).$

So, for a given $K$, by taking $\alpha$ sufficiently large, we can make this
number smaller than $\gamma$. We then take $x_{i}=\alpha^{i-1}$,
$i=1,\dots,M$. So, for $j=0,\dots,\rho$,

${\mathbf{P}}(\text{all $h(x_{jK+1}),\dots,h(x_{jK+K})$ have the same
sign})\leq 2^{-K+2}.$

And so, for sufficiently large $n$,

${\mathbf{P}}(\text{$h_{n}$ does not have any real roots in
$[1+x_{jK+1}/n,1+x_{jK+K}/n]$})\leq 2^{-K+3}.$

By the union bound, the probability that $h_{n}$ has less than $\rho+1$ real
roots in $[1+x_{1}/n,1+x_{M}/n]$ (which is a union of $\rho+1$ such intervals
above) is at most $(\rho+1)2^{-K+3}$. By choosing $K$ sufficiently large so
that this number is smaller than $\varepsilon$, we obtain the desired tail
probability.

##  4\. Proof of Theorem 1.2: upper bound for $\zeta\neq\pm 1$

We want to show that there exists a constant $C$ such that with probability at
least $1-\varepsilon$, there is at least one root of $f_{\rho,n}$ in the ball
$B(\zeta,2C/n)$. For a sufficiently small constant $\delta>0$ depending only
on $\rho$, we consider the strip $\zeta+\zeta U_{C}/n$ that goes along the
unit circle where $U_{C}=(-\delta,\delta)\times(-C,C)\subset B(0,2C)$. Since
this strip is a subset of $B(\zeta,2C/n)$, it suffices to show that with
probability at least $1-\varepsilon$, there is at least one root of
$f_{\rho,n}$ in $U_{C}$. The use of $U_{C}$ in place of the ball allows us to
derive the upper bound using much simpler arguments because for
$z_{1},z_{2}\in U_{C}$, they are either very close or
$|z_{1}-z_{2}|\approx|y_{1}-y_{2}|$ where $y_{i}$ is the imaginary part of
$z_{i}$ and is a real number!

Since the upper bound at $\pm 1$ has been proved in Theorem 1.1, it suffices
to assume that $\zeta\neq\pm 1$.

###  4.1. The setup

Consider the following rescaled version of $f_{n}$, centered around $\zeta$

$g_{n}(z)=\frac{1}{n^{\rho+1/2}}f_{\rho,n}\left(\zeta+\frac{1}{n}\zeta
z\right),z\in U_{C}.$ (3)

The proof consists of the following steps.

  1. (1)

Construct a Gaussian process $g_{\infty}$ that shall be the limit of $g_{n}$.

  2. (2)

Show that for all $\varepsilon\geq 0$, there exists $C$ such that

${\mathbf{P}}(N_{g_{\infty}}(U_{C})=0)\leq\varepsilon.$ (4)

  3. (3)

When the random variables $\xi_{i}$ are iid standard Gaussian, show that on
$U_{C}$,

$g_{n}\xrightarrow{w}g.$ (5)

  4. (4)

Show that

${\mathbf{E}}N^{k}_{g_{n}}(U_{C})\to{\mathbf{E}}N^{k}_{g_{\infty}}(U_{C})$ (6)

for general $\xi_{i}$ (not necessarily Gaussian).

  5. (5)

Show that this implies

${\mathbf{P}}(N_{g_{n}}(U_{C})=0)\to{\mathbf{P}}(N_{g_{\infty}}(U_{C})=0).$
(7)

These steps are carried out in Sections 4.2, 4.3, 4.5, 4.6, 4.7, respectively.

###  4.2. Construct $g_{\infty}$

We have for all $z,w\in U_{C}$,

$\displaystyle{\mathbf{E}}g_{n}(z)\overline{g_{n}(w)}=\frac{1}{n^{\rho+1}}\sum_{k=0}^{n}|\zeta|^{2k}a_{k,\rho,n}^{2}(1+\frac{1}{n}z)^{k}(1+\frac{1}{n}\bar{w})^{k}$
$\displaystyle=$
$\displaystyle(1+o_{n}(1))\frac{1}{n}(1+\frac{1}{n}z)^{\rho}(1+\frac{1}{n}\bar{w})^{\rho}\frac{\partial^{2\rho}}{\partial
z^{\rho}\partial\bar{w}^{\rho}}\sum_{k=0}^{n}(1+\frac{1}{n}z)^{k}(1+\frac{1}{n}\bar{w})^{k}$
$\displaystyle=$ $\displaystyle(1+o_{n}(1))\frac{\partial^{2\rho}}{\partial
z^{\rho}\partial\bar{w}^{\rho}}\frac{(1+\frac{1}{n}z)^{n+1}(1+\frac{1}{n}\bar{w})^{n+1}-1}{n\left((1+\frac{1}{n}z)(1+\frac{1}{n}\bar{w})-1\right)}.$

where we used $|\zeta|=1$. When $\rho=0$, we have

$\displaystyle{\mathbf{E}}g_{n}(z)\overline{g_{n}(w)}$
$\displaystyle\xrightarrow{n\to\infty}$
$\displaystyle\frac{1}{z+\bar{w}}(\exp(z+\bar{w})-1)=\int_{0}^{1}e^{(z+\bar{w})t}dt=F(z+\bar{w}).$

where

$F(u):=\int_{0}^{1}e^{tu}dt=\frac{e^{u}-1}{u}.$ (8)

Similarly, for all $\rho$, we have

$\displaystyle{\mathbf{E}}g_{n}(z)\overline{g_{n}(w)}$
$\displaystyle\xrightarrow{n\to\infty}$
$\displaystyle\frac{\partial^{2\rho}}{\partial
z^{\rho}\partial\bar{w}^{\rho}}F(z+\bar{w})=F^{(2\rho)}(z+\bar{w})=\int_{0}^{1}t^{2\rho}e^{t(z+\bar{w})}dt.$

Moreover, since $\zeta\neq\pm 1$,

$\displaystyle{\mathbf{E}}g_{n}(z){g_{n}(w)}$ $\displaystyle=$
$\displaystyle(1+o(1))\frac{\partial^{2\rho}}{\partial z^{\rho}\partial
w^{\rho}}\frac{\zeta^{2n+2}(1+\frac{1}{n}z)^{n+1}(1+\frac{1}{n}w)^{n+1}-1}{n\left(\zeta^{2}(1+\frac{1}{n}z)(1+\frac{1}{n}w)-1\right)}\xrightarrow{n\to\infty}0$

as the denominator blows up (with its derivatives bounded) and the numerator
is bounded. So, it is logical to define the tentative limit $g_{\infty}$ to be

$g_{\infty}(z)=\int_{0}^{1}t^{\rho}e^{zt}dB(t)$

where $B$ is the standard complex Brownian motion
$B(t)=\frac{1}{\sqrt{2}}(B_{1}(t)+\sqrt{-1}B_{2}(t))$ with $B_{1},B_{2}$ being
independent standard real Brownian motions. We shall prove in Section 4.5 that
$g_{n}$ indeed converges to $g_{\infty}$ when $\zeta\neq\pm 1$.

###  4.3. Upper bound the hole radius for $g_{\infty}$

In this section, we want to show that there exists at least one root of
$g_{\infty}$ in $U_{C}$ with high probability (as $C\to\infty$). By
Chebyshev’s inequality, it suffices to show the following

$\mathbf{Var}N_{g_{\infty}}(U_{C})=o_{C\to\infty}({\mathbf{E}}N_{g_{\infty}}(U_{C}))^{2}.$
(9)

Indeed, we have

${\mathbf{P}}(N_{g_{\infty}}(U_{C})=0)\leq\frac{\mathbf{Var}N_{g_{\infty}}(U_{C})}{({\mathbf{E}}N_{g_{\infty}}(U_{C}))^{2}}=o(1).$

So, by choosing $C$ to be sufficiently large, the probability of $g_{\infty}$
having no roots can be arbitrarily small.

To prove (9), let us evaluate ${\mathbf{E}}N_{g_{\infty}}(U_{C})$. Since
$g_{\infty}$ is a Gaussian analytic function ([7]), we can use Kac-Rice
formula for Gaussian case. By [7, Corollary 3.4.2], we have

$\displaystyle{\mathbf{E}}N_{g_{\infty}}(U_{C})$ $\displaystyle=$
$\displaystyle\int_{U_{C}}\rho_{1}(z):=\int_{U_{C}}\frac{V-T^{2}S^{-1}}{\pi
S}dz=\int_{U_{C}}\frac{VS-T^{2}}{\pi S^{2}}dz$

where

$S(u)={\mathbf{E}}|g_{\infty}(z)|^{2},\quad
T(u)={\mathbf{E}}g_{\infty}^{\prime}(z)\overline{g_{\infty}(z)},\quad
V(u)={\mathbf{E}}|g_{\infty}^{\prime}(z)|^{2}.$

Let $u=z+\bar{z}\in{\mathbb{R}}$. By the definition (8) of $F$, we have

$S=F^{(2\rho)}(u).$

Taking derivative gives

$T=F^{(2\rho+1)}(u),V=F^{(2\rho+2)}(u).$

Since

$F(u)=\frac{e^{u}-1}{u}=\sum_{k=0}^{\infty}\frac{u^{k}}{(k+1)!},$

we have

$F^{(k)}(0)=\frac{1}{k+1}.$

Therefore,

$F^{(\rho+2)}(0)F^{(\rho)}(0)-\left(F^{(\rho+1)}(0)\right)^{2}\neq 0.$ (10)

And so,
$V(0)S(0)-T^{2}(0)=\frac{1}{(2\rho+1)(2\rho+3)}-\frac{1}{(2\rho+2)^{2}}>0.$ By
choosing $\delta$ to be sufficiently small, it holds that for all $|u|\leq
2\delta$, we have

$V(u)S(u)-T^{2}(u)=\Theta_{\rho}(1)\quad\mbox{and}\quad
S(u)=\Theta_{\rho}(1).$ (11)

By the definition of $U_{C}$, $u=z+\bar{z}\in[-2\delta,2\delta]$ for all $z\in
U_{C}$. Hence, $\rho_{1}(z)=\Theta_{\rho}(1)$ which implies

$\displaystyle{\mathbf{E}}N_{g_{\infty}}(U_{C})=\Theta(C).$

It remains to show the following.

######  Lemma 4.1.

We have

$\displaystyle\mathbf{Var}N_{g_{\infty}}(U_{C})$ $\displaystyle\ll$
$\displaystyle C\log C.$ (12)

###### Proof of Lemma 4.1.

Let $N=N_{g_{\infty}}(U_{C})$. We have

$\displaystyle\mathbf{Var}(N_{g_{\infty}}(U_{C}))$ $\displaystyle=$
$\displaystyle{\mathbf{E}}N(N-1)+{\mathbf{E}}N={\mathbf{E}}N(N-1)+O(C).$

So, it suffices to show that ${\mathbf{E}}N(N-1)\ll C\log C$. By [7, Corollary
3.4.2],

$\displaystyle{\mathbf{E}}N(N-1)$ $\displaystyle=$
$\displaystyle\int_{U_{C}}\int_{U_{C}}\rho_{2}(z_{1},z_{2})-\rho_{1}(z_{1})\rho_{1}(z_{2})dz_{1}dz_{2}$
(13)

where

$\displaystyle\pi^{2}\rho_{2}(z_{1},z_{2})=\frac{\text{per}(\mathcal{V}-\mathcal{T}\mathcal{S}^{-1}\mathcal{T}^{*})}{\det(\mathcal{S})},$
(14)

with $\mathcal{S},\mathcal{T},\mathcal{V}$ being $2\times 2$ matrices defined
by

$\mathcal{S}_{i,j}={\mathbf{E}}g(z_{i})\bar{g}(z_{j})=\int_{0}^{1}t^{2\rho}e^{tu_{ij}}dt=F^{(2\rho)}(u_{ij}),\quad
u_{ij}=z_{i}+\bar{z}_{j}$
$\mathcal{T}_{ij}={\mathbf{E}}g^{\prime}(z_{i})\bar{g}(z_{j})=\int_{0}^{1}t^{2\rho+1}e^{tu_{ij}}dt=F^{(2\rho+1)}(u_{ij})$
$V_{ij}={\mathbf{E}}g^{\prime}(z_{i})\overline{g^{\prime}}(z_{j})=\int_{0}^{1}t^{2\rho+2}e^{tu_{ij}}dt=F^{(2\rho+2)}(u_{ij}).$

Since $\mathcal{T}^{*}=\mathcal{T}$, by letting
$\mathcal{W}=\mathcal{T}\mathcal{S}^{-1}\mathcal{T}^{*}$, we have

$\displaystyle\mathcal{W}_{ij}$ $\displaystyle=$
$\displaystyle=\sum_{k=1}^{2}\sum_{\ell=1}^{2}\mathcal{T}_{ik}(\mathcal{S}^{-1})_{k\ell}\mathcal{T}_{\ell
j}$ (15)

where

$\det(\mathcal{S})\mathcal{S}^{-1}=\begin{pmatrix}\mathcal{S}_{22}&-\mathcal{S}_{12}\\\
-\mathcal{S}_{21}&\mathcal{S}_{11}\end{pmatrix}.$

We present a straightforward observation from the forms of $\rho_{1}$ and
$\rho_{2}$.

######  Lemma 4.2.

If $\mathcal{S},\mathcal{T},\mathcal{V},\mathcal{W}$ were diagonal (namely,
setting the off-diagonal entries to $0$), then
$\rho_{2}(z_{1},z_{2})-\rho_{1}(z_{1})\rho_{1}(z_{2})=0.$ .

Motivated by this, we will show that for $z_{1}$ and $z_{2}$ far away, the
above matrices are indeed diagonally dominated, and hence
$\rho_{2}(z_{1},z_{2})-\rho_{1}(z_{1})\rho_{1}(z_{2})$ is small. In
particular, we show the following.

######  Lemma 4.3 (Off diagonal).

For all $D_{0}$ satisfying $10\leq D_{0}\leq C$, let
$\mathcal{D}_{C,D_{0}}=\\{(z_{1},z_{2})\in U_{C}^{2}:|z_{1}-z_{2}|\geq
D_{0}\\}$. For all $(z_{1},z_{2})\in U_{C}$, it holds that

$|\rho_{2}(z_{1},z_{2})-\rho_{1}(z_{1})\rho_{1}(z_{2})|\ll|z_{1}-z_{2}|^{-1}.$

In particular, we have

$\iint_{\mathcal{D}_{C,D_{0}}}\rho_{2}(z_{1},z_{2})-\rho_{1}(z_{1})\rho_{1}(z_{2})dz_{1}dz_{2}\ll
C\log(C/D_{0}).$

When $z_{1}$ and $z_{2}$ are close, we show that
$\rho_{2}(z_{1},z_{2})-\rho_{1}(z_{1})\rho_{1}(z_{2})$ is bounded and hence
the contribution from the diagonal region is negligible.

######  Lemma 4.4 (Diagonal).

There exists a constant $M$, independent of $C$, such that for all
$(z_{1},z_{2})\in U_{C}^{2}$, we have

$|\rho_{2}(z_{1},z_{2})-\rho_{1}(z_{1})\rho_{1}(z_{2})|\ll M.$

This implies

$\iint_{U_{C}^{2}\setminus\mathcal{D}_{C,D_{0}}}\left(\rho_{2}(z_{1},z_{2})-\rho_{1}(z_{1})\rho_{1}(z_{2})\right)dz_{1}dz_{2}\ll
CD_{0}.$

Assuming these lemmas, letting $D_{0}=\log C$, we get that
$\mathbf{Var}N(U_{C})\ll C\log C$ as desired. This finishes the proof of Lemma
4.1. ∎

###### Proof of Lemma 4.2.

We write $x_{i}={\operatorname{Re}}(z_{i})$ and
$y_{i}={\operatorname{Im}}(z_{i})$ for $i=1,2$. Note that $|x_{i}|\leq\delta$
for all $i$. By (11),

$|\mathcal{S}_{11}|=\Theta(1),|\mathcal{S}_{22}|=\Theta(1).$ (16)

Let $|z_{1}-z_{2}|=D$, we have $D\geq D_{0}\geq 10$. So, $|y_{1}-y_{2}|\gg D$
and hence

$|u_{12}|=|u_{21}|\gg D.$ (17)

For all $u\in{\mathbb{C}}$, since $(uF(u))^{(k)}=kF^{(k-1)}+uF^{(k)}$ and
since the left-hand side equals $e^{u}$ for all $k$, we get

$F^{(k)}(u)=\frac{e^{u}-kF^{(k-1)}(u)}{u}.$

For $u_{12}=z_{1}+\bar{z}_{2}$, we have $|{\operatorname{Re}}(u_{12})|\leq
2\delta$ and so,

$|e^{u_{12}}|=e^{{\operatorname{Re}}(u_{12})}=O(1).$

Hence, by (17) and induction in $k$, it holds for all $k\leq\rho$ that

$|F^{(k)}(u_{12})|\ll\frac{1}{|u_{12}|}$

which gives

$|\mathcal{S}_{12}|=O(D^{-1}).$ (18)

Similarly, $|\mathcal{S}_{21}|=O(D^{-1})$.

Thus,

$\det(\mathcal{S})=\mathcal{S}_{11}\mathcal{S}_{22}-O(D^{-2})=(1+O(D^{-2}))\mathcal{S}_{11}\mathcal{S}_{22}=\Theta(1).$

And so,

$\mathcal{S}^{-1}=\begin{pmatrix}\frac{1+O(D^{-2})}{\mathcal{S}_{11}}&O(D^{-1})\\\
O(D^{-1})&\frac{1+O(D^{-2})}{\mathcal{S}_{22}}\end{pmatrix}.$

Similarly, the same bounds as in (16) and (18) hold for $\mathcal{T}$ and
$\mathcal{V}$ in place of $\mathcal{S}$.

Using these bounds, we get that $\mathcal{W}_{12}$ is the sum of 4 terms each
of which is of order $O(D^{-1})$. So, $|\mathcal{W}_{12}|=O(D^{-1})$.
Likewise, $|W_{21}|=O(D^{-1})$.

Similarly,

$\mathcal{W}_{11}=O(D^{-1})+\mathcal{T}_{11}^{2}(\mathcal{S}^{-1})_{11}=O(D^{-1})+\frac{(1+O(D^{-2}))\mathcal{T}_{11}^{2}}{\mathcal{S}_{11}}=\frac{\mathcal{T}_{11}^{2}}{\mathcal{S}_{11}}+O(D^{-1})=\Theta(1)$

and

$\mathcal{W}_{22}=\frac{\mathcal{T}_{22}^{2}}{\mathcal{S}_{22}}+O(D^{-1})=\Theta(1).$

Therefore,

$|(\mathcal{V}-\mathcal{W})_{12}|=O(D^{-1}),|(\mathcal{V}-\mathcal{W})_{12}|=O(D^{-1}).$

And

$(\mathcal{V}-\mathcal{W})_{11}=\mathcal{V}_{11}-\frac{\mathcal{T}_{11}^{2}}{\mathcal{S}_{11}}+O(D^{-1})=\Theta(1)$

where in the last equality, we used (11). So,

$\text{per}(\mathcal{V}-\mathcal{W})=\left(\mathcal{V}_{11}-\frac{\mathcal{T}_{11}^{2}}{\mathcal{S}_{11}}\right)\left(\mathcal{V}_{22}-\frac{\mathcal{T}_{22}^{2}}{\mathcal{S}_{22}}\right)+O(D^{-1}).$

All in all, we get

$\displaystyle\pi^{2}\left|\rho_{2}(z_{1},z_{2})-\rho_{1}(z_{1})\rho_{1}(z_{2})\right|$
$\displaystyle=$
$\displaystyle\frac{\left(V_{11}-\frac{T_{11}^{2}}{\mathcal{S}_{11}}\right)\left(V_{22}-\frac{T_{22}^{2}}{\mathcal{S}_{22}}\right)+O(D^{-1})}{(1+O(D^{-2}))\mathcal{S}_{11}\mathcal{S}_{22}}-\frac{\left(V_{11}-\frac{T_{11}^{2}}{\mathcal{S}_{11}}\right)\left(V_{22}-\frac{T_{22}^{2}}{\mathcal{S}_{22}}\right)}{\mathcal{S}_{11}\mathcal{S}_{22}}$
$\displaystyle=$ $\displaystyle\frac{O(D^{-1})}{\Theta(1)}=O(D^{-1}).$

Integrating this over $\mathcal{D}_{C,D_{0}}$ we get

$\displaystyle\iint_{\mathcal{D}_{C,D_{0}}}\rho_{2}(z_{1},z_{2})-\rho_{1}(z_{1})\rho_{1}(z_{2})dz_{1}dz_{2}$
$\displaystyle\ll$
$\displaystyle\iint_{\mathcal{D}_{C,D_{0}}}|z_{1}-z_{2}|^{-1}dz_{1}dz_{2}$
$\displaystyle\ll$ $\displaystyle\int_{U_{C}}dz_{1}\int_{w\in[-2\delta,2\delta
2]\times[-2C,2C],|w|\geq D_{0}}|w|^{-1}dw$ $\displaystyle=$
$\displaystyle\Theta(C\int_{D_{0}}^{C}y^{-1}dy)=\Theta(C\log(C/D_{0})).$

This finishes the proof of Lemma 4.2. ∎

###### Proof of Lemma 4.4.

Since $\rho_{1}(z_{1})\rho_{1}(z_{2})$ is bounded over $U_{C}^{2}$, we only
need to show the boundedness of $\rho_{2}$. By the first part of Lemma 4.2, we
can reduce to the region

$U_{\text{diag}}:=\\{(z_{1},z_{2})\in U_{C}^{2}:|z_{1}-z_{2}|\ll 1\\}.$

Since $\rho_{2}$ can be written as a function of $(u_{ij})_{i,j=1,2}$, it is
also a function of $x_{1},x_{2}$ and $\Delta_{y}:=y_{1}-y_{2}$ where we recall
$x_{i}={\operatorname{Re}}(z_{i})$ and $y_{i}={\operatorname{Im}}(z_{i})$.
Note that $U_{\text{diag}}$ is a subset of
$\\{(z_{1},z_{2}):|x_{1}|\leq\delta,|x_{2}|\leq\delta,\Delta_{y}\ll 1\\}$
which is a compact set. So, if we can show that $\rho_{2}$ is in fact a
continuous function of $x_{1},x_{2}$ and $\Delta_{y}:=y_{1}-y_{2}$, we
conclude that it is bounded $U_{\text{diag}}$. To show continuity, note that
the only possible singularities of $\rho_{2}$ occur when
$\det(\mathcal{S})=0$. Thus, it suffices to show the following.

######  Lemma 4.5.

If $\det(\mathcal{S}(z_{1},z_{2}))=0$ then $z_{1}=z_{2}$.

######  Lemma 4.6.

For all $z\in U_{C}$, $\rho_{2}(z_{1},z_{2})$ is continuous at
$(z_{1},z_{2})=(z,z)$.

Proving these lemmas will complete the proof of Lemma 4.4. ∎

###### Proof of Lemma 4.5.

Assume that $\det(\mathcal{S}(z_{1},z_{2}))=0$. Since
$\mathcal{S}(z_{1},z_{2})$ is a $2\times 2$ (complex) matrix, there exist
deterministic complex numbers $w_{1},w_{2}$ such that

$(w_{1}\ w_{2})\mathcal{S}(w_{1}\ w_{2})^{\mathcal{T}}=0.$

In other words,

${\mathbf{E}}|w_{1}g(z_{1})+w_{2}g(z_{2})|^{2}=0.$

Since the left-hand side equals
$\int_{0}^{1}|w_{1}e^{tz_{1}}+w_{2}e^{tz_{2}}|^{2}dt$, we conclude that the
integrand is 0 for almost all $t$ (and hence for all $t$ by continuity).
Therefore, it is necessary that $z_{1}=z_{2}$. ∎

###### Proof of Lemma 4.6.

We need to show that for all $z\in{\mathbb{C}}$ (or jut $U_{C}$ if needed),

$\lim_{(\varepsilon,\delta)\to(0,0)}\rho_{2}(z,z+\varepsilon+i\delta)\text{
exists}.$

We shall perform Taylor expansion to the order 2 of the functions appearing in
(14). Here $z_{1}=z=x+iy,z_{2}=z+\varepsilon+i\delta$. Then

$u_{11}=2x=:u,u_{12}=(2x+\varepsilon)-i\delta,u_{22}=2x+2\varepsilon.$

Let

$a=F^{(2\rho)}(u),b=F^{(2\rho+1)}(u),c=F^{(2\rho+2)}(u),d=F^{(2\rho+3)}(u),e=F^{(2\rho+4)}(u).$

So,

$\mathcal{S}_{11}=F^{(2\rho)}(u)=a,$
$\mathcal{S}_{22}=F^{(2\rho)}(u+2\varepsilon)\approx a+2\varepsilon
b+2\varepsilon^{2}c,$
$\mathcal{S}_{12}=\overline{\mathcal{S}_{21}}=F^{(2\rho)}(u+\varepsilon-i\delta)\approx
a+(\varepsilon-i\delta)b+\frac{(\varepsilon-i\delta)^{2}}{2}c=\left(a+\varepsilon
b+\frac{\varepsilon^{2}-\delta^{2}}{2}c\right)-i\left(\delta
b+\varepsilon\delta c\right).$

That is

$\mathcal{S}\approx\begin{pmatrix}a&\left(a+\varepsilon
b+\frac{\varepsilon^{2}-\delta^{2}}{2}c\right)-i\left(\delta
b+\varepsilon\delta c\right)\\\ \left(a+\varepsilon
b+\frac{\varepsilon^{2}-\delta^{2}}{2}c\right)+i\left(\delta
b+\varepsilon\delta c\right)&a+2\varepsilon b+2\varepsilon^{2}c\end{pmatrix}.$
(19)

The denominator of $\rho_{2}$ is

$\displaystyle\det(\mathcal{S})$ $\displaystyle=$
$\displaystyle\mathcal{S}_{11}\mathcal{S}_{22}-\mathcal{S}_{12}\mathcal{S}_{21}\approx
a(a+2\varepsilon b+2\varepsilon^{2}c)-\left(a+\varepsilon
b+\frac{\varepsilon^{2}-\delta^{2}}{2}c\right)^{2}-\left(\delta
b+\varepsilon\delta c\right)^{2}$ $\displaystyle=$ $\displaystyle
2\varepsilon^{2}ac-\varepsilon^{2}b^{2}-(\varepsilon^{2}-\delta^{2})ac-\delta^{2}b^{2}+o(\varepsilon^{2}+\delta^{2})\approx(\varepsilon^{2}+\delta^{2})(ac-b^{2}).$

Note that by (11), $ac-b^{2}\neq 0$.

Since $\mathcal{T}$ and $\mathcal{V}$ are similar to $\mathcal{S}$, we get

$\mathcal{T}\approx\begin{pmatrix}b&\left(b+\varepsilon
c+\frac{\varepsilon^{2}-\delta^{2}}{2}d\right)-i\left(\delta
c+\varepsilon\delta d\right)\\\ \left(b+\varepsilon
c+\frac{\varepsilon^{2}-\delta^{2}}{2}d\right)+i\left(\delta
c+\varepsilon\delta d\right)&b+2\varepsilon c+2\varepsilon^{2}d\end{pmatrix}$

and

$\mathcal{V}\approx\begin{pmatrix}c&\left(c+\varepsilon
d+\frac{\varepsilon^{2}-\delta^{2}}{2}e\right)-i\left(\delta
d+\varepsilon\delta e\right)\\\ \left(c+\varepsilon
d+\frac{\varepsilon^{2}-\delta^{2}}{2}e\right)+i\left(\delta
d+\varepsilon\delta e\right)&c+2\varepsilon d+2\varepsilon^{2}e\end{pmatrix}.$

Note that the numerator of (14) for $\rho_{2}$ equals

$\text{per}(\mathcal{V}-\mathcal{W})=\frac{1}{\det(\mathcal{S})^{2}}\text{per}((\det\mathcal{S})\mathcal{V}-(\det\mathcal{S})\mathcal{W}).$

So far, (14) becomes

$\pi^{2}\rho_{2}\approx\frac{\text{per}((\det\mathcal{S})\mathcal{V}-(\det\mathcal{S})\mathcal{W})}{(\varepsilon^{2}+\delta^{2})^{3}(ac-b^{2})^{3}}=:\frac{\text{per}(\mathcal{Y})}{(\varepsilon^{2}+\delta^{2})^{3}(ac-b^{2})^{3}}.$

Next, we write down $(\det\mathcal{S})\mathcal{W}$. We have

$(\det\mathcal{S})\mathcal{W}_{11}=\mathcal{T}_{11}\mathcal{S}_{22}\mathcal{T}_{11}-\mathcal{T}_{12}\mathcal{S}_{21}\mathcal{T}_{11}-\mathcal{T}_{11}\mathcal{S}_{12}\mathcal{T}_{21}+\mathcal{T}_{12}\mathcal{S}_{11}\mathcal{T}_{21}$
$(\det\mathcal{S})\mathcal{W}_{12}=\mathcal{T}_{11}\mathcal{S}_{22}\mathcal{T}_{12}-\mathcal{T}_{12}\mathcal{S}_{21}\mathcal{T}_{12}-\mathcal{T}_{11}\mathcal{S}_{12}\mathcal{T}_{22}+\mathcal{T}_{12}\mathcal{S}_{11}\mathcal{T}_{22}$
$(\det\mathcal{S})\mathcal{W}_{21}=\mathcal{T}_{21}\mathcal{S}_{22}\mathcal{T}_{11}-\mathcal{T}_{22}\mathcal{S}_{21}\mathcal{T}_{11}-\mathcal{T}_{21}\mathcal{S}_{12}\mathcal{T}_{21}+\mathcal{T}_{22}\mathcal{S}_{11}\mathcal{T}_{21}$
$(\det\mathcal{S})\mathcal{W}_{22}=\mathcal{T}_{21}\mathcal{S}_{22}\mathcal{T}_{12}-\mathcal{T}_{22}\mathcal{S}_{21}\mathcal{T}_{12}-\mathcal{T}_{21}\mathcal{S}_{12}\mathcal{T}_{22}+\mathcal{T}_{22}\mathcal{S}_{11}\mathcal{T}_{22}$

So,

$\displaystyle\mathcal{Y}_{11}$ $\displaystyle=$
$\displaystyle(\varepsilon^{2}+\delta^{2})(ac-b^{2})c-(\mathcal{T}_{11}\mathcal{S}_{22}\mathcal{T}_{11}-\mathcal{T}_{12}\mathcal{S}_{21}\mathcal{T}_{11}-\mathcal{T}_{11}\mathcal{S}_{12}\mathcal{T}_{21}+\mathcal{T}_{12}\mathcal{S}_{11}\mathcal{T}_{21})$
$\displaystyle=$
$\displaystyle{(\varepsilon^{2}+\delta^{2})(ac-b^{2})c}-b^{2}(a+2\varepsilon
b+2\varepsilon^{2}c)$
$\displaystyle{+2b{\operatorname{Re}}\left(\left(b+\varepsilon
c+\frac{\varepsilon^{2}-\delta^{2}}{2}d\right)+i\left(\delta
c+\varepsilon\delta d\right)\right)\left(\left(a+\varepsilon
b+\frac{\varepsilon^{2}-\delta^{2}}{2}c\right)-i\left(\delta
b+\varepsilon\delta c\right)\right)}$ $\displaystyle-a\left(b+\varepsilon
c+\frac{\varepsilon^{2}-\delta^{2}}{2}d\right)^{2}-a\left(\delta
c+\varepsilon\delta d\right)^{2}$ $\displaystyle=$
$\displaystyle{(\varepsilon^{2}+\delta^{2})(ac-b^{2})c}-b^{2}(a+2\varepsilon
b+2\varepsilon^{2}c)$ $\displaystyle{+2b\left(\left(b+\varepsilon
c+\frac{\varepsilon^{2}-\delta^{2}}{2}d\right)\left(a+\varepsilon
b+\frac{\varepsilon^{2}-\delta^{2}}{2}c\right)+\left(\delta
c+\varepsilon\delta d\right)\left(\delta b+\varepsilon\delta c\right)\right)}$
$\displaystyle-a\left(b+\varepsilon
c+\frac{\varepsilon^{2}-\delta^{2}}{2}d\right)^{2}-a\left(\delta
c+\varepsilon\delta d\right)^{2}.$

Grouping the terms with $\varepsilon$, $\delta$, $\varepsilon\delta$,
$\varepsilon^{2}$, $\delta^{2}$ and smaller order terms together, we get

$\displaystyle\mathcal{Y}_{11}$ $\displaystyle=$
$\displaystyle\varepsilon\left(-2b^{3}+2b^{3}+2abc-2abc\right)+\varepsilon^{2}\left((ac-b^{2})c+b^{2}c+abd-
ac^{2}-abd\right)$
$\displaystyle+\delta^{2}\left((ac-b^{2})c-b^{2}c-abd+2b^{2}c+abd-
ac^{2}\right)+O(\varepsilon^{3}+\varepsilon\delta^{2}+\varepsilon^{2}\delta+\delta^{3})$
$\displaystyle=$ $\displaystyle
O(\varepsilon^{3}+\varepsilon\delta^{2}+\varepsilon^{2}\delta+\delta^{3}).$

We now try to accomplish the same estimate for the other three entries of
$\mathcal{Y}$. We have

$\displaystyle\mathcal{Y}_{22}$ $\displaystyle=$
$\displaystyle(\varepsilon^{2}+\delta^{2})(ac-b^{2})\left(c+2\varepsilon
d+2\varepsilon^{2}e\right)-(\mathcal{T}_{22}\mathcal{S}_{11}\mathcal{T}_{22}-\mathcal{T}_{22}\mathcal{S}_{21}\mathcal{T}_{12}-\mathcal{T}_{21}\mathcal{S}_{12}\mathcal{T}_{22}+\mathcal{T}_{21}\mathcal{S}_{22}\mathcal{T}_{12})$
$\displaystyle=$
$\displaystyle{(\varepsilon^{2}+\delta^{2})(ac-b^{2})\left(c+2\varepsilon
d+2\varepsilon^{2}e\right)}-a\left(b+2\varepsilon
c+2\varepsilon^{2}d\right)^{2}$ $\displaystyle{+2\left(b+2\varepsilon
c+2\varepsilon^{2}d\right){\operatorname{Re}}\left(\left(\left(b+\varepsilon
c+\frac{\varepsilon^{2}-\delta^{2}}{2}d\right)+i\left(\delta
c+\varepsilon\delta d\right)\right)\right.}$ $\displaystyle{\hskip
122.34685pt\left.\left(\left(a+\varepsilon
b+\frac{\varepsilon^{2}-\delta^{2}}{2}c\right)-i\left(\delta
b+\varepsilon\delta c\right)\right)\right]}$
$\displaystyle-\left(b+2\varepsilon
c+2\varepsilon^{2}d\right)\left(b+\varepsilon
c+\frac{\varepsilon^{2}-\delta^{2}}{2}d\right)^{2}-\left(b+2\varepsilon
c+2\varepsilon^{2}d\right)\left(\delta c+\varepsilon\delta d\right)^{2}.$

So,

$\displaystyle\mathcal{Y}_{22}$ $\displaystyle=$
$\displaystyle{(\varepsilon^{2}+\delta^{2})(ac-b^{2})\left(c+2\varepsilon
d+2\varepsilon^{2}e\right)}-a\left(b+2\varepsilon
c+2\varepsilon^{2}d\right)^{2}$ $\displaystyle{+2\left(b+2\varepsilon
c+2\varepsilon^{2}d\right)\left(\left(b+\varepsilon
c+\frac{\varepsilon^{2}-\delta^{2}}{2}d\right)\left(a+\varepsilon
b+\frac{\varepsilon^{2}-\delta^{2}}{2}c\right)+\delta^{2}\left(c+\varepsilon
d\right)\left(b+\varepsilon c\right)\right)}$
$\displaystyle-{\left(a+2\varepsilon
b+2\varepsilon^{2}c\right)\left(b+\varepsilon
c+\frac{\varepsilon^{2}-\delta^{2}}{2}d\right)^{2}-\left(a+2\varepsilon
b+2\varepsilon^{2}c\right)\delta^{2}\left(c+\varepsilon d\right)^{2}}$

Comparing this with $\mathcal{Y}_{11}$, we get

$\displaystyle\mathcal{Y}_{22}$ $\displaystyle=$
$\displaystyle\mathcal{Y}_{11}+{(\varepsilon^{2}+\delta^{2})(ac-b^{2})\left(2\varepsilon
d+2\varepsilon^{2}e\right)}+(2\varepsilon
b^{3}+2\varepsilon^{2}b^{2}c-4\varepsilon^{2}ac^{2}-4\varepsilon^{2}abd-4\varepsilon
abc)$ $\displaystyle{+2\left(2\varepsilon
c+2\varepsilon^{2}d\right)\left(\left(b+\varepsilon
c+\frac{\varepsilon^{2}-\delta^{2}}{2}d\right)\left(a+\varepsilon
b+\frac{\varepsilon^{2}-\delta^{2}}{2}c\right)+\left(\delta
c+\varepsilon\delta d\right)\left(\delta b+\varepsilon\delta c\right)\right)}$
$\displaystyle-{\left(2\varepsilon
b+2\varepsilon^{2}c\right)\left(b+\varepsilon
c+\frac{\varepsilon^{2}-\delta^{2}}{2}d\right)^{2}-\left(2\varepsilon
b+2\varepsilon^{2}c\right)\delta^{2}\left(c+\varepsilon
d\right)^{2}}+O(\varepsilon^{3}+\varepsilon\delta^{2}+\varepsilon^{2}\delta+\delta^{3})$
$\displaystyle=$ $\displaystyle
O(\varepsilon^{3}+\varepsilon\delta^{2}+\varepsilon^{2}\delta+\delta^{3})+(2\varepsilon
b^{3}+2\varepsilon^{2}b^{2}c-4\varepsilon^{2}ac^{2}-4\varepsilon^{2}abd-4\varepsilon
abc)$ $\displaystyle{+2\left(2\varepsilon
c+2\varepsilon^{2}d\right)\left(b+\varepsilon c\right)\left(a+\varepsilon
b\right)}-{\left(2\varepsilon b+2\varepsilon^{2}c\right)\left(b+\varepsilon
c\right)^{2}}$

which gives

$\displaystyle\mathcal{Y}_{22}$ $\displaystyle=$ $\displaystyle
O(\varepsilon^{3}+\varepsilon\delta^{2}+\varepsilon^{2}\delta+\delta^{3})+(2\varepsilon
b^{3}+2\varepsilon^{2}b^{2}c-4\varepsilon^{2}ac^{2}-4\varepsilon^{2}abd-4\varepsilon
abc)$ $\displaystyle{+2\left(2\varepsilon
c+2\varepsilon^{2}d\right)\left(ab+\varepsilon b^{2}+\varepsilon
ac\right)}-{\left(2\varepsilon
b+2\varepsilon^{2}c\right)\left(b^{2}+2\varepsilon bc\right)}$
$\displaystyle=$ $\displaystyle
O(\varepsilon^{3}+\varepsilon\delta^{2}+\varepsilon^{2}\delta+\delta^{3})+(2\varepsilon
b^{3}+2\varepsilon^{2}b^{2}c-4\varepsilon^{2}ac^{2}-4\varepsilon^{2}abd-4\varepsilon
abc)$ $\displaystyle{+\left(4\varepsilon
abc+4\varepsilon^{2}b^{2}c+4\varepsilon^{2}ac^{2}+4\varepsilon^{2}abd\right)}-{\left(2\varepsilon
b^{3}+6\varepsilon^{2}b^{2}c\right)}$ $\displaystyle=$ $\displaystyle
O(\varepsilon^{3}+\varepsilon\delta^{2}+\varepsilon^{2}\delta+\delta^{3}).$

Finally,

$\displaystyle\mathcal{Y}_{12}$ $\displaystyle=$
$\displaystyle(\varepsilon^{2}+\delta^{2})(ac-b^{2})\left(\left(c+\varepsilon
d+\frac{\varepsilon^{2}-\delta^{2}}{2}e\right)-i\left(\delta
d+\varepsilon\delta e\right)\right)$
$\displaystyle-(T_{11}\mathcal{S}_{22}T_{12}-T_{12}\mathcal{S}_{21}T_{12}-T_{11}\mathcal{S}_{12}T_{22}+T_{12}\mathcal{S}_{11}T_{22})$
$\displaystyle=$ $\displaystyle
O(\varepsilon^{3}+\varepsilon\delta^{2}+\varepsilon^{2}\delta+\delta^{3})+(\varepsilon^{2}+\delta^{2})(ac-b^{2})c$
$\displaystyle{-b\left(a+2\varepsilon
b+2\varepsilon^{2}c\right)\left(\left(b+\varepsilon
c+\frac{\varepsilon^{2}-\delta^{2}}{2}d\right)-i\left(\delta
c+\varepsilon\delta d\right)\right)}$
$\displaystyle{+\left(\left(b+\varepsilon
c+\frac{\varepsilon^{2}-\delta^{2}}{2}d\right)-i\left(\delta
c+\varepsilon\delta d\right)\right)^{2}\left(\left(a+\varepsilon
b+\frac{\varepsilon^{2}-\delta^{2}}{2}c\right)+i\left(\delta
b+\varepsilon\delta c\right)\right)}$ $\displaystyle{+b\left(b+2\varepsilon
c+2\varepsilon^{2}d\right)\left(\left(a+\varepsilon
b+\frac{\varepsilon^{2}-\delta^{2}}{2}c\right)-i\left(\delta
b+\varepsilon\delta c\right)\right)}$ $\displaystyle-a\left(b+2\varepsilon
c+2\varepsilon^{2}d\right)\left(\left(b+\varepsilon
c+\frac{\varepsilon^{2}-\delta^{2}}{2}d\right)-i\left(\delta
c+\varepsilon\delta d\right)\right).$

So,

$\displaystyle\mathcal{Y}_{12}$ $\displaystyle=$ $\displaystyle
O(\varepsilon^{3}+\varepsilon\delta^{2}+\varepsilon^{2}\delta+\delta^{3})+(\varepsilon^{2}+\delta^{2})(ac-b^{2})c$
$\displaystyle{-ab\left(b+\varepsilon
c+\frac{\varepsilon^{2}-\delta^{2}}{2}d-i\left(\delta c+\varepsilon\delta
d\right)\right)-2\varepsilon b^{2}(b+\varepsilon c-i\delta
c)-2\varepsilon^{2}b^{2}c}$ $\displaystyle{+b\left(b+\varepsilon
c+\frac{\varepsilon^{2}-\delta^{2}}{2}d-i\left(\delta c+\varepsilon\delta
d\right)\right)\left(a+\varepsilon
b+\frac{\varepsilon^{2}-\delta^{2}}{2}c+i\left(\delta b+\varepsilon\delta
c\right)\right)}$ $\displaystyle{+\varepsilon c\left(\left(b+\varepsilon
c\right)-i\delta c\right)\left(a+\varepsilon b+i\delta
b\right)+\frac{\varepsilon^{2}-\delta^{2}}{2}dab}$ $\displaystyle{-i\delta
c\left(\left(b+\varepsilon c\right)-i\delta c\right)\left(a+\varepsilon
b+i\delta b\right)-i\varepsilon\delta abd}$
$\displaystyle{+b^{2}\left(a+\varepsilon
b+\frac{\varepsilon^{2}-\delta^{2}}{2}c-i\left(\delta b+\varepsilon\delta
c\right)\right)+2\varepsilon bc\left(a+\varepsilon b-i\delta
b\right)+2\varepsilon^{2}abd}$ $\displaystyle-ab\left(b+\varepsilon
c+\frac{\varepsilon^{2}-\delta^{2}}{2}d-i\left(\delta c+\varepsilon\delta
d\right)\right)-2\varepsilon ac\left(b+\varepsilon c-i\delta
c\right)+2\varepsilon^{2}abd$

giving

$\displaystyle\mathcal{Y}_{12}$ $\displaystyle=$ $\displaystyle
O(\varepsilon^{3}+\varepsilon\delta^{2}+\varepsilon^{2}\delta+\delta^{3})+(\varepsilon^{2}+\delta^{2})(ac-b^{2})c$
$\displaystyle{-ab\left(b+\varepsilon
c+\frac{\varepsilon^{2}-\delta^{2}}{2}d-i\left(\delta c+\varepsilon\delta
d\right)\right)-2\varepsilon b^{2}(b+\varepsilon c-i\delta
c)-2\varepsilon^{2}b^{2}c}$ $\displaystyle{+b^{2}\left(a+\varepsilon
b+\frac{\varepsilon^{2}-\delta^{2}}{2}c+i\left(\delta b+\varepsilon\delta
c\right)\right)+\varepsilon bc\left(a+\varepsilon b+i\delta
b\right)+\frac{\varepsilon^{2}-\delta^{2}}{2}abd}$ $\displaystyle{-i\delta
bc\left(a+\varepsilon b+i\delta b\right)-i\varepsilon\delta abd}$
$\displaystyle{+\varepsilon bc\left(a+\varepsilon b+i\delta
b\right)+\varepsilon^{2}c^{2}\left(a+\varepsilon b+i\delta
b\right)-i\varepsilon\delta ac^{2}+\frac{\varepsilon^{2}-\delta^{2}}{2}dab}$
$\displaystyle{-i\delta bc\left(a+\varepsilon b+i\delta
b\right)-i\varepsilon\delta ac^{2}-\delta^{2}ac^{2}-i\varepsilon\delta abd}$
$\displaystyle{+b^{2}\left(a+\varepsilon
b+\frac{\varepsilon^{2}-\delta^{2}}{2}c-i\left(\delta b+\varepsilon\delta
c\right)\right)+2\varepsilon bc\left(a+\varepsilon b-i\delta
b\right)+2\varepsilon^{2}abd}$ $\displaystyle-ab\left(b+\varepsilon
c+\frac{\varepsilon^{2}-\delta^{2}}{2}d-i\left(\delta c+\varepsilon\delta
d\right)\right)-2\varepsilon ac\left(b+\varepsilon c-i\delta
c\right)+2\varepsilon^{2}abd$ $\displaystyle=$ $\displaystyle
O(\varepsilon^{3}+\varepsilon\delta^{2}+\varepsilon^{2}\delta+\delta^{3}).$

Since the product of any two terms in
$\\{\varepsilon^{3},\varepsilon\delta^{2},\varepsilon^{2}\delta,\delta^{3}\\}$
is bounded by $(\varepsilon^{2}+\delta^{2})^{3}$, we yield the continuity of
$\rho_{2}$. ∎

###  4.4. More on $N_{g_{n}}(U_{C})$ and $N_{g_{\infty}}(U_{C})$

Before moving on the next section to show the convergence of $g_{n}$ to
$g_{\infty}$ and their number of roots, we will first show that
$N_{g_{n}}(U_{C})$ have uniformly bounded higher moments.

######  Lemma 4.7.

There exists a constant $A$ such that the following holds. For any $\ell\geq
0$, we have

${\mathbf{E}}N_{g_{n}}^{\ell}(U_{C})\leq(A\ell)^{\ell}.$

###### Proof.

Let $k\geq A\ell$ for some large constant $A$. We want to bound the
probability that $N_{g_{n}}(U_{C})\geq k$. We divide $U_{C}$ into
$O(C\eta^{-2})$ (possibly overlapping) open balls $B_{i}=B(c_{i},\eta)$
centered at $c_{i}$ of radius $\eta$, which is chosen to be sufficiently
small. Then there exists $i$ such that $B_{i}$ contains at least
$s=k\eta^{2}/C$ roots. Then by Hermite interpolation, as $g_{n}$ is analytic
with probability one, we have

$|g_{n}(c_{i})|\leq\frac{1}{s!}\eta^{s}\sup_{z\in
B(c_{i},\eta)}|g_{n}^{(s)}(z)|.$ (20)

By Taylor expanding $g_{n}^{s}(z)$ around $c_{i}$, we obtain for any $m\geq 0$
(we later choose $m=\log n$),

$\displaystyle|g_{n}^{(s)}(z)|\leq\sum_{j=s}^{s+m-1}\frac{|g_{n}^{(j)}(c_{i})|}{(j-s)!}\eta^{j-s}+\sup_{w\in
B(c_{i},\eta)}\frac{|g_{n}^{(s+m)}(w)|}{m!}\eta^{m}.$ (21)

For each $j$, $g_{n}^{(j)}(c_{i})$ is a Gaussian random variable with mean 0
and variance equals that of
$\frac{1}{n^{j+\rho+1/2}}f_{n}^{(j)}(\zeta+\frac{1}{n}\zeta c_{i})$, which is
of order

$(1+o(1))\frac{1}{n^{2j+2\rho+1}}\sum_{h=j}^{n}h^{2}(h-1)^{2}\dots(h-j+1)^{2}a_{h,\rho,n}^{2}|\zeta+\frac{1}{n}\zeta
c_{i}|^{2h-2j}=O_{C}(1)$

where we used $|\zeta+\frac{1}{n}\zeta c_{i}|^{j}\leq(1+2C/n)^{n}\leq
e^{2C}=O_{C}(1)$.
So, by Gaussianity, for all $M_{j}\geq 1$,

$\displaystyle{\mathbf{P}}(|g_{n}^{(j)}(c_{i})|\gg M_{j})\leq e^{-M_{j}}.$

Finally, for the supremum term, we observe

$\displaystyle\sup_{w\in B(c_{i},\eta)}|g_{n}^{(s+m)}(w)|$ $\displaystyle\ll$
$\displaystyle\frac{1}{n^{s+m+\rho+1/2}}\sum_{h=s+m}^{n}h(h-1)\dots(h-s-m+1)|a_{h,\rho,n}||\xi_{h}|(|c_{i}|+\eta)^{h}$
$\displaystyle\ll$ $\displaystyle\frac{1}{n^{1/2}}\sum_{h=1}^{n}|\xi_{h}|.$

Note that if we hadn’t used another round of Taylor expansion in (21) and just
applied the above bound to $|g_{n}^{(s)}(z)|$ and take supremum, the term
$n^{-1/2}\sum_{h=1}^{n}|\xi_{h}|$, which can be as large as $\sqrt{n}$, would
be too big to handle. Here, we performed (21) so that the extra term
$\eta^{m}/m!$ would swallow the $n^{-1/2}\sum_{h=1}^{n}|\xi_{h}|$. Indeed, for
an $M_{0}$ to be chosen,

${\mathbf{P}}(\sum_{h=1}^{n}|\xi_{h}|\geq M_{0}n)\leq
n{\mathbf{P}}(|\xi_{h}|\geq M_{0})\leq ne^{-M_{0}}.$

Thus,

${\mathbf{P}}(\sup_{w\in B(c_{i},\eta)}|g_{n}^{(s+m)}(w)|\gg n^{1/2}M_{0})\ll
ne^{-M_{0}}.$

Combining all of these events, we conclude that with probability at least
$1-ne^{-M_{0}}-\sum_{j=s}^{s+m-1}e^{-M_{j}}$, we have

$\sup_{z\in
B(c_{i},\eta)}|g_{n}^{(s)}(z)|\ll\sum_{j=s}^{s+m-1}\frac{M_{j}}{(j-s)!}\eta^{j-s}+\frac{n^{1/2}M_{0}\eta^{m}}{m!}.$

On this event,

$|g_{n}(c_{i})|\ll\frac{1}{s!}\eta^{s}\left(\sum_{j=s}^{s+m-1}\frac{M_{j}}{(j-s)!}\eta^{j-s}+\frac{n^{1/2}M_{0}\eta^{m}}{m!}\right)$

which only happens with probability at most

$O_{C}(1)\frac{1}{s!}\eta^{s}\left(\sum_{j=s}^{s+m-1}\frac{M_{j}}{(j-s)!}\eta^{j-s}+\frac{n^{1/2}M_{0}\eta^{m}}{m!}\right)$

since $g_{n}(c_{i})$ is a Gaussian random variable with variance
$\Theta_{C}(1)$. All together, we get that the probability that
$N_{g_{n}}^{\ell}(U_{C})\geq k$ is at most (up to a constant depending on
$C$),

$\frac{C}{\eta^{2}}\left[\frac{1}{s!}\eta^{s}\left(\sum_{j=s}^{s+m-1}\frac{M_{j}}{(j-s)!}\eta^{j-s}+\frac{n^{1/2}M_{0}\eta^{m}}{m!}\right)+\sum_{j=s}^{s+m-1}e^{-M_{j}}\right]+ne^{-M_{0}}$

for any choice of $\eta$, $M_{0},M_{j}$, with $s=k\eta^{2}/C$. For instance,
we choose $s=8\ell$, we get
$\eta=\frac{\sqrt{8C\ell}}{\sqrt{k}}\leq\sqrt{\frac{8C}{A}}$. By setting

$M_{j}=s\log\frac{1}{\eta}+\eta^{-(j-s)/2},\quad M_{0}=2\ell\log k+\log n,$

we obtain the tail probability of

$\frac{C}{\eta^{2}}\left[\frac{1}{s!}\eta^{s}\left(e^{\sqrt{\eta}}+\frac{n^{1/2}M_{0}\eta^{m}}{m!}\right)+\eta^{s}\right]+k^{-2\ell}.$

Sending $m\to\infty$, the term with $m$ goes to $0$, so we end up with

$\eta^{s-2}+k^{-2\ell}\ll(C\ell)^{4\ell}k^{-4\ell+1}+k^{-2\ell}.$

So,

$\displaystyle{\mathbf{E}}N_{g_{n}}^{\ell}(U_{C})$ $\displaystyle\ll$
$\displaystyle(A\ell)^{\ell}+\ell\sum_{k=A\ell}^{\infty}k^{\ell-1}{\mathbf{P}}(N_{g_{n}}(U_{C})\geq
k)$ $\displaystyle\ll$
$\displaystyle(A\ell)^{\ell}+\ell\sum_{k=A\ell}^{\infty}\left((C\ell)^{4\ell}k^{-3\ell}+k^{-\ell-1}\right)\ll(A\ell)^{\ell}$

as desired.

∎

###  4.5. Convergence of $g_{n}$ to $g_{\infty}$ when $\xi_{i}$ are iid
$\mathcal{N}(0,1)$

Now, we prove (5). We first start with two simple results for the Gaussian
models.

######  Lemma 4.8.

With probability one, $g_{n}$ and $g_{\infty}$ do not have double roots in
$U_{C}$.

###### Proof.

For $g_{n}$, if it has a double root then $f_{n}$ also has a double root. As
this is a polynomial of degree $n$, if $f_{n}(z)$ and $f_{n}^{\prime}(z)$ have
common roots then the resultant must have zero determinant. But the resultant
is a non-degenerate multivariate function of the Gaussians, so it is zero with
probability zero.

For $g_{\infty}(z)$, for any $\alpha>0$, we divide $U_{C}$ into
$O(C\alpha^{-2})$ balls $B_{i}$ of radius $\alpha$. We will show that the
probability there exists $i$ such that $N_{i}$, the number of zeros in
$B_{i}$, is greater than $2$ is of order $O(\alpha^{4})$, from which we see
that the given probability will be bounded by $O(\alpha^{2})$ after taking
union bound. Indeed, using the boundedness of $\rho_{2}$ in Lemma 4.4,

${\mathbf{P}}(N_{i}\geq 2)\leq{\mathbf{E}}(N_{i}(N_{i}-1))=\int_{B_{i}\times
B_{i}}\rho_{2}(z_{1},z_{2})dz_{1}dz_{2}\leq
O(|B_{i}|\times|B_{i}|)=O(\alpha^{4}).$

Sending $\alpha$ to $0$, we conclude that the probability that $g_{\infty}$
has double roots is 0. ∎

Our next simple result is the following.

######  Claim 4.9.

With probability one, $g_{n}(z)$ and $g_{\infty}$ do not have roots on the
boundary $\partial U_{C}$ of $U_{C}$.

###### Proof.

We will show for $g_{\infty}$ as the treatment for $g_{n}$ is similar. From
(11), we saw that for all $\alpha>0$ sufficiently small, $\rho_{1}(z)=O(1)$
for all $z\in U_{C}+B(0,\alpha)$. Let $N$ be the number of roots in $\partial
U_{C}+B(0,\alpha)$, then

${\mathbf{P}}(N\geq 1)\leq{\mathbf{E}}N=\int_{\partial
U_{C}+B(0,\alpha)}\rho_{1}(z)dz\leq O_{C}(\alpha).$

Sending $\alpha$ to $0$, we obtain the claim. ∎

Our treatment below is similar to [8, Section 4] where instead of real roots,
we consider complex roots. First, let ${\mathcal{H}}$ be the set of all
analytic function on the entire complex plane. We endow ${\mathcal{H}}$ with
the topology of uniform convergence on the compact sets, which can be
generated by the complete separable metric

$d(f,g)=\sum_{k}\frac{1}{2^{k}}\frac{\|f-g\|_{\bar{D}_{k}}}{1+\|f-g\|_{\bar{D}_{k}}},$

where $\bar{D}_{k}=\\{z\in{\mathbb{C}}:|z|\leq k\\}$ and $\|f\|_{K}=\sup_{z\in
K}|f(z)|$.

######  Lemma 4.10.

Let $A_{C}$ be the set of all $f\in{\mathcal{H}}$ which do not have multiple
roots in $U_{C}$ and do not have roots over the boundary of $U_{C}$. Then the
set $A_{C}$ is open.

###### Proof.

This follows from Hurwitz’s theorem. Indeed, consider a sequence
$(f_{n})_{n\in{\mathbb{N}}}$ in ${\mathcal{H}}$, which converges to some $f\in
A_{C}$ locally uniformly. We will show that $f_{n}\in A_{C}$ for sufficiently
large $n$. Let $R>0$ be large such that $U_{C}\subset D_{R}=\\{z:|z|<R\\}$.
Let $z_{1},\dots,z_{d}$ be the collections of all zeros of $f$ in $D_{R}$ with
multiplicities $m_{1},\dots,m_{d}$. Let $\alpha>0$ be sufficiently small such
that the open disks $z_{i}+D_{\alpha}$ are disjoint, and do not intersect the
boundary of the open sets $D_{R}$ and of $U_{C}$, except when $z_{i}$ are on
one of these boundaries. By Hurwitz’s theorem for sequence of (locally
convergent) analytic functions, there exists $n_{0}$ such that for all $n\geq
n_{0}$, $f_{n}$ has exactly $m_{k}$ zeros in $z_{k}+D_{\alpha}$. Now if
$z_{i}\in U_{C}$, then as $f\in A_{C}$, we must have $m_{i}=1$, and $f_{n}$
has exactly one zero in $z_{i}+D_{\alpha}$. Thus, $f_{n}\in A_{C}$ for all
$n\geq n_{0}$. ∎

######  Lemma 4.11.

The mapping $f\to\mathcal{Z}_{U_{C}}(f)=\\{z\in U_{C}:f(z)=0\\}$ to the space
of locally finite point measures on $U_{C}$ endowed with the vague topology is
continuous on $A_{C}$.

###### Proof.

This also follows from Hurwitz’s theorem with the same argument as in the
proof of Lemma 4.10, by letting the radius $\alpha$ tend to zero. ∎

######  Lemma 4.12.

We have the following weak convergence (of random elements with values in the
metric space ${\mathcal{H}}$)

$g_{n}\xrightarrow{w}g_{\infty}.$

###### Proof.

By Prokhorov’s theorem, it suffices to verify convergence in finite
dimensional and tightness. Let $z_{1},\dots,z_{k}$ be complex numbers. We
first observe that the convergence in distribution of the Gaussian vector
$(g_{n}(z_{1}),\dots,g_{n}(z_{k}))$ to the Gaussian vector
$(g_{\infty}(z_{1}),\dots,g_{\infty}(z_{k}))$ already follows from our
previous computations verifying the convergences of
${\mathbf{E}}g_{n}(z_{i})\overline{g_{n}(z_{j})}$ and
${\mathbf{E}}g_{n}(z_{i})g_{n}(z_{j})$ to
${\mathbf{E}}g_{\infty}(z_{i})\overline{g_{\infty}(z_{j})}$ and
${\mathbf{E}}g_{\infty}(z_{i})g_{\infty}(z_{j})$, respectively.

We need to verify tightness, for this, it suffices to show that for any $R>0$,
there exists $C_{R}<\infty$ such that

$\sup_{n}\sup_{|z|\leq R}{\mathbf{E}}|g_{n}(z)|^{2}<C_{R}.$

However, this is clear as

${\mathbf{E}}|g_{n}(z)|^{2}=(1+o_{n}(1))\frac{\partial^{2\rho}}{\partial
z^{\rho}\partial\bar{w}^{\rho}}\frac{(1+\frac{1}{n}z)^{n+1}(1+\frac{1}{n}\bar{w})^{n+1}-1}{n\left((1+\frac{1}{n}z)(1+\frac{1}{n}\bar{w})-1\right)}\Big{|}_{w=z}.$

∎

######  Theorem 4.13.

We have that $N_{g_{n}}(U_{C})\to N_{g_{\infty}}(U_{C})$ in distribution and
for each $k\in{\mathbb{N}}$,
$\lim_{n\to\infty}{\mathbf{E}}N_{g_{n}}^{k}(U_{C})={\mathbf{E}}N^{k}_{g_{\infty}}(U_{C})$.

###### Proof.

We have that $g_{n}\to g_{\infty}$ weakly, they are analytic and with
probability one, they all belong to $A_{C}$. By Lemma 4.11, the point process
$\mathcal{Z}_{U_{C}}(g_{n})$ converges to $\mathcal{Z}_{U_{C}}(g_{\infty})$
weakly, and hence the number of zeros $N_{g_{n}}(U_{C})$ converges to
$N_{g_{\infty}}(U_{C})$ in distribution. In particular, for all
$m\in{\mathbb{N}}$,
$p_{n,m}:={\mathbf{P}}(N_{g_{n}}(U_{C})=m)\to{\mathbf{P}}(N_{g_{\infty}}(U_{C})=m):=p_{m}$
as $n\to\infty$. By Fatou’s lemma and Lemma 4.7, it holds for all
$\ell\in{\mathbb{N}}$ that

${\mathbf{E}}N_{g_{\infty}}^{\ell}(U_{C})\leq\liminf_{n}{\mathbf{E}}N_{g_{n}}^{\ell}(U_{C})\leq(A\ell)^{\ell}.$
(22)

Fix $k\in{\mathbb{N}}$, we have for a large constant $M$,

$\displaystyle\left|{\mathbf{E}}N_{g_{n}}^{k}(U_{C})-{\mathbf{E}}N^{k}_{g_{\infty}}(U_{C})\right|$
$\displaystyle\leq$
$\displaystyle\sum_{m=0}^{M-1}m^{k}|p_{nm}-p_{m}|+{\mathbf{E}}N_{g_{n}}^{k}(U_{C})\textbf{1}_{N_{g_{n}}(U_{C})\geq
M}+{\mathbf{E}}N_{g_{\infty}}(U_{C})\textbf{1}_{N_{g_{\infty}}^{k}(U_{C})\geq
M}$ $\displaystyle\leq$
$\displaystyle\sum_{m=0}^{M-1}m^{k}|p_{nm}-p_{m}|+2\sup_{\hat{n}}\left({\mathbf{E}}N_{g_{\hat{n}}}^{2k}(U_{C})\right)^{1/2}{\mathbf{P}}\left(N_{g_{\hat{n}}}(U_{C})\geq
M\right)^{1/2}\text{ by Jensen's inequality}$ $\displaystyle\leq$
$\displaystyle\sum_{m=0}^{M-1}m^{k}|p_{nm}-p_{m}|+2(Ak)^{k}\sqrt{\frac{A}{M}}\text{
by \eqref{eq:moment} and Markov's inequality.}$

Letting $M$ and $n$ go to infinity, we obtain the convergence in moments. ∎

###  4.6. Convergence for the number of real roots

In this section, we prove (6). In other words, we prove the following
convergence of the number of roots $U_{C}$. Note that the random variables are
not necessarily Gaussian here.

The following generalizes Theorem 4.13 to non-Gaussian random variables.

######  Theorem 4.14.

Let $C$ be a fixed positive number. For all $k\geq 0$, we have

${\mathbf{E}}N^{k}_{g_{n}}(U_{C})\to{\mathbf{E}}N^{k}_{g_{\infty}}(U_{C})$

as $n\to\infty$.

Let $\tilde{g}_{n}$ be the version of $g_{n}$ when the random variables
$\xi_{i}$ are iid standard Gaussian. By Theorem 4.13, we have

${\mathbf{E}}N^{k}_{\tilde{g}_{n}}(U_{C})\to{\mathbf{E}}N^{k}_{g_{\infty}}(U_{C}).$
(23)

We note that the same proof holds with $U_{C}$ replaced by
$U_{C}+B(0,\alpha)$.

###### Proof.

Note that the number of roots of $g_{n}$ in $U_{C}$ is the same as the number
of roots in the original function $f_{n}$ in the set $\zeta+\frac{1}{n}\zeta
U_{C}$, by (3). Let $\varphi$ be a test function approximating the indicator
of $(U_{C})^{k}$, in particular, we let $\varphi$ be a smooth function such
that

$\displaystyle\textbf{1}_{(U_{C})^{k}}\leq\varphi\leq\textbf{1}_{(U_{C}+B(0,\alpha))^{k}}$
(24)

and $|\triangledown^{a}\varphi(z)|\ll 1$ for all multi-indices $a$ with
$0\leq|a|\leq 2k+4$.

By [6, Theorem 2.4] 111or perhaps a slightly readable [14, Theorem 4.3] which
was written for the Kac polynomial but it holds also for the derivatives of
the Kac polynomial. applied to the function $G=\varphi$ and the centers
$z_{1}=\dots=z_{k}=\zeta$, we get

$\displaystyle\left|{\mathbf{E}}\sum_{\zeta_{i_{1}},\dots,\zeta_{i_{k}}\in\mathcal{Z}(g_{n})}\varphi(n(\zeta_{i_{1}}/\zeta-1),\dots,n(\zeta_{i_{1}}/\zeta-1))\right.$
$\displaystyle\left.-{\mathbf{E}}\sum_{\zeta_{i_{1}},\dots,\zeta_{i_{k}}\in\mathcal{Z}(\tilde{g}_{n})}\varphi(n(\zeta_{i_{1}}/\zeta-1),\dots,n(\zeta_{i_{1}}/\zeta-1))\right|\ll
n^{-c}$

where $c>0$ is a small constant. Here, we note that the transformation
$z:=n(\zeta_{i}/\zeta-1)$ is just the inverse of the rescaling map
$\zeta_{i}=\zeta+\frac{1}{n}\zeta z$ that brings the neighborhood of $\zeta$
to $U_{C}$. We note that when $\varphi$ is replaced by
$\textbf{1}_{(U_{C})^{k}}$, the term under the expectation becomes
$N^{k}(U_{C})$. So, we have from (24) that

${\mathbf{E}}N^{k}_{g_{n}}(U_{C})\leq{\mathbf{E}}N^{k}_{\tilde{g}_{n}}(U_{C}+B(0,\alpha))+O_{\alpha}(n^{-c}).$

Using (23), we obtain

$\limsup_{n\to\infty}{\mathbf{E}}N^{k}_{g_{n}}(U_{C})\leq\limsup_{n\to\infty}{\mathbf{E}}N^{k}_{\tilde{g}_{n}}(U_{C}+B(0,\alpha))={\mathbf{E}}N^{k}_{g_{\infty}}(U_{C}+B(0,\alpha)).$

Sending $\alpha$ to 0, we obtain

$\limsup_{n\to\infty}{\mathbf{E}}N^{k}_{g_{n}}(U_{C})\leq\limsup_{\alpha\to
0}{\mathbf{E}}N^{k}_{g_{\infty}}(U_{C}+B(0,\alpha))={\mathbf{E}}N^{k}_{g_{\infty}}(U_{C})$

where the last equality follows from the dominated convergence theorem,
knowing that
${\mathbf{E}}N^{k}_{\tilde{g}_{\infty}}(U_{C}+B(0,\alpha))<\infty$ for some
$\alpha>0$ (by (22)). Similarly, we get the reverse direction and conclude the
proof. ∎

###  4.7. Upper bound the hole radius for $g_{n}$

In this section, we show the following theorem.

######  Theorem 4.15.

The random variables $N_{g_{n}}(U_{C})$ converges to $N_{g_{\infty}}(U_{C})$
in distribution as $n\to\infty$. In particular, we have (7):

${\mathbf{P}}(N_{g_{n}}(U_{C})=0)\to{\mathbf{P}}(N_{g_{\infty}}(U_{C})=0).$

Here, we recall that since the random variables $N_{g_{n}}(U_{C})$ are
discrete random variables supported on ${\mathbb{N}}$, convergence in
distribution means convergence of the probability density
${\mathbf{P}}(N_{g_{n}}(U_{C})=i)$, as $i$ varies.

###### Proof.

By Theorem 4.14, $N_{g_{n}(U_{C})}$ converges to $N_{g_{\infty}}(U_{C})$ in
moments. By 22 and the Carleman’s criteria, $N_{g_{\infty}}(U_{C})$ is
uniquely determined by its moments. Thus, we infer that $N_{g_{n}(U_{C})}$
converges to $N_{g_{\infty}}(U_{C})$ in distribution. ∎

##  5\. Proof of Theorem 1.2: Lower bound

We want to show that there exists $c$ such that

${\mathbf{E}}(N_{f_{n}}(B(\zeta,c/n)))\leq\varepsilon.$ (25)

Without loss of generality, we assume that $\varepsilon<1/100$ and
$c<\varepsilon$.

The first step is to reduce to the Gaussian case, via universality results.
Consider the Gaussian version of $f_{n}$,

$\tilde{f}_{n}=\tilde{f}_{\rho,n}=\sum_{i=0}^{n}a_{i,\rho,n}\tilde{\xi}_{i}z^{i}$

where $\tilde{\xi}_{i}$ are iid standard Gaussian.

Let $G$ be a smooth function such that approximates the indicator of the ball,
or more specifically, $\textbf{1}_{B(\zeta,c/n)}\leq
G\leq\textbf{1}_{B(\zeta,2c/n)}$ and
$||\triangledown^{a}G||_{\infty}=O(n^{a})$ for all $a\leq 3$. We now apply a
universality property of $f_{n}$ established in [6, Theorem 2.3]. This theorem
applied to the function $G$ states that the linear statistics
${\mathbf{E}}\sum_{w\in\mathcal{Z}(f_{n})}G(w)$ is universal, i.e.,

${\mathbf{E}}\sum_{w\in\mathcal{Z}(f_{n})}G(w)-{\mathbf{E}}\sum_{w\in\mathcal{Z}(\tilde{f}_{n})}G(w)\ll
n^{-\gamma}$

for a constant $\gamma$ independent of $n$ and $\zeta$.

Using this, we obtain

$\displaystyle{\mathbf{E}}(N_{f_{n}}(B(\zeta,c/n)))$ $\displaystyle\leq$
$\displaystyle{\mathbf{E}}\sum_{w\in\mathcal{Z}(f_{n})}G(w)={\mathbf{E}}\sum_{w\in\mathcal{Z}(\tilde{f}_{n})}G(w)+o(1)$
$\displaystyle\leq$
$\displaystyle{\mathbf{E}}(N_{\tilde{f}_{n}}(B(\zeta,2c/n)))+o(1).$

Thus, it suffices to prove that

${\mathbf{E}}(N_{\tilde{f}_{n}}(B(\zeta,2c/n)))\leq\varepsilon/2.$ (26)

In other words, it suffices to prove for the Gaussian case. To this end, we
let $B_{c}=B(0,2c)$ and definite the functions $g_{n}$ and $g_{\infty}$ as
before. We apply the Kac-Rice formula to $g_{\infty}$ to get

$\displaystyle{\mathbf{E}}N_{g_{\infty}}(B_{c})$ $\displaystyle=$
$\displaystyle\int_{B_{c}}\rho_{1}(z).$

By (11), we have for all $z\in B_{c}$, $\rho_{1}(z)\ll 1.$ Thus,

${\mathbf{E}}N_{g_{\infty}}(B_{c})\ll c^{2}\leq\varepsilon/4.$

By the same argument as for $U_{C}$ (noting that $B_{c}\subset U_{C}$ for
small $c$ and for $C\geq c$), we obtain the same limit as in Theorem 4.13. So,
we get

$\lim_{n\to\infty}{\mathbf{E}}N_{g_{n}}(B_{c})={\mathbf{E}}N_{g_{\infty}}(B_{c})\leq\varepsilon/4.$

So, by choosing $n$ to be sufficiently small, we obtain (26) as desired.

##  6\. Acknowledgment

We thank Manjunath Krishnapur for suggesting helpful references on correlation
functions.

## References

  * [1] John Baez, Dan Christensen, and Sam Derbyshire.  The beauty of roots.  https://math.ucr.edu/home/baez/roots/beauty_web.pdf.
  * [2] Peter Borwein and Christopher Pinner.  Polynomials with $\\{0,+1,-1\\}$ coefficients and a root close to a given point.  Canad. J. Math., 49(5):887–915, 1997.
  * [3] Jeremiah Buckley, Alon Nishry, Ron Peled, and Mikhail Sodin.  Hole probability for zeroes of gaussian taylor series with finite radii of convergence.  Probability Theory and Related Fields, 171:377–430, 2018.
  * [4] Danny Calegari, Sarah Koch, and Alden Walker.  Roots, Schottky semigroups, and a proof of Bandt’s conjecture.  Ergodic Theory Dynam. Systems, 37(8):2487–2555, 2017.
  * [5] Nicholas A Cook, Hoi H Nguyen, Oren Yakir, and Ofer Zeitouni.  Universality of poisson limits for moduli of roots of kac polynomials.  International mathematics research notices, 2023(8):6648–6690, 2023\.
  * [6] Yen Do, Oanh Nguyen, and Van Vu.  Roots of random polynomials with coefficients with polynomial growth.  Annals of Probability, 46(5):2407–2494, 2018.
  * [7] John Ben Hough, Manjunath Krishnapur, Yuval Peres, and Bálint Virág.  Zeros of Gaussian analytic functions and determinantal point processes, volume 51.  American Mathematical Society Providence, RI, 2009.
  * [8] Alexander Iksanov, Zakhar Kabluchko, and Alexander Marynych.  Local universality for real roots of random trigonometric polynomials.  Electron. J. Probab., 21:Paper No. 63, 19, 2016.
  * [9] Zakhar Kabluchko and Dmitry Zaporozhets.  Asymptotic distribution of complex zeros of random analytic functions.  The Annals of Probability, 42(4):1374–1395, 2014.
  * [10] Manjunath Krishnapur, Erik Lundberg, and Oanh Nguyen.  The number of limit cycles bifurcating from a randomly perturbed center.  arXiv preprint arXiv:2112.05672, 2021.
  * [11] Marcus Michelen.  Real roots near the unit circle of random polynomials.  Transaction of AMS, 2020.
  * [12] Marcus Michelen and Julian Sahasrabudhe.  Random polynomials: the closest roots to the unit circle.  arXiv preprint arXiv:2010.10869, 2020.
  * [13] Oanh Nguyen and Van Vu.  Random polynomials: central limit theorems for the real roots.  Duke Mathematical Journal, 170(17):3745–3813, 2021.
  * [14] Oanh Nguyen and Van Vu.  Roots of random functions: A framework for local universality.  American Journal of Mathematics, 144(1):1–747, 2022.
  * [15] A. M. Odlyzko and B. Poonen.  Zeros of polynomials with $0,1$ coefficients.  Enseign. Math. (2), 39(3-4):317–348, 1993.
  * [16] Larry A Shepp and Robert J Vanderbei.  The complex zeros of random polynomials.  Transactions of the American Mathematical Society, pages 4365–4384, 1995.


# Operational approach to bilocality with joint probability distributions

Kelvin Onggadinata  Centre for Quantum Technologies, National University of
Singapore, 3 Science Drive 2, 117543 Singapore, Singapore  Department of
Physics, National University of Singapore, 3 Science Drive 2, 117543
Singapore, Singapore  Paweł Kurzyński  Institute of Spintronics and Quantum
Information, Faculty of Physics, Adam Mickiewicz University, 61-614 Poznań,
Poland  Centre for Quantum Technologies, National University of Singapore, 3
Science Drive 2, 117543 Singapore, Singapore  Dagomir Kaszlikowski  Centre for
Quantum Technologies, National University of Singapore, 3 Science Drive 2,
117543 Singapore, Singapore  Department of Physics, National University of
Singapore, 3 Science Drive 2, 117543 Singapore, Singapore

###### Abstract

We show an operational approach to bilocality with quasi-probability
distributions and quasi-stochastic processes. This approach clearly
demonstrates that negative probabilities are necessary to violate bilocality.
It also highlights a subtle interplay between bilocal and local correlations
and it can be easily extended to study N-locality.

##  I Introduction

Bilocality, discovered in [1], is a special case of local hidden variables
(LHV) model where three observers, Alice, Bob and Charlie share two sources
distributing LHVs (Alice-Bob and Bob-Charlie), instead of a common source,
feeding three of them at the same time. This scenario is relevant for the
fundamental research into non-local correlations, both quantum and post-
quantum (for instance, PR-boxes [2]), but it also has a more pragmatic side
related to quantum networks [3, 4, 5] where $N$ observers share non-local
correlations ($N$-locality), distributed in a variety of ways. Although,
$N$-locality is a subset of LHV hypothesis that has been around since 1935 [6]
it is a highly non-trivial one because of its non-linearity.

Let us now describe bilocality using mathematical language, most commonly used
in the literature on the topic (see [7]). In the subsequent sections we will
change this language, gaining new techniques and new qualitative insights into
the problem at hand.

Consider three parties, Alice, Bob, and Charlie performing measurements with
random inputs $x$, $y$, and $z$ and the corresponding outcomes $a$, $b$, and
$c$. A tripartite collection of probability distributions is called local if
the probability of observing the outcomes $a,b,c$ given the inputs $x,y,z$ can
be written as

$P(a,b,c|x,y,z)=\sum_{\lambda}\rho(\lambda)P(a|x,\lambda)P(b|y,\lambda)P(c|z,\lambda)\,.$
(1)

$\lambda$’s are called local hidden variables (LHVs) and they are
deterministic instructions of what output to produce for a given measurement
setting. For instance, if the inputs are binary ($0,1$) and outcomes are
dichotomic ($\pm 1$), a hidden variable can say: “output $+1$ for Alice if she
chooses the input $0$, $-1$ if she chooses $1$, output $-1$ for Bob if he
chooses $0$ and $-1$ if he chooses $1$ and $-1$ for Charlie choosing $0$ and
$+1$ if it is $1$.” This can be conveniently denoted as a string
$(+1,-1;-1,-1;-1,+1)$. These instructions are distributed by the source with
some probability $\rho$ and, when they reach measuring apparatus, they are
further locally processed with conditional probabilities
$P(a|x,\lambda),P(b|y,\lambda)$ and $P(c|z,\lambda)$.

We always insist that $P(a,b,c|x,y,z)$ is no-signalling, i.e., the
distributions of local outcomes depend only on local inputs. For instance,
$P(a|x)=\sum_{b,c}P(a,b,c|x,y,z)$ depends only on $x$, not on $y$ or/and $z$.

This scenario assumes that the instructions encoded in $\lambda$’s are
distributed to all parties simultaneously. This can be relaxed. For instance,
in Fig. 1, there are two independent sources $S_{1}$ and $S_{2}$, a scenario
we call bilocal [8]. This assumption further decomposes Eq. (1) into

$\displaystyle P(a,b,c|x,y,z)$ $\displaystyle=$
$\displaystyle\sum_{\lambda_{1}\lambda_{2}}\rho_{1}(\lambda_{1})\rho_{2}(\lambda_{2})P(a|x,\lambda_{1})$
(2) $\displaystyle\quad\times
P(b|y,\lambda_{1},\lambda_{2})P(c|z,\lambda_{2})\,.$

This relaxation has lead to new, interesting scenarios, like the one discussed
in [1] and has been subsequently generalised to a more complex setups [9].

Figure 1: Illustration of a typical bilocal hidden variable model. Two
independent sources, $S_{1}$ and $S_{2}$, emit hidden variables $\lambda_{1}$
and $\lambda_{2}$ to the pair Alice-Bob(B) and Bob(B’)-Charlie. Alice, Bob,
and Charlie perform random measurements $x,y,z$ with the corresponding
outcomes $a,b,c$. Note that Bob measures jointly both subsystems he receives.

To proceed, we first introduce a notation, which we will use interchangeably
with the previous one:

$P(a,b,c|x,y,z)\leftrightarrow P(a_{x},b_{y},c_{z})\,.$

Now, a joint probability distribution is a probability distribution $q$ over
all measurement outcomes $\vec{a}=(a_{0},...,a_{|X|-1})$,
$\vec{b}=(b_{0},...,b_{|Y|-1})$, $\vec{c}=(c_{0},...,c_{|Z|-1})$ with the
following properties:

$\displaystyle\sum_{\vec{a},\vec{b},\vec{c}}q(\vec{a},\vec{b},\vec{c})$
$\displaystyle=$ $\displaystyle 1$ (3) $\displaystyle\sum_{\vec{a}\setminus
a_{x},\vec{b}\setminus b_{y},\vec{c}\setminus
c_{z}}q(\vec{a},\vec{b},\vec{c})$ $\displaystyle=$ $\displaystyle
P(a_{x},b_{y},c_{z}).$ (4)

Here, $|X|,|Y|$ and $|Z|$ denote the number of measurement settings for Alice,
Bob and Charlie, and $\vec{a}\setminus a_{x}$ indicates all the variables in
$\vec{a}$, excluding $a_{x}$. In this paper we consider positive joint
probability distributions (JPD), i.e., all probabilities are positive, and
negative joint quasi-probability distributions (JQD) where some probabilities
can be negative. Note that both JPD and JQD are automatically no-signalling.

Fine [10] showed that one can find a JPD with properties (3) and (4) if and
only if the marginal distribution $P$ in (4) is local. If $P$ is nonlocal,
some of $q$’s entries must be negative. Of course, (4) guarantees that all
observable statistics are semi-positive or else we are forced to find an
operational interpretation of negative probabilities, which goes beyond the
scope of this paper.

When all outputs are binary, any joint probability distribution (JPD and JQD)
can be written in a compact form by labelling the outcomes as $\\{+1,-1\\}$.
In the simplest Bell-CHSH scenario with two parties, two inputs each, we can
change $(a_{0}a_{1}b_{0}b_{1})\to(s_{0}s_{1}s_{2}s_{3})$ for convenience,
giving us

$\displaystyle q(s_{0}s_{1}s_{2}s_{3})$ $\displaystyle=$
$\displaystyle\frac{1}{2^{4}}\Big{[}1+\sum_{j=0}^{3}s_{j}E^{(1)}_{j}+\sum_{j<k}s_{j}s_{k}E^{(2)}_{jk}$
(5)
$\displaystyle\,+\sum_{j<k<l}s_{j}s_{k}s_{l}E^{(3)}_{jkl}+s_{0}s_{1}s_{2}s_{3}E^{(4)}\Big{]}$

where $E^{(1)}_{j}=\langle s_{j}\rangle_{q}$, $E^{(2)}_{jk}=\langle
s_{j}s_{k}\rangle_{q}$ etc. The observable statistics drawn from $P$ determine
the correlation functions (correlators) $E^{(1)}_{j}$’s, $E^{(2)}_{02}$,
$E^{(2)}_{03}$, $E^{(2)}_{12}$ and $E^{(2)}_{13}$. The remaining, unobservable
correlators, should satisfy $|E|\leq 1$ to be consistent with observable
outcomes $\pm 1$ but are arbitrary otherwise. The extension to more inputs and
more parties is straightforward.

##  II The bilocal scenario with joint probability distributions

We illustrate the bilocal scenario with JPDs in the Fig. 2. It can be
mathematically written as

$r(\vec{a},\vec{b},\vec{c})=\sum_{\vec{\beta},\vec{\beta}^{\prime}}S(\vec{b}|\vec{\beta},\vec{\beta}^{\prime})r_{0}(\vec{a},\vec{\beta},\vec{\beta}^{\prime},\vec{c})\,,$
(6)

where $\vec{a}=(a_{0},a_{1},\dots,a_{N-1})$ is the vector of measurement
outcomes $a_{k}=\pm 1$ ($k=0,1,\dots,N-1$) for $N$ settings (inputs). Similar
definition applies to $\vec{b},\vec{c},\vec{\beta},\vec{\beta}^{\prime}$ but
these vectors do not necessarily have the same number of measurement settings.
$S(\vec{b}|\vec{\beta},\vec{\beta}^{\prime})$ represents Bob’s data
processing. Bilocality imposes the following constraints:

  1. (i)

Independent preparation — the initial JPD takes the factorization form

$r_{0}(\vec{a},\vec{\beta},\vec{\beta}^{\prime},\vec{c})=q_{AB}(\vec{a},\vec{\beta})q_{B^{\prime}C}(\vec{\beta}^{\prime},\vec{c})\,,$
(7)

where
$q_{AB}(\vec{a},\vec{\beta}),q_{B^{\prime}C}(\vec{\beta}^{\prime},\vec{c})\geq
0$ are generated by two completely independent local sources.

  2. (ii)

Local processing — Mirroring the assumptions made in (2), Bob has a device
that stochastically processes data received from both sources. This device is
mathematically a stochastic process
$S(\vec{b}|\vec{\beta},\vec{\beta}^{\prime})$ such that
$S(\vec{b}|\vec{\beta},\vec{\beta}^{\prime})\geq 0$ and
$\sum_{\vec{b}}S(\vec{b}|\vec{\beta},\vec{\beta}^{\prime})=1$.

Figure 2: Alternative description of the bilocal model in terms of joint
probability distributions with each boxes having binary inputs and outputs.
Initially, two independent sources distribute local hidden variables with
joint probability distributions $q_{AB}(a_{0}a_{1}\beta_{0}\beta_{1})$ and
$q_{B^{\prime}C}(\beta_{0}^{\prime}\beta_{1}^{\prime}c_{0}c_{1})$ to the pair
Alice-Bob(B) and Bob(B’)-Charlie, respectively. The joint measurement on Bob’s
part is replaced by a local processing $S$ acting on the previously obtained
outcomes
$\\{\beta_{0},\beta_{1},\beta_{0}^{\prime},\beta_{1}^{\prime}\\}\rightarrow\\{b_{0},b_{1}\\}$.

###  II.1 Bilocal inequality

To warm up and see the joint probability distribution formalism at work, let
us re-derive the bilocal inequality from [8]. In the scenario considered
there, Alice and Charlie have two inputs and two dichotomic outputs, i.e.,
$\vec{a}=(a_{0},a_{1})$ and $\vec{c}=(c_{0},c_{1})$ with
$a_{x},c_{z}\in\\{+1,-1\\}$.

Let us define

$\displaystyle\left\langle I\right\rangle_{r}$ $\displaystyle=$
$\displaystyle\frac{1}{4}\sum_{i,j=0,1}\left\langle
A_{i}B_{0}C_{j}\right\rangle_{r}\,,$ (8) $\displaystyle\left\langle
J\right\rangle_{r}$ $\displaystyle=$
$\displaystyle\frac{1}{4}\sum_{i,j=0,1}(-1)^{i+j}\left\langle
A_{i}B_{1}C_{j}\right\rangle_{r}\,,$ (9)

where

$\left\langle
A_{i}B_{j}C_{k}\right\rangle_{r}=\sum_{a_{0}a_{1}}\sum_{b_{0}b_{1}}\sum_{c_{0}c_{1}}a_{i}b_{j}c_{k}r(a_{0}a_{1}b_{0}b_{1}c_{0}c_{1})\,.$
(10)

We now show that any bilocal joint probability distribution given in (6)
together with assumptions (i)-(ii) imply

$\sqrt{|\left\langle I\right\rangle_{r}|}+\sqrt{|\left\langle
J\right\rangle_{r}|}\leq 1\,.$ (11)

First, let us take the absolute value of (8) and use the inequality
$|\sum_{k}x_{k}|\leq\sum_{k}|x_{k}|$ to obtain

$\displaystyle|\left\langle I\right\rangle_{r}|$ $\displaystyle\leq$
$\displaystyle\sum_{\vec{a},\vec{b},\vec{c}}\sum_{\vec{\beta},\vec{\beta}^{\prime}}\frac{1}{4}|b_{0}(a_{0}+a_{1})(c_{0}+c_{1})|$
(12) $\displaystyle\quad\times
S(\vec{b}|\vec{\beta},\vec{\beta}^{\prime})q_{AB}(\vec{a},\vec{\beta})q_{B^{\prime}C}(\vec{\beta}^{\prime},\vec{c})$
$\displaystyle=$
$\displaystyle\sum_{\vec{a},\vec{\beta}}\frac{|a_{0}+a_{1}|}{2}q_{AB}(\vec{a},\vec{\beta})$
$\displaystyle\times\sum_{\vec{\beta}^{\prime},\vec{c}}\frac{|c_{0}+c_{1}|}{2}q_{B^{\prime}C}(\vec{\beta}^{\prime},\vec{c})$
$\displaystyle\times\sum_{\vec{b}}|b_{0}|S(\vec{b}|\vec{\beta},\vec{\beta}^{\prime})\,.$

Since $b_{0}$ can only take values $\pm 1$, we have
$\sum_{\vec{b}}|b_{0}|S(\vec{b}|\vec{\beta},\vec{\beta}^{\prime})=\sum_{\vec{b}}S(\vec{b}|\vec{\beta},\vec{\beta}^{\prime})=1$.
The above simplifies to

$|\left\langle
I\right\rangle_{r}|\leq\sum_{a_{0}a_{1}}\frac{|a_{0}+a_{1}|}{2}p_{A}(\vec{a})\sum_{c_{0}c_{1}}\frac{|c_{0}+c_{1}|}{2}p_{C}(\vec{c})$
(13)

with $p_{A}(\vec{a})=\sum_{\vec{\beta}}q_{AB}(\vec{a},\vec{\beta})$,
$p_{C}(\vec{c})=\sum_{\vec{\beta}^{\prime}}q_{B^{\prime}C}(\vec{\beta}^{\prime},\vec{c})$.

In a similar manner we get

$|\left\langle
J\right\rangle_{r}|\leq\sum_{a_{0}a_{1}}\frac{|a_{0}-a_{1}|}{2}p_{A}(\vec{a})\sum_{c_{0}c_{1}}\frac{|c_{0}-c_{1}|}{2}p_{C}(\vec{c})\,.$
(14)

We use the inequality
$\sqrt{uv}+\sqrt{u^{\prime}v^{\prime}}\leq\sqrt{u+u^{\prime}}\sqrt{v+v^{\prime}}$
for $u,v,u^{\prime},v^{\prime}\geq 0$ and obtain

$\displaystyle\sqrt{|\left\langle I\right\rangle_{r}|}+\sqrt{|\left\langle
J\right\rangle_{r}|}$ (15)
$\displaystyle\quad\leq\sqrt{\sum_{a_{0}a_{1}}\left(\frac{|a_{0}+a_{1}|}{2}+\frac{|a_{0}-a_{1}|}{2}\right)p_{A}(a_{0}a_{1})}$
$\displaystyle\quad\times\sqrt{\sum_{c_{0}c_{1}}\left(\frac{|c_{0}+c_{1}|}{2}+\frac{|c_{0}-c_{1}|}{2}\right)p_{C}(c_{0}c_{1})}=1\,.$

The last equality comes from
$|a_{0}+a_{1}|/2+|a_{0}-a_{1}|/2=\max\\{|a_{0}|,|a_{1}|\\}=1$ and
$|c_{0}+c_{1}|/2+|c_{0}-c_{1}|/2=\max\\{|c_{0}|,|c_{1}|\\}=1$. The remaining
summation over $p_{A}(a_{0}a_{1})$ and $p_{C}(c_{0}c_{1})$ simply returns 1.

##  III Violations of bilocality

In this section we use the JPD/JQD formalism to explain how violation of
bilocality reported in [8] happens in simple, operational terms. To be clear,
we focus on this particular scenario but whatever we present here can be
applied to any $N$-locality scenario with appropriate modifications.

In [8] the authors model Bob’s data processing on entanglement swapping. They
also assume that Alice and Bob(B), as well as, Bob(B’) and Charlie, share
maximally entangled states. Their measurement settings are carefully chosen to
avoid revealing any quantum non-locality. These assumptions can be now
directly mapped into our formalism: (1) we need to choose Alice-Bob(B) and
Bob(B’)-Charlie joint probability distributions and (2) Bob(BB’) data
processing $S$.

What follows is an extension of the joint probability distribution approach
considered in [11], now applied to bilocality instead of just locality.
Locality is equivalent to JPD and nonlocality requires JQD. On the other hand,
bilocal correlations are included in local correlations, therefore there are
some local correlations that are not bilocal and should violate the bilocality
inequality. Such correlations must be described by a JPD, hence non-bilocality
must be manifested by some more subtle properties of JPDs. Below we find what
these properties are.

To this end, consider the following joint probability distribution

$\displaystyle q_{AB}(a_{0}a_{1}\beta_{0}\beta_{1})$
$\displaystyle=\frac{1}{16}\Big{[}1+\mu_{1}(a_{0}\beta_{0}+a_{0}\beta_{1}+a_{1}\beta_{0}-a_{1}\beta_{1})\Big{]}$
(16a) $\displaystyle
q_{B^{\prime}C}(\beta_{0}^{\prime}\beta_{1}^{\prime}c_{0}c_{1})$
$\displaystyle=\frac{1}{16}\Big{[}1+\mu_{2}(\beta_{0}^{\prime}c_{0}+\beta_{0}^{\prime}c_{1}+\beta_{1}^{\prime}c_{0}-\beta_{1}^{\prime}c_{1})\Big{]}$
(16b)

where $\mu_{1,2}\in[0,1]$. The above $q$’s have the following physical
interpretations: (i) for $0\leq\mu_{1,2}\leq\frac{1}{2}$ we have JPDs that
saturate the LHV bound for the Bell-CHSH inequality; (ii) for
$\frac{1}{2}<\mu_{1,2}\leq\frac{1}{\sqrt{2}}$ we have JQDs corresponding to
the quantum mechanical violations of the Bell-CHSH inequality; (iii) for
$\frac{1}{2}<\mu_{1,2}\leq 1$ we deal with JQDs representing super quantum
correlations generated by noisy PR-boxes; (iv) the region for $\mu_{1,2}>1$
has no known physical meaning because one would observe negative probabilities
in observable probability marginals.

Next we chose Bob’s stochastic process as

$\displaystyle
S_{\eta}(b_{0}b_{1}|\beta_{0}\beta_{1}\beta_{0}^{\prime}\beta_{1}^{\prime})$
$\displaystyle\qquad=S_{\eta}(b_{0}|\beta_{0}\beta_{0}^{\prime})S_{\eta}(b_{1}|\beta_{1}\beta_{1}^{\prime})$
$\displaystyle\qquad=\frac{1}{2}\Big{[}1+\eta
b_{0}\beta_{0}\beta_{0}^{\prime}\Big{]}\times\frac{1}{2}\Big{[}1+\eta
b_{1}\beta_{1}\beta_{1}^{\prime}\Big{]}\,,$ (17)

where $S_{\eta}$ is positive for $\eta\leq 1$. Inserting these expressions
into (6), we get

$\displaystyle r(a_{0}a_{1}b_{0}b_{1}c_{0}c_{1})$ $\displaystyle=$
$\displaystyle\frac{1}{2^{6}}\Big{[}1+\eta\mu_{1}\mu_{2}\big{\\{}b_{0}(a_{0}+a_{1})(c_{0}+c_{1})$
(18) $\displaystyle\qquad+b_{1}(a_{0}-a_{1})(c_{0}-c_{1})\big{\\}}\Big{]}\,,$

whence

$\sqrt{|\left\langle I\right\rangle_{r}|}+\sqrt{|\left\langle
J\right\rangle_{r}|}=2\sqrt{\eta\mu_{1}\mu_{2}}\,.$ (19)

Thus violations of the bilocal inequality happen when
$\eta\mu_{1}\mu_{2}>\frac{1}{4}$.

Let us summarize the bounds on $\eta\mu_{1}\mu_{2}$ for the different regions

$\eta\mu_{1}\mu_{2}\stackrel{{\scriptstyle\text{BL}}}{{\leq}}\frac{1}{4}\stackrel{{\scriptstyle\text{BQ}}}{{\leq}}\frac{1}{2}\stackrel{{\scriptstyle\text{L,BNS}}}{{\leq}}1$
(20)

where BQ and BNS refers to bi-quantum and bi-no-signaling, respectively. Note
that the upper bound is imposed such that negative probabilities are not
observed at the level of observable marginals
$r_{\eta\mu_{1}\mu_{2}}(a_{j}b_{k}c_{l})$.

The violation of the bilocal inequality depends on the values of $\eta$ and
$\mu_{1,2}$. In Fig. 3 we show different violation regions that correspond to
different origins of negativity. Here we assume that $\mu_{1}=\mu_{2}=\mu$, so
$q_{AB}$ and $q_{B^{\prime}C}$ are symmetric in their correlations. We group
the violations as follows:

  * •

JQD and positive $S_{\eta}$ (orange shaded region): $\mu>\frac{1}{2}$ and
$\eta\leq 1$.

  * •

JPD and negative $S_{\eta}$ (yellow shaded region): $\mu\leq\frac{1}{2}$ and
$\eta>1$.

  * •

JQD and negative $S_{\eta}$ (gray shaded region): $\mu>\frac{1}{2}$ and
$\eta>1$.

Therefore, negative probabilities are required either in the initial joint
probability distributions shared between Alice-Bob(B) and Bob(B’)-Charlie, or
in the local processing performed by Bob(BB’). This is also true for
$\mu_{1}\neq\mu_{2}$.

Figure 3: Strength of the nonbilocal correlations as a function of $\eta$ and
$\mu$. The regions below the blue, red and black lines corresponds to bilocal,
bi-quantum, and bi-no-signaling correlations, respectively. The three shaded
regions correspond to different origins of negativity. The details are
described in the main text.

The above classification of different mechanisms responsible for violations of
bilocality is one of the main results of this paper. We can now see that
violations of bilocality require negative probabilities either in the initial
probability distributions shared by the observers or in Bob’s (BB’) data
processing or both. We discussed a similar mechanism in the context of quantum
nonlocality in [12]. This brings locality and bilocality violations into a
universal mathematical and conceptual framework.

Let us make some remarks and observations. Note that the JQD

$\displaystyle r_{\text{BQ}}(a_{0}a_{1}b_{0}b_{1}c_{0}c_{1})$ $\displaystyle=$
$\displaystyle\frac{1}{2^{6}}\Big{[}1+\frac{1}{2}\Big{\\{}b_{0}(a_{0}+a_{1})(c_{0}+c_{1})$
(21) $\displaystyle\quad+b_{1}(a_{0}-a_{1})(c_{0}-c_{1})\Big{\\}}\Big{]}$

saturates the maximal BQ bound. If we mix $r_{\text{BQ}}$ with white noise, we
get a noisy JQD

$\displaystyle r_{\text{BQ}}^{V}(a_{0}a_{1}b_{0}b_{1}c_{0}c_{1})$
$\displaystyle=$ $\displaystyle
Vr_{\text{BQ}}(a_{0}a_{1}b_{0}b_{1}c_{0}c_{1})+\frac{1-V}{2^{6}}$ (22)
$\displaystyle=$
$\displaystyle\frac{1}{2^{6}}\Big{[}1+\frac{V}{2}\Big{\\{}b_{0}(a_{0}+a_{1})(c_{0}+c_{1})$
$\displaystyle\quad+b_{1}(a_{0}-a_{1})(c_{0}-c_{1})\Big{\\}}\Big{]}\,,$

where $V\in[0,1]$ is customarily called visibility. The above JQD becomes
bilocal when $V\leq\frac{1}{2}$, in which case it also becomes a JPD. This is
in agreement with [8]. However, note that the above is a global noise
admixture. If each source that produces the initial correlations introduces
their own white noise with probability $1-v_{p}$ and $1-v_{q}$, then
$V=v_{p}v_{q}$.

In addition, $r_{\text{BQ}}$ is not unique. One could choose a different JQD
that mimics $r_{\text{BQ}}$, i.e., having the same observable probabilities.
In fact, there exists a JPD that is equivalent to $r_{\text{BQ}}$ and still
violates the bilocal inequality. For instance, one could consider

$\displaystyle r_{\text{BQ}}^{+}(a_{0}a_{1}b_{0}b_{1}c_{0}c_{1})$
$\displaystyle=\frac{1}{2^{6}}\Big{[}1+\frac{1}{2}\Big{\\{}b_{0}(a_{0}+a_{1})(c_{0}+c_{1})$
$\displaystyle\quad+b_{1}(a_{0}-a_{1})(c_{0}-c_{1})+a_{0}a_{1}c_{0}c_{1}\Big{\\}}\Big{]}\,.$
(23)

It is easy to check that it is nonnegative. However, it is also inseparable

$\sum_{b_{0}b_{1}}r_{\text{BQ}}^{+}(a_{0}a_{1}b_{0}b_{1}c_{0}c_{1})\neq
p(a_{0}a_{1})q(c_{0}c_{1})\,.$ (24)

Hence, it does not come from two independent sources, it is not interesting to
us. Moreover, the JQD in Eq. (21) is optimal in the sense that it has the
minimal amount of negativity while staying separable. Here, we use the
$L_{1}$-norm measure of negativity [13] defined as

$N(q)=\sum_{x}|q(x)|\geq 1$ (25)

with equality saturated iff $q$ is completely positive.

Interestingly, since $b_{0}$ and $b_{1}$ are co-measurable in the entanglement
swapping scenario utilised in [8], we expect that the marginals
$r(a_{x}b_{0}b_{1}c_{z})$ are semi-positive. Indeed, they behave properly up
to the maximal BQ bound and break down once we exceed it. This is not
surprising because above this bound we deal with noisy PR-boxes that do not
have a quantum mechanical counterpart. We designed Bobs (BB’)
(quasi)stochastic process to simulate quantum mechanical measurement, not an
imaginary PR-box measurement. Perhaps this is an argument one can use to
reject PR-boxes as non physical but more research is needed, which goes beyond
the scope of this paper.

Finally, let us discuss Bob’s transformation $S_{\eta}$. When $\eta=1$,
$S_{\eta=1}$ is the classical AND operation on the bits
$(\beta_{0},\beta_{0}^{\prime})\to b_{0}$ and
$(\beta_{1},\beta_{1}^{\prime})\to b_{1}$. $S_{\eta=1}$ is a positive
stochastic matrix, consistent with the requirement that for any JQD $p$ and
$q$ the final JQD $r$ has positive observable marginals. Note that any
stochastic matrix $S_{\eta}$, acting on JQD, does not increase its negativity.
When $\eta>1$, one can think of $S_{\eta}$ as a nebit (negative bit)
controlled operation where AND is performed with probability $(1+\eta)/2$ and
NAND with probability $(1-\eta)/2$. Nebit, introduced in [12], is a catalyst
turning stochastic operations into quasi-stochastic ones and it is useful to
quantify non-classical behaviour or in this case, non-bilocality.

##  IV Final remarks

Bilocality is a relatively new phenomenon in the topic of quantum and post-
quantum nonlocality that has so far generated a host of interesting and still
unsolved problems, see for instance [4, 9]. It has also practical implications
with the advent of quantum internet [14, 15] and distributed quantum computing
[16, 17].

In this paper, we present a certain way to analyse bilocality and more
generally, $N$-locality, based on a well established, albeit still niche,
research into quasi-probabilities and quasi-stochastic processes. We
demonstrated that our approach recovers known results (the bilocality
inequality from [8]) but it also offers operational insights (stochastic and
quasi-stochastic data processing) into the mechanisms responsible for
violations of bilocality that have not been yet discussed in the literature.

Extensions to $N$-locality are straightforward as adding more independent JPDs
and/or JQDs and more quasi-stochastic and/or stochastic processes does not
change the formal structure of our approach.

## Acknowledgements

We thank Valerio Scarani for suggestions and discussions. This research is
supported by the National Research Foundation, Singapore and A*STAR under its
CQT Bridging Grant. PK is supported by the Polish National Science Centre
(NCN) under the Maestro Grant no. DEC-2019/34/A/ST2/00081.

## References

  * Branciard  _et al._ [2010] C. Branciard, N. Gisin, and S. Pironio, Phys. Rev. Lett. 104, 170401 (2010).
  * Popescu and Rohrlich [1994] S. Popescu and D. Rohrlich, Foundations of Physics 24, 379 (1994).
  * Rosset  _et al._ [2016] D. Rosset, C. Branciard, T. J. Barnea, G. Pütz, N. Brunner, and N. Gisin, Phys. Rev. Lett. 116, 010403 (2016).
  * Renou  _et al._ [2019] M.-O. Renou, E. Bäumer, S. Boreiri, N. Brunner, N. Gisin, and S. Beigi, Phys. Rev. Lett. 123, 140401 (2019).
  * Hansenne  _et al._ [2022] K. Hansenne, Z.-P. Xu, T. Kraft, and O. Gühne, Nature communications 13, 1 (2022).
  * Einstein  _et al._ [1935] A. Einstein, B. Podolsky, and N. Rosen, Phys. Rev. 47, 777 (1935).
  * Brunner  _et al._ [2014] N. Brunner, D. Cavalcanti, S. Pironio, V. Scarani, and S. Wehner, Rev. Mod. Phys. 86, 419 (2014).
  * Branciard  _et al._ [2012] C. Branciard, D. Rosset, N. Gisin, and S. Pironio, Phys. Rev. A 85, 032119 (2012).
  * Tavakoli  _et al._ [2022] A. Tavakoli, A. Pozas-Kerstjens, M.-X. Luo, and M.-O. Renou, Reports on Progress in Physics 85, 056001 (2022).
  * Fine [1982] A. Fine, Phys. Rev. Lett. 48, 291 (1982).
  * Abramsky and Brandenburger [2014] S. Abramsky and A. Brandenburger, in  _Horizons of the mind. A tribute to Prakash Panangaden_ (Springer, 2014) pp. 59–75.
  * Onggadinata  _et al._ [2021] K. Onggadinata, P. Kurzynski, and D. Kaszlikowski, arXiv preprint arXiv:2106.07945 https://doi.org/10.48550/arXiv.2106.07945 (2021).
  * Oas  _et al._ [2014] G. Oas, J. A. de Barros, and C. Carvalhaes, Physica Scripta T163, 014034 (2014).
  * Kimble [2008] H. J. Kimble, Nature 453, 1023 (2008).
  * Wehner  _et al._ [2018] S. Wehner, D. Elkouss, and R. Hanson, Science 362, eaam9288 (2018).
  * Acín  _et al._ [2006] A. Acín, N. Gisin, and L. Masanes, Phys. Rev. Lett. 97, 120405 (2006).
  * Acín  _et al._ [2007] A. Acín, N. Brunner, N. Gisin, S. Massar, S. Pironio, and V. Scarani, Phys. Rev. Lett. 98, 230501 (2007).