# Generating function of $SO(N)$ random matrix

I am interested in the generating function of $SO(N)$ random matrix, that is, I want to compute $$Z_N[J]=\int dM e^{{\rm Tr} (J^T M)},$$ where $dM$ is the $SO(N)$ Haar measure, and $J$ is an arbitrary $N\times N$ matrix. From this generating function, I can generate all correlations $\langle M_{ij}M_{kl}\cdots\rangle$ by taking derivatives with respect to the elements of $J$.

Due to the invariance of the measure, one sees that $Z[J]=Z[U^TJV]$ with $U,V\in SO(N)$, and thus $Z$ only depends on the singular values of $J$. (Stated otherwise, $Z$ only depends on the $N$ invariants ${\rm Tr}((J^TJ)^n)$, $n=1,...,N$).

Finally, at least for $N=2$ and $N=3$, one can show that $Z$ is also invariant under permutations of the singular values of $J$ (maybe it can be generalized for all $N$ ?).

It is not too hard to compute explicitly $Z_2[J]$, which is given in terms of a Bessel function of the sum of the two singular values of $J$.

Is there a way (or has it been done in the literature) to compute $Z_N$ for any $N$ ? I would already be happy with $Z_3$, which I cannot manage to compute explicitly.

EDIT : Here is an attempt, which kind of works for $N=2$, but for which I am stuck for $N=3$. If we define a Laplacian $\Delta=\sum_{ij}\frac{\partial^2}{\partial J_{ij}^2}$ (with $J_{ij}$ the elements of $J$), one shows easily that $$\Delta Z[J]=N Z[J]. \tag{1}$$ If we call $\lambda_i$ the singular values of $J$ (with $\lambda_1>\lambda_2>\ldots$), using the fact that $Z[J]=Z[\lambda_1,\lambda_2,\ldots]$, one shows (at least for $N=2$ and $N=3$, but it might be generalizable to $N\geq4$) that $$\Delta Z=\frac{1}{D}\sum_{i}\frac{\partial}{\partial \lambda_i}\left(D\frac{\partial}{\partial \lambda_i}Z\right),$$ where $D=\prod_{i< j}(\lambda_i^2-\lambda_j^2)$ is related to the Jacobian to go from $J_{ij}$ to $\lambda_i$. This equation looks nice enough, so my hope is that a solution exists, I am not quite sure how to find it for $N=3$.

In the case $N=2$, we can compute $Z_2$ exactly via its definition, and it reads $Z_2[\lambda_1,\lambda_2]=I_0(\lambda_1+\lambda_2)$, with $I_\nu$ the modified Bessel function of the first kind. One checks that this is indeed a solution of Eq. (1).

Unfortunately, even in that case, it is not clear to me how to find this solution starting from Eq. (1) only. Given that $D=\lambda_1^2-\lambda_2^2$, it is tempting to defined $u=(\lambda_1+\lambda_2)/2$ and $v=(\lambda_1-\lambda_2)/2$. Then Eq. (1) is solved by separation of variables and we find a family of solution $Z_{2,\mu}$ (where I have already use the fact that $Z_N[0]=1$) : $$Z_{2,\mu}[u,v]=I_0(\sqrt{\mu}u)I_0(\sqrt{4-\mu}v).$$ Clearly the solution to my problem corresponds to $\mu=4$, but it is not clear to me what is the rigorous argument to pick this value of $\mu$ (since $Z_N>0$ $\forall J$, we must have $\mu\geq 4$ as $I_0$ can be negative for imaginary variables; but how to select $\mu=4$ as the only viable solution ?).

One way to solve this issue of $\mu$ is to use the fact that $Z_2[J=\lambda Id_2]$ can be computed explicitly: $Z_2[J=\lambda Id_2]=I_0(2\lambda)$, which unambiguously selects $\mu=4$.

Since we can always compute the expansion of $Z_N[\lambda Id_N]$ explicitly at least for small $\lambda$, this kind of argument might be enough to fix the constant also for $N>2$. For $N=3$, a few special cases can be computed explicitly, which might help too.

Any insight for the solution of Eq. (1) would be greatly appreciated.

• If you expand $e^{{\rm Tr}(J^TM)}=\sum_n \frac{1}{n!}[{\rm Tr}(J^TM)]^n$, you can use the results discussed in this question Jun 7, 2018 at 11:29
• @Marcel : Thank you for the reference. In addition of being quite hard to decipher for a physicist like me, I was hopping to get an explicit solution (for $N=3$, say). This solution exists at least for $N=2$: $Z_2=I_0\left(\sqrt{{\rm Tr}(J J^T)+2\det J}\right)$, with $I_\nu$ the modified Bessel function of the first kind. It is quite unclear to me how to recover this results from the answers to the question you linked to.
Jun 7, 2018 at 12:04

Singular value decomposition: The real $N\times N$ matrix $J$ has singular values $\sigma_1,\sigma_2\ldots\sigma_N\geq 0$ in the decomposition $J=U\,{\rm diag}\,(\sigma_1,\sigma_2\ldots\sigma_N)V$ with $U,V\in{\rm O}(N)$. We need to restrict $U,V$ to ${\rm SO}(N)$, which we can do by allowing for one of the singular values to become negative: $$J=U\,{\rm diag}\,(\lambda_1,\lambda_2\ldots\lambda_N)V,\;\;U,V\in{\rm SO}(N),$$ $$\lambda_n=\sigma_n\geq 0,\;\;{\rm for}\;\;n=1,2,\ldots N-1,$$ $$\lambda_N=\begin{cases} \sigma_N\;{\rm if}\;{\rm det}\,M\geq 0,\\ -\sigma_N\;{\rm if}\;{\rm det}\,M<0. \end{cases}$$ Note that this need to consider the $\lambda_n$'s implies that the integral $Z_N$ is not only a function of the invariants ${\rm Tr}((J^TJ)^n)$.

Without loss of generality we may now take a diagonal $J={\rm diag}\,(\lambda_1,\lambda_2,\ldots \lambda_N)$ and then we need the integral $$Z_N=\int_{{\rm O}(N)} dM \exp\left(\sum_{n=1}^N \lambda_n M_{nn}\right).$$
Large-$N$ limit: I do not think there is a closed form expression for arbitrary $N$, but for $N\gg 1$ the diagonal matrix elements of $M$ are independent Gaussians (mean zero, variance $1/N$), resulting in $$Z_N\rightarrow \exp\left(\frac{1}{2N}\sum_{n=1}^N \lambda_n^2\right),\;\;N\gg 1.$$ This result for $Z_N$ in the large-$N$ limit depends only on the singular values (because $\lambda_n$ enters only quadratically), but for small-$N$ this is no longer true.

Small-$N$ results: For small $N$ we can use the parametrization of the orthogonal group given in Appendix B of arXiv:1405.3115.

• $N=2$ $$Z_2=\frac{1}{\pi}\int_0^{\pi} \exp\bigl[(\lambda_1+\lambda_2)\cos\theta\bigr]\,d\theta=I_0(\lambda_1+\lambda_2)$$ Note that in terms of the singular values considered in the OP, this would read $Z_2=I_0(\sigma_1\pm\sigma_2)$, where the $\pm$ sign is the sign of ${\rm det}\,J$.

• $N=3$ $$Z_3=\frac{1}{8\pi^2}\int_0^{2\pi} d\alpha\int_0^{2\pi} d\alpha'\int_0^\pi \sin\theta\, d\theta$$ $$\qquad\qquad\times\exp\left[\cos \alpha \cos \alpha' (\lambda_1+\lambda_2 \cos \theta)-\sin \alpha \sin \alpha' (\lambda_1 \cos \theta+\lambda_2)+\lambda_3\cos \theta\right]$$ The two integrals over $\alpha,\alpha'$ can be carried out (by rewriting them as integrals over $\alpha\pm\alpha'$), with the result $$Z_3=\frac{1}{2}\int_{-1}^1 I_0\left[\tfrac{1}{2}(\lambda_1-\lambda_2) (1-x)\right] I_0\left[\tfrac{1}{2} (\lambda_1+\lambda_2)(1+x)\right]\,e^{\lambda_3 x}\,dx.$$ As shown by Paul Enta, this integral can be written as the inverse Laplace transform ${\cal L}^{-1}[F(s)](t)$ for $t=1$ of the function $$F_3(s)=\frac{1}{\sqrt{(s-\lambda_1-\lambda_2-\lambda_3) (s+\lambda_1+\lambda_2-\lambda_3) (s+\lambda_1-\lambda_2+\lambda_3) (s-\lambda_1+\lambda_2+\lambda_3)}},$$ to demonstrate that it is invariant under permutation of the $\lambda_n$'s.
Notice also that $F_3(s)$ is invariant under a sign change of two of the $\lambda_n$'s, but not under a sign change of one single $\lambda_n$. This is consistent with the definition of the $\lambda_n$'s in terms of the singular values, given above.

A closed-form expression of either the Bessel-function integral or the inverse Laplace transform seems not forthcoming, except for some special cases (see the comments by Adam).

Average over the full orthogonal group: If we integrate over the full group ${\rm O}(N)$, instead of only over ${\rm SO}(N)$, we should combine the results for $\pm\lambda_N$, $$Z_{{\rm O}(N)}(\lambda_1,\lambda_2,\ldots \lambda_{N-1},\lambda_N)= \tfrac{1}{2}Z_{{\rm SO}\,(N)}(\lambda_1,\lambda_2,\ldots \lambda_{N-1},\lambda_N)+ \tfrac{1}{2}Z_{{\rm SO}(N)}(\lambda_1,\lambda_2,\ldots \lambda_{N-1},-\lambda_N).$$ There is now no need to distinguish the $\lambda_n$'s from the singular values $\sigma_n$. In particular, for $N=2$ one has $$Z_{{\rm O}(2)}=\tfrac{1}{2}I_0(\sigma_1+\sigma_2)+\tfrac{1}{2}I_0(\sigma_1-\sigma_2).$$

• There is at least a closed form for $N=2$, see the edit part of my question, and my hope is that it can be compute for some specific $N$ (say $N=3$). Thanks for the large $N$ limit, though.
Jun 7, 2018 at 14:07
• Sure, it is not hard to find an integral form of $Z_3$, but I'm looking for a more explicit form. One of the main issue with this representation is that the symmetry under permutation of $\lambda_i$ is not explicit (this would be needed for further calculation), see also math.stackexchange.com/questions/2808873/…
• Your integral for $Z_3$ seems to give $Z_3(0.3,0.4,0.5)\neq Z_3(0.4,0.3,0.5)$, but it seems clear that the true $Z_3$ should be symmetric in $\lambda$. Jun 7, 2018 at 16:30
• @NeilStrickland --- I corrected a typo (the integrals over $\alpha$ and $\alpha'$ should run from $0$ to $2\pi$, not from $0$ to $\pi$), and now it is symmetric; thanks for catching this. Jun 7, 2018 at 20:47
• Using this representation, one can also compute the case $\lambda_3=0$ and $\lambda_1=\lambda_2=\lambda$, which reads $I_0(2\lambda)+\frac\pi2 L_1(2\lambda) I_0(2\lambda)-\frac\pi2 L_0(2\lambda) I_1(2\lambda)$ with $L_\nu(x)$ the Struve L function.
For the $n=3$ case, it is best to use quaternionic methods to express $SO(3)$ as a quotient of $S^3$. Explicitly, for $a\in S^3$ one can check that the matrix $$A = \left[\begin{matrix} a_1^2-a_2^2-a_3^2+a_4^2 & 2a_1a_2-2a_3a_4 & 2a_1a_3+2a_2a_4 \\ 2a_1a_2+2a_3a_4 & -a_1^2+a_2^2-a_3^2+a_4^2 & -2a_1a_4+2a_2a_3 \\ 2a_1a_3-2a_2a_4 & 2a_1a_4+2a_2a_3 & -a_1^2-a_2^2+a_3^2+a_4^2 \end{matrix}\right]$$ lies in $SO(3)$. (This is based on an identification $\mathbb{H}\simeq\mathbb{R^4}$ with $1\in\mathbb{H}$ corresponding to $(0,0,0,1)\in\mathbb{R}^4$.) This construction comes from a group homomorphism, so the Haar measure on $SO(3)$ corresponds to the natural round measure on $S^3$, suitably rescaled. Thus, we just want to integrate the function $\exp(f)$ over $S^3$, where $$f(a) = (a_1^2-a_2^2-a_3^2+a_4^2)\lambda_1 + (-a_1^2+a_2^2-a_3^2+a_4^2)\lambda_2 + (-a_1^2-a_2^2+a_3^2+a_4^2)\lambda_3.$$ Using $\sum_ia_i^2=1$, we get $f(a)=\sum_i\lambda_i-2g(a)$, where $$g(a) = (a_2^2+a_3^2)\lambda_1 + (a_1^2+a_3^2)\lambda_2 + (a_1^2+a_2^2)\lambda_3.$$ This is at least visibly invariant under permutations of $\lambda$. Next, because $a_4$ does not appear in $g(a)$, it is natural to use a stereographic parametrisation $\sigma\colon\mathbb{R}^3\to S^3$, as follows: $$\sigma(x_1,x_2,x_3) = \frac{(2x_1,2x_2,2x_3,\|x\|^2-1)}{\|x\|^2+1}$$ The Jacobian determinant of $\sigma$ turns out to be $8/(\|x\|^2+1)^3$. This leads to an expression for $Z_3$ as an integral over $\mathbb{R}^3$ that retains all the natural symmetries, but it does not seem easy to evaluate either numerically or symbolically.
• Yes, written in this way, the integral can be seen to be invariant under permutation. But unfortunately, it does not help to write $Z$ as a function explicitely invariant (i.e. without using a change of variable in the integral). I was hoping that maybe a clever parametrization of SO(3) might do the job (see also this question : math.stackexchange.com/questions/2808873/…), but I haven't got any luck about that yet.