Rank probabilities for real random $N\times N\times 2$ tensors

We prove that the probability $P_N$ for a real random Gaussian $N\times N\times 2$ tensor to be of real rank $N$ is $P_N=(\Gamma((N+1)/2))^N/G(N+1)$, where $\Gamma(x)$, $G(x)$ denote the gamma and Barnes $G$-functions respectively. This is a rational number for $N$ odd and a rational number multiplied by $\pi^{N/2}$ for $N$ even. The probability to be of rank $N+1$ is $1-P_N$. The proof makes use of recent results on the probability of having $k$ real generalized eigenvalues for real random Gaussian $N\times N$ matrices. We also prove that $\log P_N= (N^2/4)\log (e/4)+(\log N-1)/12-\zeta '(-1)+{\rm O}(1/N)$ for large $N$, where $\zeta$ is the Riemann zeta function.


Introduction
The (real) rank of a real m × n × p 3-tensor or 3-way array T is the well defined minimal possible value of r in an expansion where ⊗ denotes the tensor (or outer) product [1,3,4,8].
If the elements of T are choosen randomly according to a continuous probability distribution, there is in general (for general m, n and p) no generic rank, i.e., a rank which occurs with probability 1. Ranks which occur with strictly positive probabilities are called typical ranks. We assume that all elements are independent and from a standard normal (Gaussian) distribution (mean 0, variance 1). Until now, the only analytically known probabilities for typical ranks were for 2 × 2 × 2 and 3 × 3 × 2 tensors [2,7]. Thus in the 2 × 2 × 2 case the probability that r = 2 is π/4 and the probability that r = 3 is 1 − π/4, while in the 3 × 3 × 2 case the probability of the rank equaling 3 is the same as the probability of it equaling 4 which is 1/2. Before these analytic results the first numerical simulations were performed by Kruskal in 1989, for 2 × 2 × 2 tensors [8], and the approximate values 0.79 and 0.21 obtained for the probability of ranks r = 2 and r = 3 respectively. For N × N × 2 tensors ten Berge and Kiers [10] have shown that the only typical ranks are N and N + 1. From ten Berge [9], it follows that the probability P N for an N × N × 2 tensor to be of rank N is equal to the probability that a pair of real random Gaussian N × N matrices T 1 and T 2 (the two slices of T ) has N real generalized eigenvalues, i.e., the probability that det(T 1 − λT 2 ) = 0 has only real solutions λ [2,9]. Knowledge about the expected number of real solutions to det(T 1 − λT 2 ) = 0 obtained by Edelman et al. [5] led to the analytical results for N = 2 and N = 3 in [2]. Forrester and Mays [7] have recently determined the probabilities p N,k that det(T 1 − λT 2 ) = 0 has k real solutions, and we here apply the results to P N = p N,N to obtain explicit expressions for the probabilities for all typical ranks of N ×N ×2 tensors for arbitrary N , hence settling this open problem for tensor decompositions. We also determine the precise asymptotic decay of P N for large N and give some recursion formulas for P N .
2 Probabilities for typical ranks of N × N × 2 tensors As above, assume that T 1 and T 2 are real random Gaussian N × N matrices and let p N,k be the probability that det(T 1 − λT 2 ) = 0 has k real solutions. Then Forrester and Mays [7] prove: where the asterisk indicates that the sum is over k values of the same parity as N . For N even we have while for N odd and The expressions for β l and β l+1/2 are given in [7], but are not needed here, and ⌈·⌉ denotes the ceiling function.
The method used in [7] relies on first obtaining the explicit form of the element probability density function for A real Schur decomposition is used to introduce k real and (N − k)/2 complex eigenvalues, with the imaginary part of the latter required to be positive (the remaining (N − k)/2 eigenvalues are the complex conjugate of these), for k = 0, 2, . . . , N (N even) and k = 1, 3, . . . , N (N odd). The variables not depending on the eigenvalues can be integrated out to give the eigenvalue probability density function, in the event that there are k real eigenvalues. And integrating this over all allowed values of the real and positive imaginary part complex eigenvalues gives P N,k . From Theorem 1 we derive our main result: Theorem 2. Let P N denote the probability that a real N × N × 2 tensor whose elements are independent and normally distributed with mean 0 and variance 1 has rank N . We have where is the Barnes G-function and Γ(x) denotes the gamma function. More explicitly P 2 = π/4, and for N ≥ 4 even while for N odd Hence P N for N odd is a rational number but for N even it is a rational number multiplied by π N/2 . The probability for rank N + 1 is 1 − P N .
Proof. From [2] we know that P N = p N,N . Hence, by Theorem 1 Since the values of β l and β l+1/2 are not needed for the determination of P N . By (3) we immediately find if N is even. For N odd we use (4) to get Substituting the expressions for α l and α l+1/2 into these formulas we obtain, after simplifying, for N even and for N odd where to obtain the second equality use has been made of the duplication formula for the gamma function, and to obtain the third equality the expression (9) for the Barnes G-function has been used. Furthermore, for each N even where to obtain the final equation use is made of the fundamental gamma function recurrence and for N odd , N = 3, 7, 11, . . .

Recursion formulas and asymptotic decay
By Theorem 2 it is straightforward to calculate P N +1 /P N from either (8) or (10) and (11), and P N +2 /P N from either (8) or (10) and (11).
Numerically, it is clear that P N → 0 as N → ∞. Some qualitative insight into the rate of decay can be obtained by recalling P N = p N,N and considering the behaviour of p N,k as a function of k. Thus we know from [5] that for large N , the mean number of real eigenvalues E N := k p N,k is to leading order equal to πN/2, and from [7] that the corresponding variance √ 2π e −x 2 /2 , and is thus p N,k is a standard Gaussian distribution after centering and scaling in k by appropriate multiples of √ N . It follows that p N,N is, for large N , in the large deviation regime of p N,k . We remark that this is similarly true of p N,N in the case of eigenvalues of N × N real random Gaussian matrices (i.e. the individual matrices T 1 , T 2 of (7)), for which it is known p N,N = 2 −N (N −1)/4 [5], [6,Section 15.10].
In fact from the exact expression (8) the explicit asymptotic large N form of P N can readily be calculated. For this, let A = e −ζ ′ (−1)+1/12 = 1.28242712...
Theorem 4. For large N , or equivalently Proof. We require the x → ∞ asymptotic expansions of the Barnes G-function [12] and the gamma function Γ(x + 1) = √ 2πx x e For future purposes, we note that a corollary of (32), and the elementary large x expansion is the asymptotic formula To make use of these expansions, we rewrite (8) Furthermore, in the notation (36) it follows from (31) and (32) and further use of (33) (only the explicit form of the leading term is now required) that Γ(N/2 + 1) N G(N + 1) = e −y 2 log(4/e) e y log y+ 1 12 log 2y e 1/6−ζ ′ (−1) 1 + O 1 y .
This corollary follows trivially from Theorem 4. It can however also be derived directly from the recursion formulas in Corollary 3, without use of Theorem 4.