Press "Enter" to skip to content

Libres pensées d'un mathématicien ordinaire Posts

Accéder aux articles plus facilement


L’extension Click & Read pour navigateur mise au point l’été dernier par l’INIST du CNRS réécrit les liens DOI contenus dans les pages web et ajoute une petite icône verte qui permet si possible d’accéder à la ressource, en donnant priorité aux versions ouvertes sur arXiv et HAL lorsqu’elles contiennent la meta-donnée DOI. Malheureusement, l’extension C&R ne permet pas de choisir entre la version éditeur et la version arXiv/HAL, ce qui fait qu’elle ne remplace pas vraiment Panist et Istex comme on pourrait s’y attendre, et en pratique ne dispense pas d’un copier-coller du DOI dans Sci-Hub (*). À ce propos, il existe également une extension Sci-Hub pour navigateur, dont le fonctionnement est un peu différent et assez complémentaire. Mais répétons-le, Sci-Hub est un service malheureusement illégal.

(*) Service illégal de mutualisation d’accès aux abonnements, très en avance sur son temps !

Leave a Comment

Dynamics of a planar Coulomb gas

This short post is about a joint work with François Bolley and Joaquín Fontbona on a Dyson Ornstein-Uhlenbeck process, published a few years ago. The numerics and graphics are new.

Above, two numerical experiments with $N=66$ particles and initial conditions on a uniform grid on a circle, and on a segment. Data is displayed as moving particles as well as trajectories. The one started from a segment is run ten times slower than the one started from a circle. These interacting particles on $\mathbb{R}^2=\mathbb{C}$ follow the Dyson Ornstein-Uhlenbeck dynamics
=\sqrt{\frac{2}{N}}\mathbb{d} B^{i,N}_t
-2X^{i,N}_t\,\mathbb{d} t
-\frac{\beta}{N}\sum_{j\neq i}\frac{X^{j,N}_t-X^{i,N}_t}{|X^{i,N}_t-X^{j,N}_t|^2}\mathbb{d}t,
\quad 1\leq i\leq N,
where $\beta\geq0$ is a parameter and $B$ is a standard Brownian motion on $(\mathbb{R}^2)^N=\mathbb{C}^N$.

The long time behavior of the particle system is
where $X^N_\infty$ follows the $\beta$ Ginibre planar Coulomb gas
\prod_{1\leq i < j\leq N}|x_i-x_j|^\beta \mathrm{d}x_1\cdots\mathrm{d}x_N, $$ while the high dimensional limit at equilibrium is the uniform law on a disc $$ \frac{1}{N}\sum_{i=1}^N\delta_{X^{i,N}_\infty} \overset{\mathrm{weak}}{\underset{N\to\infty}{\longrightarrow}} \frac{\mathbf{1}_{z\in\mathbb{C}:|z|\leq\sqrt{\frac{\beta}{2}}}}{\pi}. $$

Barycenter (first moment) of the cloud of particles. The stochastic process ${(m_t)}_{t\geq0}={(\frac{1}{N}\sum_{i=1}^nX^{i,N}_t)}_{t\geq0}$ is an Orsntein-Uhlenbeck process solving
where $W$ is a standard Brownian motion on $\mathbb{R}^2=\mathbb{C}$. In particular, for $t\geq0$,

This follows from the Itô formula and the fact that the real and imaginary parts of the first moment are eigenfunctions of the Dyson Ornstein-Uhlenbeck operator.

Second moment of the cloud of particles. The process ${(R_t)}_{t\geq0}={(\frac{1}{N}\sum_{i=1}^n|X^{i,N}_t|^2)}_{t\geq0}$ is a Cox-Ingersoll-Ross process ${(R_t)}_{t\geq0}$ solving the stochastic differential equation
\mathrm{d}R_t = \sqrt{\frac{8}{N^2}R_t}\mathrm{d}W_t
+4\left[\frac{1}{N}+ \frac{\beta}{4}\frac{N-1}{N}-R_t \right]\mathrm{d}t,
where ${(W_t)}_{t\geq0}$ is a real standard Brownian motion. In particular, its invariant distribution is the Gamma law $\Gamma_N=\mathrm{Gamma}(a,\lambda)$ with shape parameter $a=N+\frac{\beta}{4}N(N-1)$ and scale parameter $\lambda=N^2$. Moreover, for any initial condition $R_0$ and $t \geq 0$,
&\leq \mathrm{e}^{-4t}\mathrm{Wasserstein}_1(\mathrm{Law}(R_0),\Gamma_N)\\
\mathbb{E}[R_t\mid R_0]
&= R_0\mathrm{e}^{-4t}+\Bigr(\frac{1}{N}+\frac{\beta}{4}\frac{N-1}{N}\Bigr)(1-\mathrm{e}^{-4t}).
Note that $\frac{\beta}{4}$ is the second moment of the uniform distribution on the disc of radius $\sqrt{\frac{\beta}{2}}$.

This follows from Itô’s formula, Lévy’s characterization of BM, and the fact that the second moment is, up to a constant, an eigenfunction of the Dyson Ornstein-Uhlenbeck operator.

Open problem. Compute the spectral gap of the dynamics, or find how it depends over $N$.

Julia code.

using Plots

# Define the 2D Dyson Ornstein-Uhlenbeck process
Base.@kwdef mutable struct DOU2D
    n::Int = 2         # number of particles
    β::Float64 = 2.    # repulsion parameter  
    T::Float64 = 10.   # terminal time
    dt::Float64 = 1E-3 # time increment
    r::Float64  = 1.   # initial conditions uniform on segment [-r,r]+3*im
    m = floor(Int, T / dt) # number of times
    x::Array{Complex{Float64},2} = zeros(Complex{Float64},m,n)
end # end struct

function compute!(X::DOU2D)
    # initial condition on grid on centered circle of radius r
    #Theta = range(1, X.n, length = X.n) * 2 * pi / X.n
    #X.x[1,:] = X.r * (cos.(Theta) + im * sin.(Theta)) # unif grid on circle
    # initial condition on grid on segment [-r,r] - 3*im
    X.x[1,:] = range(-X.r, X.r, length = X.n) .- 3 * im
    # trajectories
    for i = 2:X.m
     dB = (randn(1,X.n) + im * randn(1,X.n))/sqrt(2)
     X.x[i,:] = copy(X.x[i-1,:])
     for j in 1:X.n # particles
         X.x[i,j] += sqrt(2/X.n) * dB[j] * sqrt(X.dt)
         X.x[i,j] += -2 * X.x[i-1,j] * X.dt
         for k in 1:X.n
             if (j == k) continue end
             X.x[i,j] += X.β * X.dt / (X.n * (conj(X.x[i-1,j] - X.x[i-1,k])))
         end # for 
     end # for
    end # for
end # function

# process and graphics
dou = DOU2D(n = 66, r = 3., T = 5., dt = 1E-4)
pdou2d = plot()
for j in 1:dou.n # particles
    plot!(pdou2d, real(dou.x[:,j]), imag(dou.x[:,j]), aspect_ratio =:equal, legend = false)
end # for
r = dou.r + 1
anim = @animate for i = 1:100:dou.m
    scatter(real(dou.x[i,:]), imag(dou.x[i,:]),
            xlims = (-r,r), ylims = (-r,r),
            aspect_ratio =:equal, legend = false)
gif(anim, "dou2danim.gif", fps = 10)

Further reading.

Leave a Comment

Science ouverte (mathématiques)


Quelques bonnes nouvelles sur le front de la science ouverte du côté des sociétés savantes :

Leave a Comment

Spectral radius from characteristic polynomial

Phase portrait of reciprocal of characteristic polynomial of a high dimensional Ginibre matrix
Phase portrait of reciprocal characteristic polynomial of a high dimensional Ginibre matrix – Julia code given below. The interior and exterior of the unit disc of the complex plane are swapped by the reciprocal transform. The spectrum lies essentially outside.

This post is about the work arXiv:2012.05602 in collaboration with Charles Bordenave and David García-Zelada on the spectral radius. Clearly this is my favorite work so far!

Consider a square random matrix with independent and identically distributed entries of mean zero and unit variance. We show that as the dimension tends to infinity, the spectral radius is equivalent to the square root of the dimension in probability. This result can also be seen as the convergence of the support in the circular law theorem under optimal moment conditions. In the proof we establish the convergence in law of the reciprocal characteristic polynomial to a random analytic function outside the unit disc, related to a hyperbolic Gaussian analytic function. The proof is short and differs from the usual approaches for the spectral radius. It relies on a tightness argument and a joint central limit phenomenon for traces of fixed powers.

Model. Let $\{a_{jk}\}_{j,k\geq1}$ be independent and identically distributed complex random variables with mean zero and unit variance, namely $$\mathbb{E}[a_{11}]=0\quad\text{and}\quad\mathbb E[|a_{11}|^2] = 1.$$ For all $n\geq1$, let us consider the Girko random matrix $$A_n={(a_{jk})}_{1 \leq j,k \leq n}.$$ When $a_{11}$ is Gaussian with independent and identically distributed real and imaginary parts then $A_n$ belongs to the complex Ginibre ensemble. We are interested in the matrix $$\frac{1}{\sqrt{n}}A_n.$$By the law of large numbers, almost surely, as $n\to\infty$, the rows (and the columns) have unit Euclidean norm and are orthogonal. Its characteristic polynomial at point $z\in\mathbb{C}$ is $$p_n(z)=\det\Bigr(z-\frac{A_n}{\sqrt{n}}\Bigr)$$ where $z$ stands for $z$ times identity matrix. The $n$ roots of $p_n$ in $\mathbb{C}$ are the eigenvalues of $\frac{1}{\sqrt{n}}A_n$. They form a multiset $\Lambda_n$ which is the spectrum of $\frac{1}{\sqrt{n}}A_n$. The spectral radius is $$\rho_n=\max_{\lambda\in\Lambda_n}|\lambda|.$$ Following Ginibre, Girko, Bai, …, and finally Tao and Vu, the circular law (of nature) states that the empirical measure of the elements of $\Lambda_n$ tends weakly as $n\to\infty$ to the uniform distribution on the closed unit disc: almost surely, for every nice Borel set $B\subset\mathbb{C}$, $$\lim_{n\to\infty}\frac{\mathrm{card}(B\cap\Lambda_n)}{n} =\frac{\mathrm{area}(B\cap\overline{\mathbb{D}})}{\pi},$$ where “$\mathrm{area}$” stands for the Lebesgue measure on $\mathbb C$, and where $\overline{\mathbb{D}}=\{z\in\mathbb{C}:|z|\leq1\}$ is the closed unit disc. The circular law, which involves weak convergence, does not provide the convergence of the spectral radius, it gives only $$\varliminf_{n\to\infty}\rho_n\geq1.$$

Our result on the spectral radius.  We have $\rho_n\overset{\mathbb{P}}{\underset{n\to\infty}{\longrightarrow}}1$, in the sense that for all $\varepsilon>0$, $$\lim_{n\to\infty}\mathbb{P}(|\rho_n-1|\geq\varepsilon)=0.$$ The moments assumptions of zero mean and unit variance are optimal. Moreover the $\frac{1}{\sqrt{n}}$ scaling is no longer adequate for entries of infinite variance.

We have $\rho_n\leq\sigma_n$ where $\sigma_n$ is the operator norm of $\frac{1}{\sqrt{n}}A_n$ (largest singular value). Bai and al have proved that the condition $\mathbb{E}[|a_{11}|^4]<\infty$ is necessary and sufficient for the convergence of $\sigma_n$ as $n\to\infty$ and that the limit is then $2$. A striking aspect of the spectral radius, in comparison with the operator norm, is that it converges without such an extra condition.

Proof. It does not involve any Hermitization or norms of powers in the spirit of Gelfand’s spectral radius formula. The idea is to show that on $\mathbb{C}\cup\{\infty\}\setminus\overline{\mathbb{D}}$, the function $$z\mapsto z^{-n}p_n(z)$$ tends as $n\to\infty$ to a random analytic function which does not vanish. The first step for mathematical convenience is to convert $\mathbb{C}\cup\{\infty\}\setminus\overline{\mathbb{D}}$ into $\mathbb D=\{z\in\mathbb{C}:|z|<1\}$ by noting that $p_n(z)=z^nq_n(1/z)$, $z\not\in\overline{\mathbb{D}}$, where for all $z\in\mathbb{D}$,
q_n(z) = \det\left(1- z\frac{A_n}{\sqrt n}\right)
\]is the reciprocal polynomial of the characteristic polynomial $p_n$. Let $\mathrm{H}(\mathbb{D})$ be the set of holomorphic or complex analytic functions on $\mathbb{D}$, equipped with the topology of uniform convergence on compact subsets, the compact-open topology, studied by Shirai. This allows to see $q_n$ as a random variable on $\mathrm{H}(\mathbb{D})$ and gives a meaning to convergence in law of $q_n$ as $n\to\infty$, namely, $q_n$ converges in law to some random element $q$ of $\mathrm{H}(\mathbb{D})$ if for every bounded real continuous function $f$ on $\mathrm{H}(\mathbb{D})$, $$\mathbb E[f(q_n)] \to \mathbb E[f(q)].$$ Namely, we prove that $$q_n
\xrightarrow[n \to \infty]{\mathrm{law}}\kappa \mathrm{e}^{-F},$$ where $F$ is the random holomorphic function on $\mathbb D$ defined by $$F(z)=\sum_{k=1}^\infty X_k \frac{z^k}{\sqrt k}$$ where $\{X_k\}_{k\geq1}$ are independent complex Gaussian random variables such that
$$\mathbb E\big[X_k\big]=0,\quad\mathbb E\left[|X_k|^2\right] = 1\quad\text{and}\quad \mathbb E\left[X_k^2 \right] = \mathbb E \left[a_{11}^2 \right]^k,$$ and where $\kappa: \mathbb D \to \mathbb C$ is the holomorphic function defined for all $z\in\mathbb{D}$ by $$\kappa(z) = \sqrt{1-z^2 \mathbb E \left[a_{11}^2 \right]}.$$ The square root defining $\kappa$ is the one such that $\kappa(0)=1$. Notice that it is a well-defined holomorphic function on the simply connected domain $\mathbb{D}$ since the function $z\mapsto1-z^2\mathbb{E}[a_{11}^2]$ does not vanish on $\mathbb{D}$, indeed $|\mathbb{E}[a_{11}^2]|\leq\mathbb{E}[|a_{11}|^2]=1$.

The proof of the convergence of $q_n$ is partially inspired by a work due to Basak and Zeitouni and relies crucially on a joint combinatorial central limit theorem for traces of fixed powers inspired from a work by Janson and Nowicki. Unlike previous arguments used in the literature for the analysis of Girko matrices, the approach does not rely on Girko Hermitization, Gelfand spectral radius formula, high order traces, resolvent method or Cauchy-Stieltjes transform. The first step consists in showing the tightness of ${(q_n)}_{n\geq1}$, by using a decomposition of the determinant into orthogonal elements related to determinants of submatrices, as in the work of Basak and Zeitouni:$$\det\Bigr(1-z\frac{A_n}{\sqrt{n}}\Bigl)=1+\sum_{k=1}^n(-z)^n\sum_{\substack{I\subset\{1,\ldots,n\}\\|I|=k}}\frac{\det((A_n)_{I,I})}{\sqrt{n}^k}.$$Knowing this tightness, the problem is reduced to show the convergence in law of these elements. A reduction step, inspired by the work of Janson and Nowicki, consists in truncating the entries, reducing the analysis to the case of bounded entries. The next step consists in a central limit theorem for product of traces of powers of fixed order$$\det\Bigr(1-z\frac{A_n}{\sqrt{n}}\Bigr)=\exp\Bigr(-\sum_{k=1}^\infty\frac{\mathrm{Trace}(A_n^k)}{\sqrt{n}^k}\frac{z^k}{k}\Bigr).$$It is important to note that we truncate with a fixed threshold with respect to $n$, and the order of the powers in the traces are fixed with respect to $n$. This is in sharp contrast with the usual Füredi-Komlós truncation-trace approach related to the Gelfand spectral radius formula used in all the previous approaches for the spectral radius.

Moment assumptions. The universality for the first order global asymptotics stated by the circular law depends only on the trace $\mathbb{E}[|a_{11}|^2]$ of the covariance matrix of $\Re a_{11}$ and $\Im a_{11}$. The universality stated by the convergence of $q_n$, just like for the central limit theorem, depends on the whole covariance matrix. Since
$$\mathbb{E}[a_{11}^2]=\mathbb{E}[(\Re a_{11})^2]-\mathbb{E}[(\Im a_{11})^2]+2\mathrm{i}\mathbb{E}[\Re a_{11}\Im a_{11}],$$ we can see that $\mathbb{E}[a_{11}^2]=0$ if and only if $$\mathbb{E}[(\Re a_{11})^2]=\mathbb{E}[(\Im a_{11})^2]\quad\text{and}\quad\mathbb{E}[\Re a_{11}\Im a_{11}]=0.$$  Moreover, we cannot in general get rid of $\mathbb{E}[a_{11}^2]$ by simply multiplying $A_n$ by a phase.

Hyperbolic Gaussian analytic function. When $\mathbb E\left[a_{11}^2\right]=0$ then $\kappa=1$ while the random analytic function $F$ which appears in the limit is a degenerate case of the well-known hyperbolic Gaussian Analytic Functions (GAFs). It can also be obtained as the antiderivative of the $L=2$ hyperbolic GAF which is $0$ at $z=0$. This $L=2$ hyperbolic GAF is related to the Bergman kernel and could be called the Bergman GAF. These GAFs appear also at various places in mathematics and physics and, in particular, in the asymptotic analysis of Haar unitary matrices.

Cauchy-Stieltjes transform. If $\mathbb{E}[a_{11}^2]=0$ then by returning to $p_n$, taking the logarithm and the derivative with respect to $z$ in the convergence of $q_n$, we obtain the convergence in law of the Cauchy-Stieltjes transform (complex conjugate of the electric field) minus $n/z$ towards $z \mapsto F'(1/z)/z^2$ which is a Gaussian analytic function on $\mathbb C \setminus \overline{\mathbb D}$ with covariance given by a Bergman kernel.

Central Limit Theorem. The convergence of $q_n$ can be seen as a central limit theorem for the log-determinant (global second order analysis). Namely for all $z\in\mathbb{D}$, we have $$|q_n(1/z)| = \exp\left[-n \left(U_n(z)-U(z)\right)\right]$$ where $$U_n(z)=-\frac{1}{n}\log|p_n(z)|\quad\text{and}\quad U(z)=-\log|z|$$ are the logarithmic potential at the point $z$ of the empirical spectral distribution of $\frac{1}{\sqrt n} A_n$ and of the uniform distribution on the unit disc $\mathbb{D}$.

Moreover, it is possible to extract from the convergence of $q_n$ a CLT for linear spectral statistics with respect to analytic functions in a neighborhood of $\overline{\mathbb{D}}$. This can be done by using the Cauchy formula for an analytic function $f$,
\begin{align*}\int f(\lambda) \mathrm \mu(\mathrm d \lambda)
&= \frac{1}{2\pi\mathrm{i}}
\int \left(\oint \frac{f(z)}{z – \lambda} \mathrm d z\right) \mu(\mathrm d
&= \frac{1}{2\pi\mathrm{i}} \oint f(z) \left(\int \frac{\mu(\mathrm d\lambda)} {z – \lambda} \right)\mathrm d z\\
&= \frac{1}{2\pi\mathrm{i}} \oint f(z)(\log \det \left(z – A \right))’ \mathrm d z
\end{align*}\]where $\mu$ is the counting measure of the eigenvalues of $A$, where the contour integral is around a centered circle of radius strictly larger than $1$, and where we have taken any branch of the logarithm. The approach is purely complex analytic. In particular, it is different from the usual approach with the logarithmic potential of $\mu$ based on the real function given by $z \mapsto \int\log|z-\lambda|\mu(\mathrm{d}\lambda) = \log|\det(z-A)|$.

Wigner case and elliptic interpolation. The finite second moment assumption is optimal. We could explore its relation with the finite fourth moment assumption for the convergence of the spectral edge of Wigner random matrices, which is also optimal. Heuristic arguments tell us that the interpolating condition on the matrix entries is $$\mathbb{E}[|a_{jk}|^2|a_{kj}|^2]<\infty,\quad j\neq k,$$ which is a finite second moment condition for Girko matrices and a finite fourth moment condition for Wigner matrices. This is work in progress.

Coupling and almost sure convergence. For simplicity, we define our random matrix $A_n$ for all $n\geq1$ by truncating from the upper left corner the infinite random matrix $\{a_{jk}\}_{j,k\geq1}$. This imposes a coupling for the matrices $\{A_n\}_{n\geq1}$. However, since our result of the spectral radius involves a convergence in probability, it remains valid for an arbitrary coupling, in the spirit of the triangular arrays assumptions used for classical CLTs. In another direction, one could ask about the upgrade to almost sure convergence. This is an open problem.

Heavy tails. An analogue of the circular law in the heavy-tailed case $\mathbb{E}[|a_{11}|^2]=\infty$ is already available but requires another scaling than $\frac{1}{\sqrt{n}}$. The spectral radius of this model tends to infinity as $n\to\infty$ but it could be possible to analyze the limiting point process at the edge as $n\to\infty$ and its universality. This is an open problem.

Julia code. Here is a Julia code illustrating the high dimensional phenomenon.

using LinearAlgebra, Plots # for eigvals(), scatter(), heatmap()

function phaseportrait(n)
 M = (randn(n,n)+im*randn(n,n))/sqrt(2*n)
 Spec = eigvals(M)
 r = 1000
 c = max.(abs.(Spec))
 c = c + c/2
 M = zeros(r,r)
 for i in 1:r
  for j in 1:r
   z = (-c + 2*c*i/r) + im * (-c + 2*c*j/r)
   M[i,j] = angle(prod(z.*Spec.-1)) # reciprocal charpoly phase
end # function

Discrete analogue. The method can be adapted to sparse Boolean matrices, replacing the Gaussian limiting regime by a Poisson limiting regime, in relation to important aspects of the high dimensional analysis of random graphs. This was done by Simon Coste.

Further reading.

Phase portrait of reciprocal characteristic polynomial of a high dimensional real Gaussian matrix
Phase portrait of reciprocal characteristic polynomial of a high dimensional real (hence the spectrum symmetry) Gaussian matrix.
Leave a Comment
Syntax · Style · .