# Deep Reinforcement Learning Control for Disturbance Rejection in a Nonlinear
Dynamic System with Parametric Uncertainty

Vincent W. Hill Tuscaloosa, AL, USA  The University of Alabama

###### Abstract

This work describes a technique for active rejection of multiple independent
and time-correlated stochastic disturbances for a nonlinear flexible inverted
pendulum with cart system with uncertain model parameters. The control law is
determined through deep reinforcement learning, specifically with a continuous
actor-critic variant of deep Q-learning known as Deep Deterministic Policy
Gradient, while the disturbance magnitudes evolve via independent stochastic
processes. Simulation results are then compared with those from a classical
control system.

###### Index Terms:

deep reinforcement learning, disturbance rejection, nonlinear control,
parametric uncertainty

##  I Introduction

One of the most challenging tasks in the design of a new system is dynamic
modeling for simulation. This is especially true for complex vehicles with
highly nonlinear dynamics. Specifically, flexible dynamics have become more
critical in recent years as a result of the drive to improve vehicle
performance by decreasing structural weight. These systems are substantially
more difficult to precisely model without undertaking costly and time-
consuming ground vibration tests, so numerical studies must be used to capture
the effects of flexibility on vehicle dynamics. While these techniques have
improved in the last two decades, they are limited in scope and are often only
used temporarily until flight data can be gathered.

Added to the challenge of flexible nonlinear dynamic system modeling is
prohibitively expensive prototype testing, especially in the case of large
aircraft and launch vehicles. Despite these difficulties, high-performance
control systems must still be developed in time for the vehicle’s initial
deployment. This means control system design must be accomplished with models
carrying a high degree of uncertainty, particularly in structural flexibility.

In addition to potential for uncertainty in vehicle mathematical modeling,
aircraft and spacecraft are subject to random excitations arising from
atmospheric conditions. These disturbances can be mitigated through passive
means, such as fins for launch vehicles or static stability for an aircraft,
but these passive disturbance alleviation methods are limited in
effectiveness. Active disturbance rejection techniques can handle a broad
class of excitation without impacting vehicle performance.

These factors place a premium on control systems that are robust to both a
high degree of dynamic model uncertainty as well as stochastic disturbances
while maintaining high performance. Deep Reinforcement Learning (DRL) is a new
technology that has emerged in the last several years is well-suited to tackle
this difficult problem. Mnih et al. [1] introduced the Deep Q Network (DQN)
algorithm in 2015 as a method for achieving human-level performance playing
Atari games such as Breakout by coupling the Q-learning algorithm with deep
neural networks for action-value function approximation. DQN operates by
recursively updating the Q-function using dynamic programming, and actions are
chosen by taking the $argmax$ of $Q$ conditioned on the current state. Over
time the DQN algorithm, known as an agent, will optimize this function to
determine the best action for every possible state.

DQN is a very effective solution for low-dimensional discrete action spaces
such as an Atari controller, but struggles in a continuous action space such
as the gimbal angle of a launch vehicle thrust vector controller or the
deflection angle of an aircraft elevator. Lillicrap et al. [2] proposed an
extension to DQN called Deep Deterministic Policy Gradient (DDPG) that
utilizes an actor-critic algorithm with deep function approximators to learn
policies in high-dimensional continuous action spaces such as those in a
control design problem.

To demonstrate the effectiveness of deep reinforcement learning control for
active disturbance rejection in nonlinear systems with parametric uncertainty,
the technique is applied to a flexible inverted pendulum with cart (FIPWC)
mathematical model with parametric uncertainty and stochastic disturbances. In
this work, DRL is used to learn the nonlinear control signal that will keep
the pendulum vertical throughout the simulation while under the influence of
time-correlated stochastic disturbances. Results from that control system are
then compared to results from the use of a conventional proportional-
derivative (PD) control law.

##  II Technical Approach

###  II-A Flexible Inverted Pendulum with Cart System Modeling

The FIPWC system is modeled as a beam with a tip mass, with flexibility
effects modeled as simple linear springs attached at the pendulum tip and
center of mass. Figures 1 and 2 provide a visual overview of the FIPWC system.

Figure 1: FIPWC system diagram.

Figure 2: FIPWC system diagram.

The FIPWC mathematical model used in this work was developed by Gorade et al.
[6] through a multi-body approach with Lagrangian equation of motion
derivation. Equation 1 provides the Lagrangian of the system, where $m_{t}$,
$m_{b}$, and $m_{c}$ are the masses of the pendulum tip, pendulum beam, and
cart, respectively, $L$ is the length of the pendulum, $l=\frac{1}{2}L$ is the
location of the beam center of mass, $k_{1}$ and $k_{2}$ are the spring
constants of the pendulum tip and pendulum beam, respectively, $z$ is the cart
motion, $\phi$ is the angular motion of the pendulum center of mass, $\theta$
is the angular motion of the pendulum tip, and $F$ is the force applied to the
cart. Table II in the Appendix provides values for the model parameters.

$\mathcal{L}=\frac{1}{2}m_{c}\dot{z}^{2}+\frac{1}{2}m_{b}(\dot{z}^{2}_{2}l\dot{z}\dot{\phi}cos\phi+l^{2}\dot{\phi})$
$+\frac{1}{2}m_{t}[\dot{z}^{2}+2\dot{z}(L\dot{\theta}cos\theta+l\dot{\phi}cos\phi+L^{2}\dot{\theta}^{2}$
$+l^{2}\dot{\phi}^{2}+2Ll\dot{\theta}\dot{\phi}cos(\theta-\phi)]+m_{t}gL(1-cos\theta)$
$-\frac{1}{2}k_{1}(Lsin\theta)^{2}+m_{b}gl(1-cos\phi)-\frac{1}{2}k_{2}(lsin\phi)^{2}$
(1)

Structural damping of the pendulum as well as friction damping of the cart are
modeled as Rayleigh damping, given in Equation 2 where $b_{1}$, $b_{2}$, and
$b_{3}$ are the damping coefficients of the tip, center of mass, and cart,
respectively.

$R_{d}=\frac{1}{2}b_{3}\dot{z}^{2}+\frac{1}{2}b_{2}[(\dot{z}+l\dot{\phi}cos\phi)^{2}+l^{2}\dot{\phi}^{2}sin^{2}\phi]+$

$\frac{1}{2}b_{1}[(\dot{z}+l\dot{\phi}cos\phi+L\dot{\theta}cos\theta)^{2}+(l\dot{\phi}sin\phi+L\dot{\theta}sin\theta)^{2}]$
(2)

Equation 3 describes the Euler-Lagrange equation that when carried out with
the Lagrangian given in Equations 1 and 2 produces the system equations of
motion, where $w=\begin{bmatrix}\theta&\phi&z\\\ \end{bmatrix}$ and
$F=\begin{bmatrix}0&0&F\\\ \end{bmatrix}$.

$\frac{\partial}{\partial
t}\frac{\partial\mathcal{L}}{\partial\dot{w}_{i}}-\frac{\partial\mathcal{L}}{\partial
w_{i}}+\frac{\partial R_{d}}{\partial\dot{w_{i}}}=F_{i}$ (3)

The system equations of motion are defined by Equations 4-6.

$m_{t}L^{2}\ddot{\theta}+m_{t}L\ddot{z}cos\theta+m_{t}Ll\ddot{\phi}cos(\theta-\phi)$
$+m_{t}Ll\dot{\phi}^{2}sin(\theta-\phi)+\frac{1}{2}k_{1}L^{2}sin(2\theta)-m_{t}gLsin\theta$
$b_{1}L\dot{z}cos\theta+b_{1}L^{2}\dot{\theta}+b_{1}Ll\dot{\phi}cos(\theta-\phi)=0$
(4)

$(m_{t}+m_{b})l^{2}\ddot{\phi}+(m_{t}+m_{b})l\ddot{z}cos\phi+m_{t}Ll\ddot{\theta}cos(\theta-\phi)$
$+(b_{1}+b_{2})l\dot{z}cos\phi+(b_{1}+b_{2})l^{2}\dot{\phi}+b_{1}Ll\dot{\theta}cos(\theta-\phi)$
$-m_{t}Ll\dot{\theta}^{2}sin(\theta-\phi)+\frac{1}{2}k_{2}l^{2}sin(2\phi)-m_{b}glsin\phi=0$
(5)

$(m_{c}+m_{b}+m_{t})\ddot{z}+(m_{t}+m_{b})l(\ddot{\phi}cos\phi-\dot{\phi}^{2})$
$m_{t}L(\ddot{\theta}cos\theta-\dot{\theta}^{2}sin\theta)+(b_{1}+b_{2}+b_{3})\dot{z}$
$(b_{1}+b_{2})l\dot{\phi}cos\phi+b_{1}L\dot{\theta}cos\theta=F$ (6)

This system of nonlinear differential equations is transformed into a matrix
format as shown in Equations 7 and 8. The matrix equations are numerically
integrated with the fourth order Runge-Kutta algorithm provided in the Python
module GNCPy [7].

$\vec{x}=\begin{bmatrix}z&\dot{z}&\phi&\dot{\phi}&\theta&\dot{\theta}\\\
\end{bmatrix}^{T}$ (7)

$\dot{\vec{x}}=\underline{A}(\vec{x})\vec{x}+\underline{B}(\vec{x})F+\underline{L}(\vec{x})$
(8)

Disturbances are injected into the system through the rate terms in $\vec{x}$
and values are determined via three independent Ornstein-Uhlenbeck (OU)
stochastic processes [8]. Equation 9 describes the OU process stochastic
differential equation, where $y(t)$ is the disturbance magnitude, $\kappa$ is
the rate of mean reversion, $\mu$ is the long-term mean of the process, $W(t)$
is the Brownian motion stochastic process, and $\sigma$ is the average
magnitude of the Brownian motion.

$dy(t)=\kappa(\mu-y(t))dt+\sigma dW(t)$ (9)

At each time step, the three independent stochastic processes are sampled and
added to the associated state value. Figures 3-4 provide sample time histories
of the OU processes for the cart velocity and pendulum angular velocities,
respectively. OU process parameters for each disturbance are provided in the
appendix.

Figure 3: Sample time history of cart disturbance stochastic process.

Figure 4: Sample time history of pendulum disturbance stochastic process.

###  II-B Parametric Uncertainty Modeling

Parametric uncertainty is applied to the model through the pendulum’s
structural stiffness and damping coefficients and the cart’s friction damping
coefficient, $k_{1,2}$ and $b_{1,2,3}$. At each instantiation of the Python
class describing the model, these values are sampled from a Gaussian
probability distribution with mean equal to the nominal value and covariance
of 50% of the nominal value. This is applied to the five uncertain parameters
prior to the start of each simulation run.

To protect against negative stiffness or damping due to a sample far along the
distribution tail, the absolute value of each parameter is used for
mathematical model construction. This makes the probability density function
of the uncertain parameters a folded normal distribution.

###  II-C Deep Reinforcement Learning Control

Deep reinforcement learning algorithms operate by controlling an agent that,
at each discrete time step $t$, performs actions $a$ that alter the state $s$
of an environment and receives rewards $r$. In this work, the agent is the
control system, the environment is the FIPWC model, the actions are the
control force applied to the cart, the states are the state vector given in
Equation 7, and the reward is defined below.

The behavior of the agent is defined by a policy $\pi$ which maps states to a
probability distribution over the actions. The environment $E$ is stochastic
and modeled at a Markov decision process with state space $\mathcal{S}$,
action space $\mathcal{A}=\mathbb{R}^{N}$, an initial state distribution
$p(s_{1})$, transition dynamics $p(s_{t+1}|s_{t},a_{t}$, and reward function
$r(s_{t},a_{t})$. The initial state distribution is assumed to be zero
$\vec{x}=\vec{0}$, the transition dynamics are given by the deterministic
equations of motion above combined with the stochastic disturbances and
parametric uncertainty, and the reward function is given by Equation 13 below.

DDPG utilizes an actor-critic formulation with two neural network function
approximators, one to choose what action to take and the other to decide how
valuable that action was. Figure 5 provides a block diagram of the DDPG agent
system.

Figure 5: Block diagram for actor-critic agent.

The critic function $Q(s,a)$ learns using the Bellman equation given in
Equation 10 where $\gamma$ is the discount factor $0<\gamma<1$ and $\alpha$ is
the learning rate. A discount factor close to zero prioritizes near-term
rewards, while a factor close to one places importance on long-term rewards.
The learning rate is a critical parameter, and is typically kept low to ensure
stable learning.

$Q(s_{t+1},a_{t+1})=Q(s_{t},a_{t})+\alpha[r_{t}$
$+\gamma\underset{a}{max}Q(s_{t+1},a)-Q(s_{t},a_{t})]$ (10)

The actor function $\mu(s|\theta^{\mu})$ specifies the current policy by
mapping states to actions and is updated by applying the chain rule to the
expected return from the initial distribution, $J=\mathbb{E}(r_{1})$, with
respect to the actor parameters as shown in Equation 11. This is the gradient
of the policy’s performance [9], also known as the policy gradient in DDPG.

$\nabla_{\theta,\mu}J=\mathbb{E}[\nabla_{a}Q(s,a|\theta^{Q})|_{s=s_{t},a=\mu(s_{t})}\nabla_{\theta_{\mu}}\mu(s|\theta^{\mu})|_{s=s_{t}}]$
(11)

However, it is not advantageous to always choose the ”best” action. At the
beginning of a training session, the neural networks will be poorly optimized
and will not know what the best actions are. Therefore
$\mu(s_{t}|\theta_{t}^{\mu})$ is augmented by additive noise sampled from a
temporally correlated stochastic process, in this work from the OU process. As
training progresses, the noise scaling term $\epsilon(t)$ is decreased to
transition the agent from ”exploration” of the state space to ”exploitation”
of what it considers to be good actions.

$\mu^{\prime}(s_{t})=\mu(s_{t}|\theta_{t}^{\mu})+\epsilon(t)\mathcal{N}$ (12)

The reward function used in this work is given by Equation 13, where $\Delta
t$ is the simulation time step, $\vec{x}_{i}$ is the $i^{th}$ term of the
state vector, $\vec{x}_{i,des}$ is the desired value for the $i^{th}$ state
(each equal to zero for this application), $F$ is the control force, and
$w_{i}$ is the reward tuning weight for the $i^{th}$ state. This function is
arbitrary, and was chosen specifically for this application. Methods exist to
learn the reward function, known as Inverse Reinforcement Learning (IRL), and
it is possible to simultaneously learn the reward and policy functions through
a coupled approach [10].

$r=-\Delta
t(\sum_{i=1}^{N_{s}}w_{i}(\vec{x}_{i}-\vec{x}_{i,des})^{2}+0.1F^{2})$ (13)

$w=\begin{bmatrix}0.1&0.5&1&1&1.2&1\\\ \end{bmatrix}$

This reward function incentivizes keeping the system state close to the
desired value, in this case with the pendulum vertical, the cart near its
initial horizontal position, and the time derivatives small. The control
effort reward is scaled by 0.1 to ensure that the agent learns that while it
is important to keep the control effort low, it is much more critical to keep
the states close to their desired values. This sum is then multiplied by the
negative of the simulation time increment to both scale the reward and
transform from a cost function (to be minimized) to a reward function (to be
maximized).

This work utilizes the DDPGAgent Python class provided by the open-source
module Keras-RL [3], neural network functionality from TensorFlow [4], and a
custom OpenAI Gym environment [5] to cohere the simulation model to Keras-RL.
The actor neural network uses twenty hidden layers with 256 neurons per layer,
while the critic neural network uses the same number of hidden layers with 512
neurons per layer. The agent was trained for 100,000 steps with 512 batches
per step. Learning parameters are given in Table II in the appendix.

##  III Simulation Results

Due to the stochastic nature of this work, a great number of simulations
utilizing the Monte Carlo approach were conducted to evaluate control system
performance. 10,000 simulations for both the classical proportional-derivative
(PD) and deep reinforcement learning control systems were performed with
rewards calculated per Equation 13, and the average reward and the standard
deviation of reward were extracted. Additionally, 10,000 simulations per
control system were performed without the cart velocity disturbance
$d_{\dot{z}}$ as the system, particularly with PD control, was notably
sensitive to this disturbance. These results are summarized in Table I.

TABLE I: Reward Statistics from 10,000 Simulations Control System | Average | Standard Deviation
---|---|---
PD no $d_{\dot{z}}$ | -111.72 | 183.52
DRL no $d_{\dot{z}}$ | -63.83 | 61.12
PD | -79,005.72 | 132,957.99
DRL | -212.35 | 105.79

Figures 6-12 present data from representative simulations for both control
systems. Figures 6 and 9 show the time history of the pendulum midpoint and
tip angles, Figures 7 and 10 plot the time history of the control force
applied to the cart, and Figures 8 and 11 display the two-dimensional motion
of the cart, pendulum midpoint, and pendulum tip. Figure 12 shows a typical
simulation of the PD system when $d_{\dot{z}}$ is present. Note that
$d_{\dot{z}}$ is present for all DRL plots.

Figure 6: Example time history of pendulum angles for PD control without
$d_{\dot{z}}$.

Figure 7: Example time history of control force for PD control without
$d_{\dot{z}}$.

Figure 8: Example FIPWC motion for PD control without $d_{\dot{z}}$.

Figure 9: Example time history of pendulum angles for DRL control with
$d_{\dot{z}}$.

Figure 10: Example time history of control force for DRL control with
$d_{\dot{z}}$.

Figure 11: Example FIPWC motion for DRL control with $d_{\dot{z}}$.

Figure 12: Example FIPWC motion for PD control with $d_{\dot{z}}$.

##  IV Discussion

It is clear from Table I that the DRL control system is able to keep the
states closer to their desired values, in this case representing the pendulum
vertical, the cart close to its initial position, and the rates low, for
longer compared to the PD control system, particularly when $d_{\dot{z}}$ is
present. The PD controller is wholly incapable of rejecting $d_{\dot{z}}$,
while the DRL controller is able to stabilize the system for all combinations
of disturbance levels as indicated by the low standard deviations in Table I.
The DRL controller achieves nearly double the reward average and half the
standard deviation of the PD controller even when $d_{\dot{z}}$ is neglected
since the system is regularly driven outside of the linear regime, commonly
considered to be $|\phi,\,\theta|<15^{\circ}$.

Figures 6-12 provide a visual insight into the superior performance of the DRL
system. The simulations that produced these data represent a typical result
close to the average performance for each system. Despite the DRL controller
having to compensate for $d_{\dot{z}}$, the pendulum angles are kept below 11
degrees for the entire simulation and below 10 degrees for two thirds of the
run. The sharp increase in the angles from 7-9 seconds can be attributed to a
particularly high-magnitude combination of disturbances, with the controller
beginning to compensate for them as the simulation ends.

Figure 12 demonstrates the PD control system’s difficulties in keeping the
pendulum vertical with $d_{\dot{z}}$ present. The standard deviation for this
set of runs is particularly high because in many simulations the angular rates
reached the DDPG agent’s observation space limit of $\pm 1\times 10^{9}$,
which returns $-1\times 10^{7}$ for the reward. This severely punishes the
agent for exploring this region of the state space, as the OpenAI Gym
environment returns $nan$ when the state is outside of that space which spoils
any training or testing runs.

##  V Conclusions

This work presented a deep reinforcement learning control approach to
stochastic disturbance rejection for a highly nonlinear dynamic system with
parametric uncertainty. Results from 10,000 Monte Carlo simulations
demonstrated superior performance for the DRL control system compared to a
classical proportional-derivative control law.

While the DRL controller produced good results, there are still some
improvements that could be made. The controller is operating at the simulation
rate of 100 Hz, and has no rate limitation as is clear from Figure 10. This is
not very realistic since real actuators would not operate at such a high
frequency, have rate limitations, and be subject to their own dynamic effects
such as rise time, settling time, and overshoot. The DDPGAgent class has a
parameter to decrease the policy update rate from the simulation rate. This
could be used to limit the controller rate to 10 Hz, for example. An actuator
dynamics model could also be incorporated to more accurately represent real-
world systems.

Future work topics of interest are implementing realistic disturbances,
nonlinear flight vehicle models, and coupling the DDPG agent with a recurrent
neural network (RNN) for uncertain parameter estimation and disturbance
observation.

## Appendix

Table II provides the simulation parameter values used in this work.

TABLE II: Simulation Parameter Values Parameter | Value (unit) | Parameter | Value (unit)
---|---|---|---
$m_{t}$ | 0.019 kg | $m_{c}$ | 0.18 kg
$m_{b}$ | 0.0215 kg | $L$ | 0.32 m
$l$ | 0.16 m | $b_{1}$ | 0.001 N-s/m
$b_{2}$ | 0.001 N-s/m | $b_{3}$ | 12 N-s/m
$k_{1}$ | 2 N/m | $k_{2}$ | 8 N/m
$\Delta t$ | 0.01 sec | $g$ | 9.81 $m/s^{2}$
$\kappa_{\dot{z}}$ | 0.01 | $\kappa_{\dot{\phi},\dot{\theta}}$ | 10
$\mu_{\dot{z}}$ | 0 | $\mu_{\dot{\phi},\dot{\theta}}$ | 0
$\sigma_{\dot{z}}$ | 0.1 m/s | $\sigma_{\dot{\phi},\dot{\theta}}$ | 1 deg/sec
$\alpha$ | 0.001 | $\gamma$ | 0.99

## References

  * [1] Mnih, V., Kavukcuoglu, K., Silver, D. et al. Human-level control through deep reinforcement learning. Nature 518, 529–533 (2015). https://doi.org/10.1038/nature14236
  * [2] Lillicrap, T., Hunt, J.J., Pritzel, A., Heess, N., Erez, T., Tassa, Y., Silver, D., and Wierstra, D. (2016). Continuous control with deep reinforcement learning. CoRR, abs/1509.02971.
  * [3] Matthias Plappert. (2016). keras-rl. https://github.com/keras-rl/keras-rl.
  * [4] Martín Abadi, et al. TensorFlow: Large-scale machine learning on heterogeneous systems, 2015\. Software available from tensorflow.org.
  * [5] Brockman, G., Cheung, V., Pettersson, L., Schneider, J., Schulman, J., Tang, J., and Zaremba, W. (2016). OpenAI Gym. ArXiv Preprint ArXiv:1606.01540.
  * [6] S. K. Gorade, P. S. Gandhi and S. R. Kurode, ”Modeling and output feedback control of Flexible Inverted Pendulum on Cart,” 2015 International Conference on Power and Advanced Control Engineering (ICPACE), 2015, pp. 436-440, doi: 10.1109/ICPACE.2015.7274987.
  * [7] Jordan D. Larson, Ryan W. Thomas, and Vincent W. Hill. (2019). GNCPy: A Python library for Guidance, Navigation, and Control algorithms, https://github.com/drjdlarson/gncpy.
  * [8] Uhlenbeck, George E and Ornstein, Leonard S. On the theory of the brownian motion. Physical review, 36(5):823, 1930.
  * [9] Silver, David, Lever, Guy, Heess, Nicolas, Degris, Thomas, Wierstra, Daan, and Riedmiller, Martin. Deterministic policy gradient algorithms. In ICML, 2014.
  * [10] R. W. Thomas and J. D. Larson, ”Inverse Reinforcement Learning for Generalized Labeled Multi-Bernoulli Multi-Target Tracking,” 2021 IEEE Aerospace Conference (50100), 2021, pp. 1-9.


# Normality of smooth statistics for planar determinantal point processes

Antti Haimi<EMAIL_ADDRESS>Faculty of Mathematics, University of
Vienna, Oskar-Morgenstern-Platz 1, A-1090 Vienna, Austria  and  José Luis
Romero<EMAIL_ADDRESS>Faculty of Mathematics, University of
Vienna, Oskar-Morgenstern-Platz 1, A-1090 Vienna, Austria
and Acoustics Research Institute, Austrian Academy of Sciences, Wohllebengasse
12-14 A-1040, Vienna, Austria

###### Abstract.

We consider smooth linear statistics of determinantal point processes on the
complex plane, and their large scale asymptotics. We prove asymptotic
normality in the finite variance case, where Soshnikov’s theorem is not
applicable. The setting is similar to that of Rider and Virág [Electron. J.
Probab., 12, no. 45, 1238–1257, (2007)] for the complex plane, but replaces
analyticity conditions by the assumption that the correlation kernel is
reproducing. Our proof is a streamlined version of that of Ameur, Hedenmalm
and Makarov [Duke Math J., 159, 31–81, (2011)] for eigenvalues of normal
random matrices. In our case, the reproducing property is brought to bear to
compensate for the lack of analyticity and radial symmetries.

###### Key words and phrases:

asymptotic normality, determinantal point process, linear statistics, Weyl-
Heisenberg DPP

###### 2020 Mathematics Subject Classification:

60F05, 60F17, 30H05, 60G55, 60B20, 82B31,46C07

A. H. and J. L. R. gratefully acknowledge support from the Austrian Science
Fund (FWF): Y 1199.

##  1\. Introduction

###  1.1. Determinantal point processes

Determinantal point processes (DPP’s) are certain interacting particle systems
which exhibit repulsion between different points [24, 6, 16, 15]. Their
central feature is that all the statistical information is encoded in a single
function, called correlation kernel. Since their introduction by Macchi in
1975 [17], DDP’s have found extensive applications in areas of knowledge as
varied as random matrix theory, number theory, quantum mechanics and machine
learning.

For a precise definition, we restrict attention to the Euclidean space and
consider a random subset $\mathcal{X}\subset\mathbb{R}^{d}$ encoded as a
random integer-valued positive Radon measure $E\mapsto\mathcal{X}(E)$, which
counts the number of points within a Borel set. We only consider _simple_
processes, that is, almost surely $\mathcal{X}$ assigns mass at most one to
each singleton $\\{x\\}$. Such a process is said to admit _correlation
functions_ (with respect to the Lebesgue measure) if for each
$n\in\mathbb{N}$, there exists a function
$\rho_{n}:(\mathbb{R}^{d})^{n}\to[0,\infty)$ such that for all mutually
disjoint Borel sets $E_{1},\ldots,E_{n}\subset\mathbb{R}^{d}$,

$\mathbb{E}\Big{[}\prod_{j=1}^{n}X(E_{j})\Big{]}=\int_{E_{1}\times\cdots\times
E_{n}}\rho_{n}(x_{1},\ldots,x_{n})\,dx_{1}\ldots dx_{n}.$

In addition, $\rho_{n}(x_{1},\ldots,x_{n})$ is required to vanish if
$x_{j}=x_{k}$ for some $j\not=k$. The correlation functions $\rho_{k}$, also
called _joint intensities_ , describe the essential statistical properties of
$\mathcal{X}$. A simple point process $\mathcal{X}$ is called _determinantal_
if there exists a function
$K:\mathbb{R}^{d}\times\mathbb{R}^{d}\to\mathbb{C}$, called _correlation
kernel_ , such that

$\rho_{n}(x_{1},\ldots,x_{n})=\det K(x_{k},x_{j})_{k,j=1,\ldots,n},\qquad
x_{1},\ldots,x_{n}\in\mathbb{R}^{d}.$

One of the most famous DDP’s is the infinite _Ginibre ensemble_ [11] on
$\mathbb{C}^{r}\simeq\mathbb{R}^{2r}$, which has correlation kernel

(1.1) $\displaystyle
K(z,w)=\frac{1}{\pi^{r}}\exp\Big{[}{\sum_{j=1}^{r}z_{j}\bar{w}_{j}-\frac{1}{2}|z_{j}|^{2}-\frac{1}{2}|w_{j}|^{2}}\Big{]},\qquad
z,w\in\mathbb{C}^{r}.$

The Macchi-Soshnikov theorem [17] provides a machinery to produce examples of
DPP’s: any Hermitian symmetric kernel
$K:\mathbb{R}^{d}\times\mathbb{R}^{d}\to\mathbb{C}$ which induces a locally
trace class integral operator with spectrum in the interval $[0,1]$ is the
correlation kernel of a uniquely determined DPP. A simple way to satisfy the
spectral assumption is to let $\mathcal{H}\subset L^{2}(\mathbb{R}^{d})$ be a
_reproducing kernel subspace_ \- that is, a closed subspace where evaluations
are continuous functionals - and let $K$ be its _reproducing kernel_ , that
is, the integral kernel representing the orthogonal projection
$L^{2}(\mathbb{R}^{d})\to\mathcal{H}$. As such, the kernel satisfies the
_reproducing formula_

(1.2) $\displaystyle K(x,y)=\int K(x,u)K(u,y)\,du.$

The Ginibre kernel (1.1) in one complex dimension $r=1$ corresponds in this
way to the weighted Bargmann-Fock space
$\mathcal{H}=\mathcal{F}_{0}(\mathbb{C})$,

(1.3)
$\displaystyle\mathcal{F}_{1}(\mathbb{C})=\Big{\\{}F(z)e^{-|z|^{2}/2}:F:{\mathbb{C}}\to{\mathbb{C}}\text{
entire, }\int_{\mathbb{C}}|F(z)|^{2}e^{-|z|^{2}}dA(z)<\infty\Big{\\}},$

where $dA$ is the Lebesgue (area) measure. Further examples, important in high
energy physics, are the polyanalytic Ginibre ensembles [13, 22] associated
with the weighted poly-Bargmann spaces [26]
$\mathcal{H}=\mathcal{F}_{k}(\mathbb{C})$,

(1.4)
$\displaystyle\mathcal{F}_{k}(\mathbb{C})=\Big{\\{}F(z)e^{-|z|^{2}/2}:\bar{\partial}^{k}F=0,\int_{\mathbb{C}}|F(z)|^{2}e^{-|z|^{2}}dA(z)<\infty\Big{\\}},$

or with the pure poly-Bargmann spaces
$\mathcal{F}_{k}(\mathbb{C})\ominus\mathcal{F}_{k-1}(\mathbb{C})$ [26].

###  1.2. Linear statistics and asymptotic normality

One of the main topics of interest around DPP’s has been the validity of
central limit theorems for various statistical observations. In this article
we will be interested in the _linear statistic_

$\displaystyle\mathrm{Tr}(f):=\sum_{x\in\mathcal{X}}f(x)$

associated with a _smooth_ function $f:\mathbb{R}^{d}\to\mathbb{C}$ and a DPP
$\mathcal{X}$. Specifically, we will consider large scale asymptotics of
linear statistics

(1.5) $\displaystyle\mathrm{Tr}_{\rho}(f)=\sum_{x\in\mathcal{X}_{\rho}}f(x)$

calculated under increasing contractions of the point process

(1.6)
$\displaystyle\mathcal{X}_{\rho}=\big{\\{}\tfrac{1}{\sqrt{\rho}}x:x\in\mathcal{X}\big{\\}},\qquad\rho\to\infty,$

or, more generally, under qualitatively similar asymptotic deformations. The
corresponding expectations are given by

$\displaystyle\mathbb{E}\big{[}\mathrm{Tr}_{\rho}(f)\big{]}=\int_{\mathbb{R}^{d}}f(x)K_{\rho}(x,x)\,dx,$

where $K_{\rho}$ is the correlation kernel of $\mathcal{X}_{\rho}$. In the
model case (1.6),

(1.7)
$\displaystyle\mathbb{E}\big{[}\mathrm{Tr}_{\rho}(f)\big{]}=\rho^{d/2}\int_{\mathbb{R}^{d}}f(x)K(\sqrt{\rho}x,\sqrt{\rho}x)\,dx,$

where $K$ is the correlation kernel of $\mathcal{X}$.

Soshnikov’s celebrated central limit theorem [25, Theorem 1], which builds on
[7], gives very general conditions under which the standardized linear
statistic

(1.8)
$\displaystyle\widetilde{\mathrm{Tr}_{\rho}}(f):=\frac{\mathrm{Tr}_{\rho}(f)-\mathbb{E}\big{[}\mathrm{Tr}_{\rho}(f)\big{]}}{\sqrt{\mathrm{Var}\big{[}\mathrm{Tr}_{\rho}(f)}\big{]}}$

converges in distribution to a normal variable. The key assumption is that the
expected value and variance of $\mathrm{Tr}_{\rho}(f)$ be related by

(1.9)
$\mathbb{E}\big{[}\mathrm{Tr}_{\rho}|f|\big{]}=\mathrm{O}\big{(}\mathrm{Var}[\mathrm{Tr}_{\rho}(f)]^{\delta}\big{)},\mbox{
as }\rho\to\infty,$

for some $\delta>0$, while $\mathrm{Var}[\mathrm{Tr}_{\rho}(f)]\to\infty$.

Soshnikov’s theorem immediately implies, for example, that the standardized
point counting statistics for the Ginibre ensemble (1.1), where $f$ is the
indicator function of a smooth domain, are asymptotically normal. For smooth
statistics, the Ginibre ensemble (1.1) satisfies the key relation (1.9) only
in complex dimension $r>1$. While standardized smooth linear statistics of the
Ginibre ensemble are also asymptotically normal in one complex dimension, the
planar case requires ad hoc arguments [20].

###  1.3. Asymptotic normality for planar DPP’s

Rider and Virág introduced an abstract model to study the normality of
standardized smooth statistics of DDP’s on certain Riemann surfaces [20]. The
central assumption is an approximate form of analyticity of the correlation
kernel, and applies to the planar Ginibre ensemble, as well as to its
analogues on the sphere and unit disk.

A second technique was introduced by Ameur, Hedenmalm and Makarov to study the
normality of fluctuations of eigenvalues of certain ensembles of random
matrices under a confining potential of increasing strength [3, 2] (see also
Berman’s work in several complex variables [5]). While the corresponding
correlation kernels are analytic, they do not have an explicit formula that
would facilitate the approach of [20, 21].

Beyond the analytic setting, normality of fluctuations of smooth statistics
was shown to hold for certain planar DDPs in [14]. Here, the correlation
kernel is _polyanalytic_ and the key insight is that the cumulants of (1.5)
can be related by precise algebraic formulas to those of the analytic Ginibre
ensemble, essentially reducing the problem to [21].

The first goal of this article is to present a general result on asymptotic
normality of smooth statistics for planar DPP’s with correlation kernels that
satisfy the _reproducing property_ (1.2). The result we derive neither
requires analyticity nor relies on algebraic relations and symmetries linking
the correlation kernel to the Ginibre one. Our main insight is that the
technique from [3] can be streamlined to substitute many of the cancellations
attributed to analyticity and radiality to instead rely on the reproducing
property. The reproducing assumption is to a great extent optimal, since under
mild assumptions non-reproducing kernels fall within the scope of Soshnikov’s
theorem (see Sections 2 and 5).

The result we present can be applied for example to infinite poly-analytic
Ginibre ensembles considered with respect to general potentials, that is,
infinite analogues of the finite ensembles introduced in [13]. These are
related to Hilbert spaces analogous to (1.4) but with respect to general
weights $Q:{\mathbb{C}}\to[0,\infty)$:

(1.10)
$\displaystyle\mathbb{H}=\Big{\\{}F(z)e^{-Q(z)/2}:\bar{\partial}^{k}F=0,\int_{\mathbb{C}}|F(z)|^{2}e^{-Q(z)}\,dA(z)<\infty\Big{\\}},$

and cannot be related to their analytic counterparts by simple algebraic
formulae, as done for standard potentials in [14].

A second motivating application is to the kernels arising from the Schrödinger
representation of the Weyl-Heisenberg group. Let $g:\mathbb{R}\to\mathbb{C}$
be a smooth and fast decaying function satisfying
$\int_{{\mathbb{C}}}|g|^{2}\,dA=1$. Then

(1.11) $\displaystyle
K_{g}(z,w)=\int_{\mathbb{R}}\overline{g(t-x_{1})}g(t-x_{2})e^{2\pi
it(y_{2}-y_{1})}\,dt,\qquad z=x_{1}+iy_{1},w=x_{2}+iy_{2}\in\mathbb{C},$

is a locally trace class reproducing kernel acting on $L^{2}({\mathbb{C}},dA)$
[10, Chapters 1 and 2]. Non-smooth linear statistics for the associated DPP
have been considered in [1, 18], providing variance asymptotics. Under mild
hypotheses, non-smooth statistics satisfy (1.9), and, as with the Ginibre
ensemble, Soshnikov’s theorem readily shows that standardized statistics
corresponding to the number of points within a growing disk are asymptotically
normal (while the investigation of precise rates for such convergence is more
challenging [9]). As an application of our main result, we shall conclude that
Weyl-Heisenberg DPP’s (1.11) also enjoy asymptotically normal _smooth_ linear
statistics. Interestingly, the kernel (1.11) may exhibit no radial symmetries.

##  2\. Results

Let $\mathcal{H}_{\rho}\subset L^{2}({\mathbb{C}},dA)$, $\rho>0$, be a family
of reproducing Hilbert spaces, where $dA$ is the Lebesgue area measure, and
denote by $K_{\rho}:{\mathbb{C}}\times{\mathbb{C}}\to\mathbb{C}$ the
corresponding reproducing kernels. We will assume that the kernels $K_{\rho}$
satisfy the following: there exist continuous _envelope functions_
$\phi_{\rho}:{\mathbb{C}}\to[0,\infty)$ such that

(2.1) $\displaystyle|K_{\rho}(z,w)|\leq\phi_{\rho}(z-w),\qquad
z,w\in{\mathbb{C}}.$

The envelops are further assumed to satisfy the following properties:

  1. $\bullet$

(_Size_)

(2.2)
$\displaystyle\sup_{\rho>0}\sup_{z\in{\mathbb{C}}}\tfrac{1}{\rho}\phi_{\rho}(z)<\infty,$

  2. $\bullet$

(_Uniform integrability_)

(2.3) $\displaystyle\sup_{\rho>0}\int_{{\mathbb{C}}}\phi_{\rho}\,dA<\infty,$

  3. $\bullet$

(_Interaction decay_)

(2.4)
$\displaystyle\lim_{\rho\to\infty}\rho\int_{{\mathbb{C}}}|z|^{3}\phi_{\rho}(z)dA(z)=0.$

The above assumptions mean that the envelops behave qualitatively like
dilations $\phi_{\rho}(z)=\rho\phi(\sqrt{\rho}z)$. The model is essentially as
in [20], but, crucially, does not assume analyticity. On the other hand, here
$K_{\rho}$ is assumed to be reproducing, while [20] requires an asymptotic
reproducing property for polynomials.

######  Remark 2.1.

A positive kernel $K_{\rho}$ with envelope $\phi_{\rho}$ satisfying (2.2) is
_locally trace class_ , that is, for every compact set $B\subset\mathbb{C}$,
the localized kernel $K^{B}_{\rho}(z,w)=1_{B}(z)K_{\rho}(z,w)1_{B}(w)$
represents a trace-class integral operator. Indeed, the operator with integral
kernel $K^{B}_{\rho}$ is positive because it is a compression of that with
kernel $K_{\rho}$, while the corresponding trace is $\lesssim\rho|B|$ by (2.1)
— if $K$ is continuous this follows from
$\int_{\mathbb{C}}K^{B}_{\rho}(z,z)\,dA(z)=\int_{B}K_{\rho}(z,z)\,dA(z)\lesssim\rho|B|$,
see, e.g., [23, Theorem 2.1], while in general one can use an approximation
argument. Hence, if $K_{\rho}$ is also reproducing, by the Soshnikov-Macchi
theorem [17] [24, Theorem 3], it is the correlation kernel of a unique
determinantal point process.

The following is our main result on asymptotic normality.

######  Theorem 2.2.

Let $\\{K_{\rho}:\rho>0\\}$ be a family of Hermitian symmetric reproducing
kernels acting on $L^{2}({\mathbb{C}},dA)$ with continuous envelopes
$\eqref{eq_env}$ satisfying (2.2), (2.3), (2.4). Let $f$ be a compactly
supported, real-valued $C^{3}$ test function. Then the cumulants of the order
$k\geq 3$ of the variable $\mathrm{Tr}_{\rho}(f)$ tend to zero as
$\rho\rightarrow\infty$. As a consequence:

  * (i)

If $\mathrm{Var}\big{[}\mathrm{Tr}_{\rho}(f)\big{]}\to\sigma^{2}$ for some
finite $\sigma\geq 0$, then
$\mathrm{Tr}_{\rho}(f)-\mathbb{E}[\mathrm{Tr}_{\rho}(f)]\to\mathcal{N}(0,\sigma^{2})$
in distribution;

  * (ii)

If $\liminf_{\rho\to\infty}\mathrm{Var}\big{[}\mathrm{Tr}_{\rho}(f)\big{]}>0$,
then $\widetilde{\mathrm{Tr}_{\rho}(}f)\to\mathcal{N}(0,1)$ in distribution.

(In (i), the degenerate case $\sigma=0$ means that the limit variable is
almost surely zero.)

While the setting of Theorem 2.2 resembles the planar case of [20], our proof
builds mainly on [3] and is presented in Section 4.

Theorem 2.2 is applicable for example to general _polyanalytic Ginibre
ensembles_ [13, 12], where $K$ is the reproducing kernel of the space (1.10)
and $Q:\mathbb{C}\to[0,\infty)$ is an adequate _potential function_. Under
mild assumptions on $Q$, a suitable adaptation of the estimates of [13, 12],
which concern finite-particle systems, shows that $K$ satisfies the
assumptions of Theorem 2.2. (While we are not aware of a citable reference for
this fact, we expect to provide the particulars in a forthcoming article.)

Concerning translation invariant processes, we obtain the following
application of Theorem 2.2.

######  Corollary 2.3.

Let $K:{\mathbb{C}}\times{\mathbb{C}}\to{\mathbb{C}}$ be a Hermitian symmetric
reproducing kernel acting on $L^{2}({\mathbb{C}},dA)$. Suppose that

(2.5) $\displaystyle|K(z,w)|=\phi(z-w),\qquad z,w\in{\mathbb{C}},$

for a bounded continuous function $\phi:{\mathbb{C}}\to[0,\infty)$ such that

(2.6) $\displaystyle\int_{\mathbb{C}}(1+|z|^{3})\phi(z)\,dA(z)<\infty.$

Define the dilated kernels

(2.7) $\displaystyle K_{\rho}(z,w)=\rho
K(\sqrt{\rho}z,\sqrt{\rho}w),\qquad\rho>0.$

For a compactly supported, real-valued $C^{3}$ test function $f$ set

(2.8) $\displaystyle\mu_{f}$
$\displaystyle:=\phi(0)\cdot\int_{{\mathbb{C}}}f(z)\,dA(z),$ (2.9)
$\displaystyle\sigma^{2}_{f}$
$\displaystyle:=\frac{1}{2}\int_{\mathbb{C}}\int_{\mathbb{C}}(\nabla f(z)\cdot
w)^{2}\phi(w)^{2}\,dA(z)dA(w).$

Then $\mathrm{Tr}_{\rho}(f)-\rho\cdot\mu_{f}\to N(0,\sigma_{f}^{2})$ in
distribution.

Note that if $\phi(0)\not=0$ Corollary 2.3 implies that (1.9) is not
satisfied. A proof of Corollary 2.3 is presented in Section 5. As an example,
the Weyl-Heisenberg kernel (1.11) satisfies the hypotheses for any Schwartz
function $g\in\mathcal{S}(\mathbb{R})$, and a suitable envelope is provided by

$\displaystyle\phi(z)=\Big{|}\int_{\mathbb{R}}{g(t)}\overline{g(t-x)}e^{-2\pi
ity}\,dt\Big{|},\qquad z=x+iy\in\mathbb{C}.$

Concerning possible generalizations of Corollary 2.3, we remark that an
analogous setting in higher dimension leads to smooth observables with
expectation and variance growing polynomially on $\rho$ (see Remark 5.1).
Similarly, as we show in Proposition 5.2, if we relax the assumptions of
Corollary 2.3 to allow for non-reproducing kernels, then for each such non-
reproducing kernel the expectation and variance of $\mathrm{Tr}_{\rho}(f)$
also grow polynomially in $\rho$. Thus, in these cases, (1.9) is satisfied for
some $\delta>0$ and Soshnikov’s asymptotic normality theorem can be invoked.

Before proceeding to the proofs, we make some final comments. An abstract
setting in which linear statistics can be shown to be asymptotically normal
was recently introduced in [8], where, moreover, quantitative convergence
rates are derived. In the planar case, however, [8, Model 1.2] is not
compatible with the envelope assumptions in Theorem 2.2, and, as in
Soshnikov’s result, $\mathrm{Var}[\mathrm{Tr}_{\rho}(f)]$ needs to diverge as
$\rho\to\infty$ [8, Proposition 4.9].

Second, we remark that Theorem 2.2 can probably be refined to cover kernels
which satisfy the assumptions only near the support of the test function,
provided that (1.2) is replaced by a suitably local approximate formula (while
in analytic settings, reproducing formulas are exactly local). Such an
extension would allow us to also treat finite particle systems in regimes
where points congregate in a compact set (droplet) as long as the test
function is supported away from its boundary. However, the expected gain in
generality did not seem to merit the inclusion of additional technicalities. A
more promising future direction is an extension of Theorem 2.2 to allow for
droplets and test functions that interact with them, paralleling what is known
in the model cases [21, 4, 14].

The reminder of the article is organized as follows: Section 3 clarifies the
notation, Theorem 2.2 is proved in Section 4, while Section 5 contains a proof
of Corollary 2.3 and a discussion on more general settings and the
applicability of Soshnikov’s theorem.

##  3\. Notation

The real and imaginary parts of $z\in\mathbb{C}$ are denoted $\mathrm{Re}(z)$
and $\mathrm{Im}(z)$. We shall also identify
${\mathbb{C}}\simeq\mathbb{R}^{2}$ and write $z=(z^{(1)},z^{(2)})$.
Derivatives with respect to real variables are denoted as follows: for
$F:{\mathbb{C}}^{n}\to\mathbb{C}$, we let
$\nabla_{rj}F(z_{1},\ldots,z_{n})\in{\mathbb{C}}$ be the partial derivative of
$F$ with respect to $z^{(j)}_{r}$, $1\leq r\leq n$, $j=1,2$. For
$z_{0}\in{\mathbb{C}}$, we denote
$\vec{z_{0}}=(z_{0},\ldots,z_{0})\in{\mathbb{C}}^{k}$ if the vector length $k$
is clear from the context.

The differential of the Lebesgue (area) measure on ${\mathbb{C}}$ is denoted
$dA$. We also let $dA$ denote the corresponding product measure on
${\mathbb{C}}^{n}$. For example, depending on convenience and readability we
write interchangeably

$\displaystyle\int_{{\mathbb{C}}^{n}}F\,dA$
$\displaystyle=\int_{{\mathbb{C}}^{n}}F(z_{1},\ldots,z_{n})\,dA(z_{1},\ldots,z_{n})=\int_{{\mathbb{C}}^{n}}F(z_{1},\ldots,z_{n})\,dA(z_{1})\ldots
dA(z_{n})$
$\displaystyle=\int_{z_{1},\ldots,z_{n}\in{\mathbb{C}}}F(z_{1},\ldots,z_{n})\,dA.$

For two non-negative functions $F,G:X\to[0,\infty)$ we write $F\lesssim G$ if
there exists a constant $C>0$ such that $F(x)\leq CG(x)$, for all $x\in X$. We
also use the following asymptotic notation with respect to the parameter
$\rho$: $F(\rho)=O(\rho^{\alpha})$ means that

(3.1) $\displaystyle\limsup_{\rho\to\infty}\rho^{-\alpha}F(\rho)<\infty$

while $F(\rho)=o(\rho^{\alpha})$ means that

(3.2) $\displaystyle\lim_{\rho\to\infty}\rho^{-\alpha}F(\rho)=0.$

For example, (2.4) means that
$\int_{\mathbb{C}}|z|^{3}\phi_{\rho}(z)dA(z)=o(1/\rho)$. When the function $F$
depends on additional parameters besides $\rho$, the limits (3.1) and (3.2)
may or may not hold uniformily on these. Such dependencies are clarified in
each case.

##  4\. Asymptotic normality

The overall strategy is as follows. To show that smooth linear statistics are
asymptotically normal we apply the cumulant method (see, e.g., [19, Appendix
A.3]). We use the exact expressions for cumulants derived in [25] to then
obtain asymptotic expansions by Taylor expanding the test function to second
order, as in [21, 3]. Elaborating on [3], we show that the reproducing
property of the correlation kernel leads to certain significant cancellations,
which allow us to conclude that higher order cumulants vanish asymptotically.

###  4.1. Cumulants

Let us fix a family of reproducing kernels $K_{\rho}$ satisfying (A1), (A2)
and (A3) and a compactly supported $C^{3}$ function
$f:\mathbb{C}\to\mathbb{R}$. Without loss of generality, we can assume that
the kernel envelops are symmetric

$\displaystyle\phi_{\rho}(z)=\phi_{\rho}(-z),\qquad z\in{\mathbb{C}},\rho>0.$

The cyclic products of the correlation kernel $K_{\rho}$ are denoted

$\displaystyle R_{\rho,k}(z_{1},\ldots,z_{k}):=K_{\rho}(z_{1},z_{2})\cdots
K_{\rho}(z_{k-1},z_{k})K_{\rho}(z_{k},z_{1}),\qquad
z_{1},\ldots,z_{k}\in\mathbb{C}.$

These are instrumental to compute the cumulants of the linear statistics
$\mathrm{Tr}_{\rho}(f)$. Indeed, let us associate with $f$ the expression

(4.1) $\displaystyle G_{k}(z_{1},\ldots
z_{k})=\sum_{j=1}^{k}\frac{(-1)^{j-1}}{j}\sum_{\begin{subarray}{c}k_{1}+k_{2}+\ldots+k_{j}=k\\\
k_{1},\ldots,k_{j}\geq 1\end{subarray}}\frac{k!}{k_{1}!k_{2}!\cdots
k_{j}!}\prod_{l=1}^{j}f(z_{l})^{k_{l}},\qquad
z_{1},\ldots,z_{k}\in\mathbb{C}.$

As shown in [25, Lemma 1], the cumulants $C_{\rho,k}(f)$ of
$\mathrm{Tr}_{\rho}(f)$ are

(4.2) $\displaystyle
C_{\rho,k}(f)=\int_{\mathbb{C}^{k}}G_{k}(z_{1},\ldots,z_{k})R_{\rho,k}(z_{1},\ldots,z_{k})\,dA(z_{1})\cdots
dA(z_{k}).$

(This is a simplification of a more general formula in [25, Lemma 1] that
concerns possibly non-reproducing kernels.) By the reproducing property of the
kernel $K_{\rho}$, we can extend integration in (4.2) to a new variable
$z_{0}$ and write

(4.3) $\displaystyle
C_{\rho,k}(f)=\int_{\mathbb{C}^{k+1}}G_{k}(z_{1},\ldots,z_{k})R_{\rho,k+1}(z_{0},z_{1},\ldots,z_{k})\,dA(z_{0})dA(z_{1})\cdots
dA(z_{k}).$

This observation is a main tool in [3] and it is similarly fundamental for our
proof.

###  4.2. Taylor expansions

As in [21, 3] we start by expanding $G_{k}$ up to second order. We write

(4.4)
$G_{k}(z_{1},\ldots,z_{k})=G_{k}(\vec{z_{0}})+{T^{1}_{k}}(z_{0},\ldots,z_{k})+{T^{2}_{k}}(z_{0},\ldots,z_{k})+{T^{*}_{k}}(z_{0},\ldots,z_{k})$

where

(4.5) $\displaystyle\vec{z_{0}}$ $\displaystyle=(z_{0},\ldots,z_{0}),$ (4.6)
$\displaystyle{T^{1}_{k}}(z_{0},\ldots,z_{k})$ $\displaystyle=\sum_{1\leq
r\leq k,1\leq j\leq
2}\nabla_{rj}G_{k}(\vec{z_{0}})\big{(}z_{r}^{(j)}-z_{0}^{(j)}\big{)},$ (4.7)
$\displaystyle{T^{2}_{k}}(z_{0},\ldots,z_{k})$
$\displaystyle=\frac{1}{2}\sum_{1\leq r_{1},r_{2}\leq k,1\leq j_{1},j_{2}\leq
2}\nabla_{r_{1}j_{1}}\nabla_{r_{2}j_{2}}G_{k}(\vec{z_{0}})\big{(}z_{r_{1}}^{(j_{1})}-z_{0}^{(j_{1})}\big{)}\big{(}z_{r_{2}}^{(j_{2})}-z_{0}^{(j_{2})}\big{)}$

and the reminder term ${T^{*}_{k}}(z_{0},\ldots,z_{k})$ satisfies

(4.8)
$\displaystyle|{T^{*}_{k}}(z_{0},\ldots,z_{k})|\lesssim\sum_{j=1}^{k}|z_{j}-z_{0}|^{3},$

with an implied constant that depends on $f$ and $k$, but is uniform in the
center of the expansion $z_{0}$ because $f$ (and hence $G_{k}$) have globally
bounded third order derivatives.

We now inspect the effect of the Taylor approximation in the cumulant
expression (4.3).

######  Lemma 4.1.

Let $B\subset\mathbb{C}$ be a compact set such that
$\mathrm{dist}(\operatorname{supp}(f),B^{c})>0$. Then, for $k\geq 2$,

(4.9)
$C_{\rho,k}(f)=\int_{B\times\mathbb{C}^{k}}\bigg{(}G_{k}(\vec{z_{0}})+{T^{1}_{k}}(z_{0},\hbox
to7.5pt{.\hss.\hss.},z_{k})+{T^{2}_{k}}(z_{0},\hbox
to7.5pt{.\hss.\hss.},z_{k})\bigg{)}\,\times\\\ \,R_{\rho,k+1}(z_{0},\hbox
to7.5pt{.\hss.\hss.},z_{k})\,dA(z_{0},\hbox
to7.5pt{.\hss.\hss.},z_{k})+\mathrm{o}(1),$

where the convergence in the implied limit depends on $f$, $B$, and $k$.

###### Proof.

Step 1. Let $\varepsilon:=\mathrm{dist}(\operatorname{supp}(f),B^{c})>0$, and
fix $k\geq 2$. We first show that

(4.10) $\displaystyle
C_{\rho,k}(f)=\int_{B\times\mathbb{C}^{k}}G_{k}(z_{1},\hbox
to7.5pt{.\hss.\hss.},z_{k})R_{\rho,k+1}(z_{0},\hbox
to7.5pt{.\hss.\hss.},z_{k})\,dA(z_{0},\hbox
to7.5pt{.\hss.\hss.},z_{k})+\mathrm{o}(1).$

Inspecting (4.1) we see that $G_{k}$ is a sum of terms of the form
$f_{1}(z_{1})\cdots f_{k}(z_{k})$, where each $f_{j}$ is a power of $f$, and
thus compactly supported or is identically $1$, and not all $f_{j}$ are $1$.
Hence, it suffices to show that for each $k\in\mathbb{N}$,

(4.11) $\displaystyle\begin{aligned}
&\int_{\mathbb{C}^{k+1}}f_{1}(z_{1})\cdots
f_{k}(z_{k})R_{\rho,k+1}(z_{0},\hbox
to7.5pt{.\hss.\hss.},z_{k})\,dA(z_{0},\hbox to7.5pt{.\hss.\hss.},z_{k})\\\
&\quad=\int_{B\times\mathbb{C}^{k}}f_{1}(z_{1})\cdots
f_{k}(z_{k})R_{\rho,k+1}(z_{0},\hbox
to7.5pt{.\hss.\hss.},z_{k})\,dA(z_{0},\hbox
to7.5pt{.\hss.\hss.},z_{k})+\mathrm{o}(1).\end{aligned}$

Moreover, we can assume that $f_{1}$ is compactly supported; otherwise, it
would be identically $1$ and we would eliminate the variable $z_{1}$ using the
reproducing property of the kernel $K_{\rho}$, and relabel the rest of the
points to obtain another instance of (4.11) with a smaller value of $k$.
(Recall that not all $f_{j}$ are $1$.)

Let us now estimate the error term corresponding to (4.11). We exploit the
support of $f_{1}$ and bound the cyclic product
$R_{\rho,k+1}(z_{0},z_{1},\ldots,z_{k})$ applying (2.2) to the factor
$K_{\rho}(z_{1},z_{2})$ and (2.3) to the others:

$\displaystyle
E_{k}:=\bigg{|}\int_{B^{c}\times{\mathbb{C}}^{k}}f_{1}(z_{1})\hbox
to7.5pt{$\cdot$\hss$\cdot$\hss$\cdot$}f_{k}(z_{k})R_{\rho,k+1}(z_{0},\hbox
to7.5pt{.\hss.\hss.},z_{k})\,dA(z_{0},\hbox
to7.5pt{.\hss.\hss.},z_{k})\bigg{|}$
$\displaystyle\,\lesssim\rho\int_{z_{0}\in B^{c}}\int_{z_{1}\in
B,\,|z_{1}-z_{0}|\geq\varepsilon}\int_{z_{2},\ldots,z_{k}\in{\mathbb{C}}}\phi_{\rho}(z_{1}-z_{0})\phi_{\rho}(z_{3}-z_{2})\hbox
to7.5pt{$\cdot$\hss$\cdot$\hss$\cdot$}\phi_{\rho}(z_{k}-z_{k-1})\phi_{\rho}(z_{0}-z_{k})\,dA$
$\displaystyle\,\lesssim\rho\int_{z_{0}\in{\mathbb{C}}}\int_{z_{1}\in
B}\int_{z_{2},\ldots,z_{k}\in{\mathbb{C}}}|z_{1}-z_{0}|^{3}\phi_{\rho}(z_{1}-z_{0})\phi_{\rho}(z_{3}-z_{2})\hbox
to7.5pt{$\cdot$\hss$\cdot$\hss$\cdot$}\phi_{\rho}(z_{k}-z_{k-1})\phi_{\rho}(z_{0}-z_{k})\,dA$
$\displaystyle\,=\rho\int_{z_{1}\in
B}\int_{z_{0}\in{\mathbb{C}}}|z_{1}-z_{0}|^{3}\phi_{\rho}(z_{1}-z_{0})\int_{z_{2},\hbox
to3.14998pt{.\hss.\hss.},z_{k}\in{\mathbb{C}}}\phi_{\rho}(z_{3}-z_{2})\hbox
to7.5pt{$\cdot$\hss$\cdot$\hss$\cdot$}\phi_{\rho}(z_{k}-z_{k-1})\phi_{\rho}(z_{0}-z_{k})\,dA\,dA(z_{0},z_{1})$

The integral on $z_{2},\ldots,z_{k}$ is seen to be $O(1)$ by integrating in
the order $z_{2},\ldots,z_{k}$ and applying (2.3). The remaining expression
can be bounded by means of (2.4):

$\displaystyle E_{k}\lesssim\rho\int_{z_{1}\in
B}\int_{z_{0}\in{\mathbb{C}}}|z_{1}-z_{0}|^{3}\phi_{\rho}(z_{1}-z_{0})\,dA(z_{0})\,dA(z_{1})\lesssim\rho\cdot\mathrm{o}(1/\rho)=\mathrm{o}(1).$

Step 2. By the Taylor expansion (4.4) and (4.8), it suffices to prove that for
all $j=1,\ldots,k$,

$\int_{B\times{\mathbb{C}}^{k}}|z_{j}-z_{0}|^{3}R_{\rho,k+1}(z_{0},\hbox
to7.5pt{.\hss.\hss.},z_{k})\,dA(z_{0},\hbox
to7.5pt{.\hss.\hss.},z_{k})=\mathrm{o}(1).$

We use the bound $|z_{j}-z_{0}|^{3}\lesssim\sum_{n=1}^{j}|z_{n}-z_{n-1}|^{3}$.
We suppose first that $j<k$ and estimate
$|K_{\rho}(z_{j},z_{j+1})|\lesssim\rho$ by (2.2), and the other factors of
$R_{\rho,k+1}(z_{0},\ldots,z_{k})$ by their envelopes to obtain

$\displaystyle E^{\prime}_{j}$
$\displaystyle:=\bigg{|}\int_{B\times{\mathbb{C}}^{k}}|z_{j}-z_{j-1}|^{3}R_{\rho,k+1}(z_{0},\hbox
to7.5pt{.\hss.\hss.},z_{k})\,dA(z_{0},\hbox
to7.5pt{.\hss.\hss.},z_{k})\bigg{|}$ $\displaystyle\lesssim\rho\int_{z_{0}\in
B}\int_{z_{j}\in\mathbb{C}}|z_{j}-z_{j-1}|^{3}\phi_{\rho}(z_{j}-z_{j-1})\int\phi_{\rho}(z_{0}-z_{k})\prod_{\stackrel{{\scriptstyle
n=1,\ldots,k}}{{n\not=j,j+1}}}\phi_{\rho}(z_{n}-z_{n-1})\,dA\,dA(z_{j})\,dA(z_{0}),$

where the innermost integral runs over all $z_{n}\in\mathbb{C}$ with $1\leq
n\leq k$ and $n\not=j$. By (2.3), if we integrate in the order
$z_{j+1},\ldots,z_{k}$ we obtain

$\displaystyle\int_{z_{k}}\cdots\int_{z_{j+1}}\phi_{\rho}(z_{0}-z_{k})\prod_{\stackrel{{\scriptstyle
n=1,\ldots,k}}{{n\not=j,j+1}}}\phi_{\rho}(z_{n}-z_{n-1})\,dA\lesssim\prod_{n=1}^{j-1}\phi_{\rho}(z_{n}-z_{n-1}).$

Second, integration on $z_{j}$ brings the factor

$\displaystyle\rho\int_{z_{j}\in\mathbb{C}}|z_{j}-z_{j-1}|^{3}\phi_{\rho}(z_{j}-z_{j-1})\,dA(z_{j})\lesssim\rho\cdot\mathrm{o}(1/\rho)=\mathrm{o}(1).$

Finally, by (2.3), integration in the order $z_{j-1},\ldots,z_{1}$ yields

$\displaystyle E^{\prime}_{j}\lesssim\mathrm{o}(1)\cdot\int_{z_{0}\in
B}\int_{z_{1}}\cdots\int_{z_{j-1}}\prod_{n=1}^{j-1}\phi_{\rho}(z_{n}-z_{n-1})\,dA\lesssim\mathrm{o}(1)\cdot|B|=\mathrm{o}(1).$

The case $j=k$ is similar: we bound $|K_{\rho}(z_{0},z_{k})|\lesssim\rho$ and
proceed as before by estimating integrals in the order $z_{k},\ldots,z_{0}$. ∎

###  4.3. Symmetry properties of cumulants

We now collect various symmetry properties of the functions $G_{k}$ in the
following lemma, which elaborates on [3, Lemmas 3.1, 3.2, 3.3] by recording
certain additional symmetries implicit in the proof of those results.

######  Lemma 4.2.

The following holds for all $z_{0}\in{\mathbb{C}}$ and $k\geq 3$:

  1. (i)

$G_{k}(\vec{z_{0}})=0$,

  2. (ii)

$\sum_{1\leq r\leq k}(\nabla_{rj}G_{k})(\vec{z_{0}})=0$ for all $1\leq j\leq
2$,

  3. (iii)

$\sum_{1\leq r_{1},r_{2}\leq
k}(\nabla_{r_{1}m}\nabla_{r_{2}l}G_{k})(\vec{z_{0}})=0$ for all $1\leq m,l\leq
2$,

  4. (iv)

$\sum_{1\leq r\leq k}(\nabla_{rm}\nabla_{rl}G_{k})(\vec{z_{0}})=0$ for all
$1\leq l,m\leq 2$,

  5. (v)

$\sum_{1\leq r_{1}\neq r_{2}\leq
k}(\nabla_{r_{1}m}\nabla_{r_{2}l}G_{k})(\vec{z_{0}})=0$, for all $1\leq
m,l\leq 2$.

###### Proof.

The proof of (i) is exactly as in [3, Lemma 3.1]. Parts (ii) and (iii) follow
by taking derivatives on (i). We now turn to (iv), that is a refined version
of [3, Lemma 3.3] which shows that
$\Delta_{\mathbb{R}^{2k}}G_{k}(\vec{z_{0}})=0$.

Let $k_{1},\ldots,k_{j}$ be positive integers such that
$k_{1}+\cdots+k_{j}=k$. Then for $1\leq r\leq j$, $1\leq m,l\leq 2$ and
$z_{1},\ldots,z_{k}\in{\mathbb{C}}$ we have

(4.12) $\displaystyle\begin{aligned}
&\nabla_{rm}\nabla_{rl}\bigg{(}\prod_{n=1}^{j}f(z_{n})^{k_{n}}\bigg{)}\\\
&\quad=\prod_{1\leq n\leq j,n\neq
r}f(z_{n})^{k_{n}}\bigg{(}k_{r}(k_{r}-1)f(z_{r})^{k_{r}-2}\nabla_{m}f(z_{r})\nabla_{l}f(z_{r})+k_{r}f(z_{r})^{k_{r}-1}\nabla_{m}\nabla_{l}f(z_{r})\bigg{)}.\end{aligned}$

This gives

$\displaystyle\sum_{1\leq r\leq
k}(\nabla_{rm}\nabla_{rl}G_{k})(\vec{z_{0}})=f(z_{0})^{k-2}\nabla_{m}f(z_{0})\nabla_{l}f(z_{0})\cdot
S_{k}+f(z_{0})^{k-1}\nabla_{m}\nabla_{l}f(z_{0})\cdot S^{\prime}_{k},$

where

$\displaystyle S_{k}$
$\displaystyle=\sum_{j=1}^{k}\frac{(-1)^{j-1}}{j}\sum_{\begin{subarray}{c}k_{1}+k_{2}+\ldots+k_{j}=k\\\
k_{1},\ldots,k_{j}\geq
1\end{subarray}}\frac{k!(k_{1}(k_{1}-1)+\cdots+k_{j}(k_{j}-1))}{k_{1}!\cdots
k_{j}!},$ $\displaystyle S^{\prime}_{k}$
$\displaystyle=\sum_{j=1}^{k}\frac{(-1)^{j-1}}{j}\sum_{\begin{subarray}{c}k_{1}+k_{2}+\ldots+k_{j}=k\\\
k_{1},\ldots,k_{j}\geq 1\end{subarray}}\frac{k\cdot k!}{k_{1}!\cdots k_{j}!}.$

$S_{k}$ is shown to be zero in the proof of [3, Lemma 3.3] (see Eqs. (3.5) and
(3.9) in [3]), while $S^{\prime}_{k}$ shown to be zero in [3, Lemma 3.1 and
Eq. (3.1)].

Finally, (v) follows by combining (iii) and (iv). ∎

###  4.4. Proof of Theorem 2.2

Step 1. We invoke Lemma 4.1 with an adequate choice of $B$ and use the
representation (4.9). First, note that, by Lemma 4.2,

(4.13) $\displaystyle G_{k}(\vec{z_{0}})=0.$

Let us show that

(4.14) $\displaystyle\int_{B\times\mathbb{C}^{k}}{T^{1}_{k}}(z_{0},\hbox
to7.5pt{.\hss.\hss.},z_{k})R_{\rho,k+1}(z_{0},\hbox
to7.5pt{.\hss.\hss.},z_{k})\,dA(z_{0},\hbox to7.5pt{.\hss.\hss.},z_{k})=0.$

We inspect (4.6), fix $j\in\\{1,2\\}$ and use the reproducing property of the
kernel $K_{\rho}$ to compute

$\displaystyle\sum_{r=1}^{k}\int_{B\times\mathbb{C}^{k}}\nabla_{rj}G_{k}(\vec{z_{0}})\big{(}z_{r}^{(j)}-z_{0}^{(j)}\big{)}R_{\rho,k+1}(z_{0},\ldots,z_{k})\,dA(z_{0})\cdots
dA(z_{k})$
$\displaystyle\quad=\sum_{r=1}^{k}\int_{B}\int_{\mathbb{C}}\nabla_{rj}G_{k}(\vec{z_{0}})\big{(}z_{r}^{(j)}-z_{0}^{(j)}\big{)}R_{\rho,2}(z_{0},z_{r})\,dA(z_{r})dA(z_{0})$
$\displaystyle\quad=\int_{B\times\mathbb{C}}\sum_{r=1}^{k}\nabla_{rj}G_{k}(\vec{z_{0}})\big{(}z^{(j)}-z_{0}^{(j)}\big{)}R_{\rho,2}(z_{0},z)\,dA(z)dA(z_{0}).$

The sum on $r$ vanishes by Lemma 4.2, which proves (4.14).

Step 2. We turn to the second order terms ${T^{2}_{k}}$ given by (4.7). We
separate the terms for which $r_{1}=r_{2}$ and classify the rest according to
$(j_{1},j_{2})$:

(4.15) $\displaystyle\int_{B\times\mathbb{C}^{k}}{T^{2}_{k}}(z_{0},\hbox
to7.5pt{.\hss.\hss.},z_{k})R_{\rho,k+1}(z_{0},\hbox
to7.5pt{.\hss.\hss.},z_{k})\,dA(z_{0},\hbox
to7.5pt{.\hss.\hss.},z_{k})=\frac{1}{2}\big{[}A_{0}+A_{1,1}+A_{1,2}+A_{2,1}+A_{2,2}\big{]},$

where

$\displaystyle A_{0}$
$\displaystyle=\int_{B\times\mathbb{C}^{k}}\sum_{\begin{subarray}{c}1\leq
r\leq k\\\ 1\leq j_{1},j_{2}\leq
2\end{subarray}}\nabla_{rj_{1}}\nabla_{rj_{2}}G_{k}(\vec{z_{0}})\big{(}z_{r}^{(j_{1})}-z_{0}^{(j_{1})}\big{)}\big{(}z_{r}^{(j_{2})}-z_{0}^{(j_{2})}\big{)}R_{\rho,k+1}(z_{0},\hbox
to7.5pt{.\hss.\hss.},z_{k})\,dA(z_{0},\hbox to7.5pt{.\hss.\hss.},z_{k}),$
$\displaystyle A_{j_{1},j_{2}}$
$\displaystyle=\int_{B\times\mathbb{C}^{k}}\sum_{\begin{subarray}{c}1\leq
r_{1},r_{2}\leq k\\\ r_{1}\neq
r_{2}\end{subarray}}\nabla_{r_{1}j_{1}}\nabla_{r_{2}j_{2}}G_{k}(\vec{z_{0}})\big{(}z_{r_{1}}^{(j_{1})}-z_{0}^{(j_{1})}\big{)}\big{(}z_{r_{2}}^{(j_{2})}-z_{0}^{(j_{2})}\big{)}R_{\rho,k+1}(z_{0},\hbox
to7.5pt{.\hss.\hss.},z_{k})\,dA(z_{0},\hbox to7.5pt{.\hss.\hss.},z_{k}).$

Step 3. We use the reproducing formula and Lemma 4.2 to compute

$\displaystyle A_{0}$ $\displaystyle=\sum_{\begin{subarray}{c}1\leq r\leq k\\\
1\leq j_{1},j_{2}\leq
2\end{subarray}}\int_{B}\int_{\mathbb{C}}\nabla_{rj_{1}}\nabla_{rj_{2}}G_{k}(\vec{z_{0}})\big{(}z_{r}^{(j_{1})}-z_{0}^{(j_{1})}\big{)}\big{(}z_{r}^{(j_{2})}-z_{0}^{(j_{2})}\big{)}R_{\rho,2}(z_{0},z_{r})\,dA(z_{r})dA(z_{0})$
$\displaystyle=\sum_{1\leq j_{1},j_{2}\leq
2}\int_{B}\int_{\mathbb{C}}\sum_{{1\leq r\leq
k}}\nabla_{rj_{1}}\nabla_{rj_{2}}G_{k}(\vec{z_{0}})\big{(}z^{(j_{1})}-z_{0}^{(j_{1})}\big{)}\big{(}z^{(j_{2})}-z_{0}^{(j_{2})}\big{)}R_{\rho,2}(z,z_{0})\,dA(z)dA(z_{0})=0.$

Step 4. Similarly, for $j\in\\{1,2\\}$,

(4.16) $\displaystyle A_{jj}$ $\displaystyle=\sum_{\begin{subarray}{c}1\leq
r_{1},r_{2}\leq k\\\
r_{1}<r_{2}\end{subarray}}\int_{B}\nabla_{r_{1}j}\nabla_{r_{2}j}G_{k}(\vec{z_{0}})\int_{\mathbb{C}^{2}}\big{(}z_{r_{1}}^{(j)}-z_{0}^{(j)}\big{)}\big{(}z_{r_{2}}^{(j)}-z_{0}^{(j)}\big{)}R_{\rho,3}(z_{0},z_{r_{1}},z_{r_{2}})\,dA(z_{r_{1}},z_{r_{2}})dA(z_{0})\,+$
(4.17) $\displaystyle\,\sum_{\begin{subarray}{c}1\leq r_{1},r_{2}\leq k\\\
r_{1}>r_{2}\end{subarray}}\int_{B}\nabla_{r_{1}j}\nabla_{r_{2}j}G_{k}(\vec{z_{0}})\int_{\mathbb{C}^{2}}\big{(}z_{r_{1}}^{(j)}-z_{0}^{(j)}\big{)}\big{(}z_{r_{2}}^{(j)}-z_{0}^{(j)}\big{)}R_{\rho,3}(z_{0},z_{r_{2}},z_{r_{1}})\,dA(z_{r_{1}},z_{r_{2}})dA(z_{0}).$

The inner integrals in (4.16) and (4.17) are equal and independent of $r_{1}$
and $r_{2}$. Thus, denoting by $I$ their common value, Lemma 4.2 gives

$\displaystyle A_{jj}=I\cdot\sum_{\begin{subarray}{c}1\leq r_{1},r_{2}\leq
k\\\ r_{1}\neq
r_{2}\end{subarray}}\nabla_{r_{1}j}\nabla_{r_{2}j}G_{k}(\vec{z_{0}})=0.$

Step 5. For the cross terms $A_{21}=A_{12}$ we invoke again the reproducing
formula,

(4.18) $\displaystyle A_{12}$ $\displaystyle=\sum_{\begin{subarray}{c}1\leq
r_{1},r_{2}\leq k\\\
r_{1}<r_{2}\end{subarray}}\int_{B}\nabla_{r_{1}1}\nabla_{r_{2}2}G_{k}(\vec{z_{0}})\int_{\mathbb{C}^{2}}\big{(}z_{r_{1}}^{(1)}-z_{0}^{(1)}\big{)}\big{(}z_{r_{2}}^{(2)}-z_{0}^{(2)}\big{)}R_{\rho,3}(z_{0},z_{r_{1}},z_{r_{2}})\,dA(z_{r_{1}},z_{r_{2}})dA(z_{0})\,+$
(4.19) $\displaystyle\,\sum_{\begin{subarray}{c}1\leq r_{1},r_{2}\leq k\\\
r_{1}>r_{2}\end{subarray}}\int_{B}\nabla_{r_{1}1}\nabla_{r_{2}2}G_{k}(\vec{z_{0}})\int_{\mathbb{C}^{2}}\big{(}z_{r_{1}}^{(1)}-z_{0}^{(1)}\big{)}\big{(}z_{r_{2}}^{(2)}-z_{0}^{(2)}\big{)}R_{\rho,3}(z_{0},z_{r_{2}},z_{r_{1}})\,dA(z_{r_{1}},z_{r_{2}})dA(z_{0}),$

and note that, by the Hermitian symmetry of the kernel $K_{\rho}$, the inner
integral in (4.18) — denoted by $I^{\prime}$ — is the complex conjugate of the
inner integral in (4.19). Thus,

(4.20)
$\displaystyle\mathrm{Re}(A_{21})=\mathrm{Re}(A_{12})=\mathrm{Re}(I^{\prime})\cdot\sum_{\begin{subarray}{c}1\leq
r_{1},r_{2}\leq k\\\ r_{1}\neq
r_{2}\end{subarray}}\nabla_{r_{1}1}\nabla_{r_{2}2}G_{k}(\vec{z_{0}})=0.$

We combine the previous calculations with (4.15) to conclude

(4.21)
$\displaystyle\mathrm{Re}\Big{[}\int_{B\times\mathbb{C}^{k}}{T^{2}_{k}}(z_{0},\ldots,z_{k})R_{\rho,k+1}(z_{0},\ldots,z_{k})\,dA(z_{0})\ldots
dA(z_{k})\Big{]}=0.$

Step 6. Since $f$ is real-valued, so is the cumulant (4.3). Hence, (4.9)
together with (4.13), (4.14), (4.20), and (4.21) yield

$\displaystyle C_{\rho,k}=\mathrm{Re}(C_{\rho,k})=o(1),$

as claimed.

Towards proving (i), suppose that
$\mathrm{Var}[\mathrm{Tr}_{\rho}(f)]\longrightarrow\sigma^{2}>0$. To conclude
the desired asymptotic normality, the cumulant method requires the uniform
integrability of the random variables
$\mathrm{Tr}_{\rho}(f)-\mathbb{E}[\mathrm{Tr}_{\rho}(f)]$ (see, e.g., [19,
Appendix A.3]). This is the case because their variances are uniformly
bounded. In the degenerate case $\sigma=0$,
$\mathrm{Tr}_{\rho}(f)-\mathbb{E}[\mathrm{Tr}_{\rho}(f)]\to 0$ in probability,
and thus in distribution.

Similarly, suppose that
$\liminf_{\rho\to\infty}\mathrm{Var}[\mathrm{Tr}_{\rho}(f)]>0$. It is easy to
see that, due to the enveloping conditions,
$\mathrm{Var}[\mathrm{Tr}_{\rho}(f)]=O(1)$; see, for example, the proof of
Corollary 2.3 below. Let
$\tilde{f}_{\rho}:=\mathrm{Var}[\mathrm{Tr}_{\rho}(f)]^{-1}f$. Then,
inspecting (4.1), we see that the cumulants $\widetilde{C}_{\rho,k}$ of
$\mathrm{Tr}_{\rho}(\tilde{f}_{\rho})$ satisfy, for $k\geq 3$,

$\displaystyle\widetilde{C}_{\rho,k}=\mathrm{Var}[\mathrm{Tr}_{\rho}(f)]^{-k}\cdot
C_{\rho,k}(f)=o(1),$

and we conclude by the cumulant method that
$\mathrm{Tr}_{\rho}(\tilde{f}_{\rho})$ is asymptotically normal. This proves
(ii). ∎

##  5\. The translation invariant case

We now specialize Theorem 2.2 to dilated kernels with translation invariant
correlations.

###### Proof of Corollary 2.3.

Let $K$ satisfy the assumptions of Corollary 2.3. Without loss of generality,
we assume that $\phi$ is symmetric. The expectation of $\mathrm{Tr}_{\rho}(f)$
is given by (1.7) and reduces to $\rho\mu_{f}$, cf. (2.8). Since $K$ is
reproducing,

$\displaystyle\mathrm{Var}[\mathrm{Tr}_{\rho}(f)]=\frac{1}{2}\int_{\mathbb{C}}\int_{\mathbb{C}}\big{[}f(z)-f(w)\big{]}^{2}|K_{\rho}(z,w)|^{2}\,dA(w)dA(z),$

as can be easily checked by inspecting (4.2) with $k=2$. Let
$B\subset{\mathbb{C}}$ be a compact set such that
$\varepsilon:=\mathrm{dist}(\operatorname{supp}(f),B^{c})>0$. We split the
variance as

(5.1) $\displaystyle\begin{aligned}
\mathrm{Var}[\mathrm{Tr}_{\rho}(f)]&=\frac{1}{2}\int_{B}\int_{\mathbb{C}}\big{[}f(z)-f(w)\big{]}^{2}|K_{\rho}(z,w)|^{2}\,dA(w)dA(z)\\\
&\qquad+\frac{1}{2}\int_{B^{c}}\int_{\mathbb{C}}f(w)^{2}|K_{\rho}(z,w)|^{2}\,dA(w)dA(z).\end{aligned}$

For the second term, note that when $w\in\operatorname{supp}(f)$ and $z\in
B^{c}$, $|z-w|\geq\varepsilon$. Thus, using (2.6), and
$\phi^{2}\leq\|\phi\|_{\infty}\cdot\phi$,

(5.2) $\displaystyle\begin{aligned}
&\int_{B^{c}}\int_{\mathbb{C}}f(w)^{2}|K_{\rho}(z,w)|^{2}\,dA(w)dA(z)\\\
&\qquad\leq\varepsilon^{-3}\|f\|^{2}_{\infty}\int_{B^{c}}\int_{\operatorname{supp}(f)}|z-w|^{3}|K_{\rho}(z,w)|^{2}\,dA(w)dA(z)\\\
&\qquad\lesssim\rho^{2}\int_{B^{c}}\int_{\operatorname{supp}(f)}|z-w|^{3}\phi(\sqrt{\rho}(z-w))\,dA(w)dA(z)\\\
&\qquad=\rho^{-1/2}\int_{\operatorname{supp}(f)}\rho\int_{B^{c}}|\sqrt{\rho}(z-w)|^{3}\phi(\sqrt{\rho}(z-w))\,dA(z)dA(w)\lesssim\rho^{-1/2}.\end{aligned}$

Taylor expanding $f$ to order 1 we see that

$\displaystyle[f(w)-f(z)]^{2}=[\nabla f(z)\cdot(w-z)]^{2}+E(z,w),$

where $|E(z,w)|\leq C|z-w|^{3}$ since $f\in C^{3}$ and compactly supported.
Thus

(5.3) $\displaystyle\begin{aligned}
&\int_{B}\int_{\mathbb{C}}|E(z,w)||K_{\rho}(z,w)|^{2}\,dA(w)dA(z)\\\
&\qquad\lesssim\rho^{2}\int_{B}\int_{\mathbb{C}}|z-w|^{3}\phi(\sqrt{\rho}(z-w))\,dA(w)dA(z)=O(\rho^{-1/2}).\end{aligned}$

Combining (5.1), (5.2) and (5.3) we obtain

$\displaystyle\mathrm{Var}[\mathrm{Tr}_{\rho}(f)]$
$\displaystyle=\frac{1}{2}\int_{B}\int_{\mathbb{C}}\big{[}\nabla
f(z)\cdot(w-z)\big{]}^{2}|K_{\rho}(z,w)|^{2}\,dA(w)dA(z)+O(\rho^{-1/2})$
$\displaystyle=\frac{1}{2}\rho^{2}\int_{B}\int_{\mathbb{C}}\big{[}\nabla
f(z)\cdot(w-z)\big{]}^{2}\phi(\sqrt{\rho}(z-w))^{2}\,dA(w)dA(z)+O(\rho^{-1/2})$
$\displaystyle=\frac{1}{2}\rho^{2}\int_{B}\int_{\mathbb{C}}\big{[}\nabla
f(z)\cdot u\big{]}^{2}\phi(\sqrt{\rho}u)^{2}\,dA(u)dA(z)+O(\rho^{-1/2})$
$\displaystyle=\frac{1}{2}\int_{B}\rho\int_{\mathbb{C}}\big{[}\nabla
f(z)\cdot(\sqrt{\rho}u)\big{]}^{2}\phi(\sqrt{\rho}u)^{2}\,dA(u)dA(z)+O(\rho^{-1/2})$
$\displaystyle=\frac{1}{2}\int_{\mathbb{C}}\int_{\mathbb{C}}\big{[}\nabla
f(z)\cdot u\big{]}^{2}\phi(u)^{2}\,dA(u)dA(z)+O(\rho^{-1/2}).$

We now invoke Theorem 2.2 (part (i)). Suitable envelopes are provided by
$\phi_{\rho}(z)=\rho\phi(\sqrt{\rho}z)$. ∎

######  Remark 5.1.

If we consider a reproducing kernel $K$ acting on $L^{2}(\mathbb{R}^{d})$,
satisfying $|K(z,w)|=\phi(z-w)$, and dilate it by
$K_{\rho}(z,w)=\rho^{d/2}K\big{(}\sqrt{\rho}(z-w)\big{)}$, then under similar
assumptions as in Corollary 2.3 reinspection of the proof of Corollary 2.3
shows

(5.4)
$\displaystyle\mathrm{Var}[\mathrm{Tr}_{\rho}(f)]=\frac{1}{2}\rho^{d/2-1}\int_{\mathbb{R}^{d}}\int_{\mathbb{R}^{d}}\big{[}\nabla
f(z)\cdot u\big{]}^{2}\phi(u)^{2}\,dA(u)dA(z)+O(\rho^{(d-3)/2}).$

Thus, for $d>2$, whenever the integral in (5.4) does not vanish,
$\mathrm{Var}[\mathrm{Tr}_{\rho}(f)]$ is of order $\rho^{d/2-1}$ for large
$\rho$. On the other hand,
$\mathbb{E}[\mathrm{Tr}_{\rho}(|f|)]=\rho^{d/2}\cdot\phi(0)\cdot\int|f|$.
Thus, (1.9) is satisfied for some $\delta>0$, and Soshnikov’s asymptotic
normality theorem is applicable.

Similarly, we now note that if we relax the assumptions of Corollary 2.3 to
allow for non-reproducing kernels, then each such kernel leads to smooth
observables with variance of order at least $\rho$ (and therefore falls within
the scope of Soshnikov’s normality theorem).

######  Proposition 5.2.

Let $K:{\mathbb{C}}\times{\mathbb{C}}\to{\mathbb{C}}$ be a continuous
Hermitian symmetric kernel acting boundedly on $L^{2}({\mathbb{C}},dA)$ and
satisfying (2.5) with $\phi$ continuous and bounded. Instead of assuming the
reproducing property, suppose only that the spectrum of $T_{K}$, the integral
operator with kernel $K$, satisfies

(5.5) $\displaystyle\sigma(T_{K})\subset[0,1],$

so that there exists a unique DPP associated with $K$. Define the dilated
kernels $K_{\rho}$ by (2.7).

Suppose that $K$ is not reproducing (i.e., $T_{K}\not=T_{K}^{2}$). Then, for
every non-zero compactly supported $C^{3}$ test function $f$,

(5.6)
$\displaystyle\liminf_{\rho\to\infty}\frac{1}{\rho}\mathrm{Var}[\mathrm{Tr}_{\rho}(f)]>0.$

###### Proof.

First note that the operator $T_{K}$ is locally trace class by assumption
(2.5). This fact along with (5.5) guarantee the existence of a unique DPP
associated with $K$ (see Remark 2.1).

Let $f$ be a non-zero compactly supported $C^{3}$ test function. Taking into
account that $K$ is not reproducing, we write
$\mathrm{Var}[\mathrm{Tr}_{\rho}(f)]$ as

$\displaystyle\frac{1}{2}\int_{{\mathbb{C}}^{2}}\big{[}f(z)-f(w)\big{]}^{2}|K_{\rho}(z,w)|^{2}\,dA(z,w)+\int_{\mathbb{C}}f(z)^{2}\Big{[}K_{\rho}(z,z)-\int_{\mathbb{C}}|K_{\rho}(z,w)|^{2}\,dA(w)\Big{]}\,dA(z)$
$\displaystyle\qquad\geq\int_{\mathbb{C}}f(z/\sqrt{\rho})^{2}\Big{[}K(z,z)-\int_{\mathbb{C}}|K(z,w)|^{2}\,dA(w)\Big{]}\,dA(z)$
$\displaystyle\qquad=\rho\cdot\int_{\mathbb{C}}f^{2}\,dA\cdot\big{[}\phi(0)-\int_{\mathbb{C}}\phi^{2}\,dA\big{]},$

where the first line is a general fact [25, Equation 5], and the last equality
follows from (2.5).

Since $f$ is non-zero, to prove (5.6), we only need to observe that for non-
reproducing $K$, $\phi(0)-\int_{\mathbb{C}}\phi^{2}\,dA\not=0$. Suppose on the
contrary that $\phi(0)=\int_{\mathbb{C}}\phi^{2}\,dA$. Then

$K(z,z)=\int_{\mathbb{C}}K(z,w)K(w,z)\,dA(w),\qquad z\in{\mathbb{C}}.$

This means that the diagonals of the integral kernels of the operators $T_{K}$
and $T_{K}^{2}$ coincide. On the other hand, by (5.5), $T_{K}-T_{k}^{2}$ is a
positive operator. Thus $\mathrm{trace}(T_{K}-T_{K}^{2})=0$ and
$T_{K}=T_{K}^{2}$ — see, e.g., [23, Theorem 2.1] — which contradicts the
assumption that $K$ is not reproducing. ∎

## References

  * [1] L. D. Abreu, J. M. Pereira, J. L. Romero, and S. Torquato.  The Weyl-Heisenberg ensemble: hyperuniformity and higher Landau levels.  J. Stat. Mech. Theory Exp., (4):043103, 16, 2017.
  * [2] Y. Ameur, H. Hedenmalm, and N. Makarov.  Berezin transform in polynomial Bergman spaces.  Comm. Pure Appl. Math., 63(12):1533–1584, 2010.
  * [3] Y. Ameur, H. Hedenmalm, and N. Makarov.  Fluctuations of eigenvalues of random normal matrices.  Duke Math. J., 159(1):31–81, 2011.
  * [4] Y. Ameur, H. Hedenmalm, and N. Makarov.  Random normal matrices and Ward identities.  Ann. Probab., 43(3):1157–1201, 2015.
  * [5] R. J. Berman.  Determinantal point processes and fermions on polarized complex manifolds: bulk universality.  In Algebraic and analytic microlocal analysis, volume 269 of Springer Proc. Math. Stat., pages 341–393. Springer, Cham, 2018.
  * [6] A. Borodin.  Determinantal point processes.  In The Oxford handbook of random matrix theory, pages 231–249. Oxford Univ. Press, Oxford, 2011.
  * [7] O. Costin and J. L. Lebowitz.  Gaussian fluctuation in random matrices.  Phys. Rev. Lett., 75(1):69–72, 1995.
  * [8] T.-C. Dinh, S. Ghosh, H.-S. Tran, and M.-H. Tran.  Quantitative Marcinkiewicz’s theorem and central limit theorems: applications to spin systems and point processes.  arXiv preprint:2107.08469v2.
  * [9] M. Fenzl and G. Lambert.  Precise deviations for disk counting statistics of invariant determinantal processes.  Int. Math. Res. Not. IMRN, (10):7420–7494, 2022.
  * [10] G. B. Folland.  Harmonic analysis in phase space, volume 122 of Annals of Mathematics Studies.  Princeton University Press, Princeton, NJ, 1989.
  * [11] J. Ginibre.  Statistical ensembles of complex, quaternion, and real matrices.  J. Mathematical Phys., 6:440–449, 1965.
  * [12] A. Haimi.  Bulk asymptotics for polyanalytic correlation kernels.  J. Funct. Anal., 266(5):3083–3133, 2014.
  * [13] A. Haimi and H. Hedenmalm.  The polyanalytic Ginibre ensembles.  J. Stat. Phys., 153(1):10–47, 2013.
  * [14] A. Haimi and A. Wennman.  A central limit theorem for fluctuations in polyanalytic Ginibre ensembles.  Int. Math. Res. Not. IMRN, (5):1350–1372, 2019.
  * [15] J. B. Hough, M. Krishnapur, Y. Peres, and B. Virág.  Zeros of Gaussian analytic functions and determinantal point processes, volume 51 of University Lecture Series.  American Mathematical Society, Providence, RI, 2009.
  * [16] K. Johansson.  Random matrices and determinantal processes.  In Mathematical statistical physics, pages 1–55. Elsevier B. V., Amsterdam, 2006.
  * [17] O. Macchi.  The coincidence approach to stochastic point processes.  Advances in Appl. Probability, 7:83–122, 1975.
  * [18] T. Matsui, M. Katori, and T. Shirai.  Local number variances and hyperuniformity of the Heisenberg family of determinantal point processes.  J. Phys. A, 54(16):Paper No. 165201, 22, 2021.
  * [19] I. Nourdin and G. Peccati.  Normal approximations with Malliavin calculus, volume 192 of Cambridge Tracts in Mathematics.  Cambridge University Press, Cambridge, 2012.  From Stein’s method to universality.
  * [20] B. Rider and B. Virág.  Complex determinantal processes and $H^{1}$ noise.  Electron. J. Probab., 12:no. 45, 1238–1257, 2007.
  * [21] B. Rider and B. Virág.  The noise in the circular law and the Gaussian free field.  Int. Math. Res. Not. IMRN, (2):Art. ID rnm006, 33, 2007.
  * [22] T. Shirai.  Ginibre-type point processes and their asymptotic behavior.  J. Math. Soc. Japan, 67(2):763–787, 2015.
  * [23] B. Simon.  Trace ideals and their applications, volume 120 of Mathematical Surveys and Monographs.  American Mathematical Society, Providence, RI, second edition, 2005.
  * [24] A. Soshnikov.  Determinantal random point fields.  Uspekhi Mat. Nauk, 55(5(335)):107–160, 2000.
  * [25] A. Soshnikov.  Gaussian limit for determinantal random point fields.  Ann. Probab., 30(1):171–187, 2002.
  * [26] N. L. Vasilevski.  Poly-Fock spaces.  In Differential operators and related topics, Vol. I (Odessa, 1997), volume 117 of Oper. Theory Adv. Appl., pages 371–386. Birkhäuser, Basel, 2000.