Matrix Whittaker processes

We study a discrete-time Markov process on triangular arrays of matrices of size \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$d\ge 1$$\end{document}d≥1, driven by inverse Wishart random matrices. The components of the right edge evolve as multiplicative random walks on positive definite matrices with one-sided interactions and can be viewed as a d-dimensional generalisation of log-gamma polymer partition functions. We establish intertwining relations to prove that, for suitable initial configurations of the triangular process, the bottom edge has an autonomous Markovian evolution with an explicit transition kernel. We then show that, for a special singular initial configuration, the fixed-time law of the bottom edge is a matrix Whittaker measure, which we define. To achieve this, we perform a Laplace approximation that requires solving a constrained minimisation problem for certain energy functions of matrix arguments on directed graphs.


Introduction
In the last few decades, we have witnessed a surge of research on stochastic integrable models, often motivated by problems in mathematical physics and enriched by deep connections with algebraic combinatorics, representation theory, symmetric functions, and integrable systems [BP14;BG16].Some of the most intensively studied models are interacting particle systems and stochastic growth processes in the Kardar-Parisi-Zhang (KPZ) universality class [Cor16;Zyg22].
From a mathematical perspective, it is natural to consider noncommutative versions of these models, which have very recently received some attention.In [O'C21] a system of interacting Brownian particles in the space of positive definite matrices was considered and shown to have an integrable structure, related to the non-Abelian Toda chain and Whittaker functions of matrix arguments (the latter introduced in that article).In the discrete-time setting, [ABO23] proved Matsumoto-Yor and Dufresne type theorems for a random walk on positive definite matrices.
On the other hand, from the theoretical physics point of view, such matrix models may find interesting applications in quantum stochastic dynamics, as set out in [GBL21].In particular, [GBL21] introduced a matrix generalisation of the classical Kesten recursion and studied a related quantum problem of interacting fermions in a Morse potential.Quoting the authors, their initial motivation was "to explore possible matrix (non-commuting) generalizations of the famous directed polymer problem (which is related to the KPZ stochastic growth equation)".
The subject of the present article is an integrable model of random walks on positive definite matrices with local interactions.This constitutes, on the one hand, a discrete-time analogue of the matrix-valued interacting diffusions studied in [O'C21] and, on the other hand, a matrix generalisation of the log-gamma polymer model.
To motivate the contributions of this article, let us first define a discrete-time exclusion process Z of N ≥ 1 ordered particles Z 1 ≤ Z 2 ≤ • • • ≤ Z N on Z moving to the right.Let (V 1 (n), . . ., V N (n)) n≥1 be a collection of independent random variables supported on Z ≥0 .At each time n, the particle positions are updated sequentially from the 1-st one to the N -th one, as follows.The 1-st particle simply evolves as a random walk on Z with time-n increment V 1 (n).Once the positions of the first i − 1 particles have been updated, if the (i − 1)-th particle has overtaken the i-th particle, then the latter is pushed forward to a temporary position to maintain the ordering; next, to complete its update, the ith particle takes V i (n) unit jumps to the right.The particle locations then satisfy the recursive relations (1.1) If one considers the initial state then the following last passage percolation formula holds: where the maximum is over all directed lattice paths π in Z 2 (i.e., at each lattice site (m, k), π is allowed to head either rightwards to (m + 1, k) or upwards to (m, k + 1)) that start from (1, 1) and end at (n, i).As a process of last passage percolation times, Z can be also associated with the corner growth process with step (or 'narrow wedge') initial configuration.Remarkable integrable versions of this model are those with geometrically and exponentially distributed jumps, first studied in [Joh00].
A positive temperature version of Z can be obtained by formally replacing the operations (max, +) with (+, ×) in the relations (1.1)-(1.2).Namely, given a collection of independent positive random variables (V 1 (n), . . ., V N (n)) n≥1 , we can consider the discrete-time Markov process Z defined by (1.4) (1.5) Considering the initial configuration we have the closed-form expression where the sum is over all directed lattice paths π in Z 2 from (1, 1) to (n, i).The variables (1.7) can be regarded as partition functions of the (1 + 1)-dimensional directed polymer, an intensively studied model of statistical mechanics.Of particular importance is the model with inverse gamma distributed weights V i (n), known as the log-gamma polymer, first considered in [Sep12].In [COSZ14] it was shown that the laws of log-gamma polymer partition functions are marginals of Whittaker measures; the latter are defined in terms of GL d (R)-Whittaker functions and were introduced in that article.
In this article, we study a noncommutative generalisation of the above Markov process of log-gamma polymer partition functions.The 'particles' of this process live in P d , the set of d × d positive definite real symmetric matrices.The random weights V i (n) are now independent inverse Wishart matrices (a matrix generalisation of inverse gamma random variables; see § 1.1).We define Z by setting (1.8) , 2 ≤ i ≤ N , (1.9) where, for a ∈ P d , a 1/2 denotes the unique b ∈ P d such that b 2 = a.The above matrix products are symmetrised to ensure that, starting from any initial configuration Z i (0) ∈ P d , each Z i (n) still belongs to P d for all n ≥ 1.The 1-st particle (1.8) evolves as a (GL dinvariant) multiplicative random walk on P d ; on the other hand, the other particles (1.9) can be viewed as analogous random walks with one-sided interactions.From this point of view, the Markov process as a whole can be also regarded as a noncommutative version of the exclusion process Z defined in (1.1)-(1.2).The natural generalisation of the initial configuration (1.6) is where I d and 0 d are the d×d identity and zero matrices, respectively.Notice that, although all but the first particle are initially zero, the process Z starting from (1.10) lives in P N d at all times n ≥ 1.
In § 3, we introduce a Markov process X = (X(n)) n≥0 , X(n) = (X i j (n)) 1≤j≤i≤N , on triangular arrays of positive definite matrices whose 'right edge', namely (X 1 1 , . . ., X N 1 ), equals Z.The evolution of X may be viewed as a noncommutative version of the dynamics on Gelfand-Tsetlin patterns with blocking and pushing interactions, studied in various contexts in [War07; WW09; Nor10; BF14; BC14; BP16].We refer to Fig. 1 for a graphical representation of such a triangular array.Moreover, as we detail in Remark 3.4, the 'left edge' of X may be regarded as a noncommutative generalisation of the strict-weak polymer studied in [OO15;CSS15].
The first main result of this article (Theorem 3.10) states that, for certain special (random) initial configurations X(0), the 'bottom edge' X N = (X N 1 , . . ., X N N ) of X also has an autonomous Markovian evolution.The transition kernel of X N is explicit and has an interpretation as a Doob h-transform with h-function given by a Whittaker function of matrix arguments.To obtain this, we prove certain intertwining relations between kernels associated to the process X and use the theory of Markov functions (reviewed in Appendix B).Another consequence of these intertwinings is that Whittaker functions are eigenfunctions of certain integral operators and possess a Feynman-Kac type interpretation.
Next, in § 4, we define matrix Whittaker measures on P N d after proving an integral identity of Whittaker functions of matrix arguments (Theorem 4.1), analogous to the wellknown Cauchy-Littlewood identity for Schur functions.The second main result of this article (Theorem 4.8) states that, for a special initial state, the fixed-time law of the bottom edge X N of X is a matrix Whittaker measure on P N d .Such an initial state, designed to match (1.10), is singular, in the sense that the particles are at the 'boundary' of P d .
Due to the singularity of the initial configuration, the proof of Theorem 4.8 will be based on a suitable limiting procedure and a careful integral approximation via Laplace's method.This will require a digression on a constrained minimisation problem for certain energy functions of matrix arguments.We chose to include this analysis in a separate section and to present it in the more general framework of directed graphs, as it may be of independent interest; see § 5.For us, the main application will be the asymptotic formula (4.8) for Whittaker functions of matrix arguments.
From our main results we deduce (see Corollary 4.10) that, under the initial configuration (1.10), the particles of the process Z defined in (1.8)-(1.9)have a fixed-time law given by the first marginal of a matrix Whittaker measure on P N d .In the scalar d = 1 case, we recover the aforementioned result of [COSZ14] for the law of the log-gamma polymer partition functions.In Corollary 4.10, we also obtain an analogous result concerning the fixed-time law of the 'left edge' of the triangular array X.
It is worth mentioning that the log-gamma polymer partition functions (1.7) were also studied in [COSZ14] as embedded in a dynamic on triangular arrays.However, such a dynamic was constructed via the combinatorial mechanism of the geometric Robinson-Schensted-Knuth correspondence; in particular, at each time step, the right edge is updated using N new (independent) random variables, whereas all the other components are updated via deterministic transformations of the current state and the newly updated right edge.It turns out that, for d = 1, the processes considered in [COSZ14] and in the present article have an identical right edge and, under the special initial configuration of Theorem 3.10, also a bottom edge process with the same Markovian evolution.However, even in the d = 1 case, the two processes, as a whole, differ.The dynamic introduced in this article is driven by random updates with N (N + 1)/2 degrees of freedom, since each particle of the triangular array is driven by an independent source of randomness (as well as by local interactions with the other particles).
Organisation of the article.In § 2, we define Whittaker functions of matrix arguments.In § 3, we introduce a Markov dynamic on triangular arrays of matrices and study the evolution of its bottom edge, using the theory of Markov functions; we also obtain a Feynman-Kac interpretation of Whittaker functions.In § 4, we define matrix Whittaker measures (through a Whittaker integral identity) and prove that they naturally arise as fixed-time laws in the aforementioned triangular process under a singular initial configuration.To do so, we need a Laplace approximation of Whittaker functions, which can be justified by solving a constrained minimisation problem for certain energy functions of matrix arguments on directed graphs: this is the content of § 5.In Appendix A, we give a proof of the Cauchy-Littlewood identity for Schur functions that resembles our proof of the Whittaker integral identity.In Appendix B, we review the theory of Markov functions for inhomogeneous discrete-time Markov processes.Finally, in Appendix C, we prove a convergence lemma related to weak convergence of probability measures.
1.1.Notation and preliminary notions.Here we introduce some notation and preliminary notions that we use throughout this work.For background and proofs, we refer to [HJ13;Ter16].The following properties hold: • if x ∈ P d and λ > 0, then λx ∈ P d ; • if x, y ∈ P d , then x + y ∈ P d (but in general xy / ∈ P d ); For x ∈ P d , there exists a unique y ∈ P d such that y 2 = x; we denote such a y by x 1/2 .For any y ∈ P d , we define the (noncommutative) 'multiplication operation' by y as (1.11)Such a symmetrised product will be used to construct a multiplicative random walk on P d (see Definition 3.1 and Remark 3.2 below).We also denote by I d and 0 d the d × d identity matrix and zero matrix, respectively.
Measure and integration on P d .Let GL d be the group of d×d invertible real matrices.Define the measure µ on P d by where dx i,j is the Lebesgue measure on R in the variable x i,j .Such a measure is the GL d -invariant measure on P d , in the sense that for all a ∈ GL d and for all suitable functions f .In other words, µ is invariant under the group action of GL d on P d Furthermore, the measure µ is preserved under the involution x → x −1 .
Wishart distributions and gamma functions.For α > The inverse of a Wishart matrix with parameter α has the distribution on P d .We will refer to the latter as the (d-variate) inverse Wishart distribution with parameter α.
Kernels and integral operators.Let (S, S) and (T, T ) be two measurable spaces.Let mS denote the set of complex-valued measurable functions on (S, S).For our purposes, a kernel from T to S will be a map L : T × S → C such that, for each t ∈ T , L(t; •) is a (complex) measure on (S, S) and, for each A ∈ S, L(•; A) is an element of mT .The kernel L can be also, alternatively, thought of as an integral operator whenever the integral is well defined.Clearly, the composition of kernels/operators yields another kernel/operator; such a composition is associative but, in general, not commutative.When the complex measure L(t; •) is a probability measure for all t ∈ T , we will talk about Markov kernels/operators.Throughout this article, the measurable spaces will be usually Cartesian powers of P d (which we denote by P k d , k ≥ 1), with their Borel sigma-algebras.Moreover, for a kernel L from P k d to P ℓ d , the measure L(t; •) will be, in most cases, absolutely continuous with respect to the reference product measure µ ⊗ℓ on P ℓ d , for any t ∈ P k d ; with a little abuse of notation, we will then also write s → L(t; s) for the corresponding density (a measurable function on P ℓ d ).

Whittaker functions
In this section we define Whittaker functions of matrix arguments following [O'C21], and then extend them to a further level of generality.Notice also that the kernels (2.8) and (2.15) defined below are matrix versions of certain kernels defined in [COSZ14, § 3.1] and [OSZ14, § 2] (see also references therein).
, read from right to left.The arrows refer to the energy function Φ N (x) in (2.3), where every summand tr[ab −1 ] corresponds to an arrow pointing from a to b in the figure .2.1.Whittaker functions of matrix arguments.We define Whittaker functions of matrix arguments as integrals over 'triangular arrays' of d × d positive definite matrices.For N ≥ 1, denote by where For a graphical representation of the array (2.1) and of the 'energy function' Φ N , see Fig. 1.For z = (z 1 , . . ., z N ) ∈ P N d , let T N d (z) ⊂ T N d be the set of all height-N triangular arrays x with N -th row x N = z.We define the Whittaker function ψ N λ (z) with argument z ∈ P N d and parameter λ ∈ C N as (2.4) Notice that, for N = 1, the expression above reduces to ψ 1 λ (z) = |z| −λ .As proved in [O'C21], the integral (2.4) is absolutely convergent for all λ ∈ C N , so that Whittaker functions are well defined.
For our purposes, it is convenient to rewrite Whittaker functions in terms of certain kernels that we now introduce.For N ≥ 1, λ ∈ C N and x ∈ T N d , define the kernel where, as always from now on, i : j denotes the tuple (i, i + 1, . . ., j − 1, j) for i ≤ j, so that x 1:(N −1) ∈ T N −1 d is the triangular array consisting of the first N − 1 rows of x.Notice that, for N = 1, (2.5) reduces to Σ 1 λ (z; ∅) = |z| −λ = ψ 1 λ (z).For z ∈ P N d , let us also define the kernel where δ is the Dirac delta kernel on P N d .Then, the Whittaker function (2.4) can be written as (2.8) We will usually regard (2.8) as a kernel by setting K N b (z; dy) := K N b (z; y)µ ⊗(N −1) (dy).We then have, for λ ∈ C N , z ∈ P N d , and (2.9) This yields a recursive definition of Whittaker functions: (2.10) 2.2.A generalisation of Whittaker functions.We now introduce a generalisation of Whittaker functions of matrix arguments, which will naturally emerge in § 4.3 and, in the scalar case d = 1, corresponds to the one considered in [OSZ14].These generalised Whittaker functions are integrals over trapezoidal arrays of positive definite matrices, similarly to how the Whittaker functions of § 2.1 are defined as integrals over triangular arrays.
Let n ≥ N ≥ 1 and denote by the set of trapezoidal arrays (2.11) with i-th row Let us now give an equivalent representation of these generalised Whittaker functions.The following kernel will play a central role in this work.For a ∈ C and z, z ∈ P N d , set (2.15) We will see P N a (z; z) as a measure in either of the two arguments, defining (2.16) We then have (2.17) We also record here two relations between the kernels (2.8) and (2.15), which follow directly from the definitions: for y = (y 1 , . . ., y N −1 ) ∈ P N −1 d , s ∈ P d , and z = (z 1 , . . ., z N ) ∈ P N d .Taking a = λ N in (2.19), multiplying both sides by ψ N −1 (λ 1 ,...,λ N−1 ) (y), integrating over P N −1 d with respect to µ ⊗(N −1) (dy), and using (2.10) and (2.17), we obtain the identity (2.20) Remark 2.1.Let us mention that we anticipate the function ψ N,n λ;s to be symmetric in the parameters λ 1 , . . ., λ n .This is not obvious from the definition, but it is suggested by an integral identity of Whittaker functions of matrix arguments that will be proven later on (see (4.1)).As argued in [O'C21, § 7.1], this symmetry is true at least in the case N = n = 2.Moreover, it is known for d = 1 and arbitrary n, N ; see, for example, [KL01], [GLO08] and [OSZ14,.

Markov dynamics
In this section, we define a Markov process X on triangular arrays, which can be viewed as a system of interacting random walks on P d .Next, we prove intertwining relations between certain transition kernels related to this process.This implies, via the theory of Markov functions, that, under certain random initial configurations, the bottom edge of the triangular process X has an autonomous stochastic evolution.A consequence of these results is that Whittaker functions of matrix arguments are eigenfunctions of certain integral operators and, thereupon, admit a Feynman-Kac interpretation.

Interacting
Remark 3.2.The random walk R of Definition 3.1 is indeed GL d -invariant, in the sense that the conjugated walk (g ⊤ R(n)g) n≥0 has the same transition kernels for any choice of g ∈ GL d (cf.[ABO23, § 3]).Instead of (3.1), one could consider a different process through the alternative symmetrisation One can check that the resulting random walk R ′ is O d -invariant, but in general not GL dinvariant.In principle, one could proceed to obtain analogous results to those presented in the present article using this alternative symmetrisation (for a similar approach in the continuous Brownian setting, see [O'C21, Prop.3.5]).However, from our point of view, the choice (3.1) is the most natural and leads to more explicit transition kernels throughout.
It is well known that the Wishart distribution (1.13) and the inverse Wishart distribution (1.14) are O d -invariant.In this article, we will focus on GL d -invariant random walks with inverse Wishart increments.
Recall from definition (2.15 Define then the renormalised kernel It is immediate to see that the (time-homogeneous) GL d -invariant random walk on P d with inverse Wishart increments of parameter a > d−1 2 has transition kernel P 1 a .We now define a discrete-time Markov process X = (X(n)) n≥0 on the set T N d of height-N triangular arrays whose components are elements of P d .
) be an inverse Wishart random matrix with parameter α(n)+β i (the same parameter across j); assume further that all these random matrices are independent of each other.We define the process X = (X(n)) n≥0 , where The i-tuple X i := (X i 1 , . . ., X i i ) will be referred to as the i-th row of X.
The fact that each X i j (n) takes values in P d follows by standard properties of positive definite matrices (cf.§ 1.1).Notice that, adopting the convention X i 0 (n) −1 = X i i+1 (n) = 0 d for all i ≥ 0 and n ≥ 0, then the last formula in (3.4) can be taken as the definition of The dynamic on T N d defined by (3.4) implies that the 'top particle' X 1 1 evolves as a GL d -invariant random walk in P d with inverse Wishart increments (W 1 1 (n)) n≥1 .Furthermore, the 'right edge' process (X 1 1 , X 2 1 , . . ., X N 1 ) equals the system (Z 1 , . . ., Z N ) of random particles in P d with one-sided interactions defined in (1.8)-(1.9),where the random weight ) also evolves as a system of particles in P d with one-sided interactions, as we now explain.
and the process L = (L 1 , . . ., L N ) satisfies the recursions Under the (singular) initial configuration one can see by induction that reduces to a sum of independent Wishart matrices: In particular, L i (i) has the Wishart distribution with parameter i j=1 (α(j) + β j ).
Remark 3.4.We make a few remarks about various specialisations of the process X and related Markov dynamics: (i) The interacting diffusion model on positive definite matrices studied in [O'C21] (see also [O'C12, § 9] for the d = 1 case) can be regarded as a continuous-time analogue of the process X defined in (3.4).
(ii) It seems that even the d = 1 case of the dynamic (3.4) has not been explicitly considered elsewhere.It is related, even though not identical, to the process constructed in [COSZ14] via the geometric Robinson-Schensted-Knuth correspondence; see the discussion in the introduction for further details.
(iii) For d = 1, under the 'step' initial configuration, the right edge can be regarded as a process of log-gamma polymer partition functions; see (1.6)-(1.7)and the discussion therein.
(iv) For d = 1, under the 'step' initial configuration (3.7), the left edge can be regarded as a process of strict-weak polymer partition functions in a gamma environment, studied in [OO15; CSS15].A strict-weak path is a lattice path π that, at each lattice site (m, k), is allowed to head either horizontally to the right to (m + 1, k) or diagonally up-right to (m + 1, k + 1).It is easily seen that the process L defined in (3.5)-(3.6), in the d = 1 case, takes the closed form expression where the sum is over all strict-weak paths π from (0, 1) to (n, i), the product is over all edges e in the path π, and d e is a weight attached to the edge e and defined as follows: . Formula (3.8) defines the strict-weak polymer partition function.
(v) The d = 1 case of (3.4) is a 'positive temperature' analogue (equivalently, a (+, ×) version) of the process defined by where W i j (n) are non-negative random variables representing jumps to the right (see e.g.[WW09]).Roughly speaking, particle X i j performs a random walk subject to certain interactions with other particles: it is pushed by X i−1 j and blocked by X i−1 j−1 .(vi) Besides [WW09], other works [War07; Nor10; BF14; BC14; BP16] studied, in various discrete and continuous settings, similar push-and-block dynamics on Gelfand-Tsetlin patterns driven by random updates with N (N +1)/2 degrees of freedom.In particular, again in the case d = 1, the process X should correspond to a certain q → 1 scaling limit of the q-Whittaker processes studied in [BC14;BP16].
Motivated to obtain the explicit Markovian evolution of X, we now introduce the following kernels.For a ∈ C, y = (y 1 , . . ., y N −1 ) with the convention (3.12) Therefore, under the above conditions on the parameters, one can renormalise these kernels, so that they integrate to 1: The following result can be easily verified using the construction of X in Definition 3.3.
Proposition 3.5.Let X as in Definition 3.3.Then, the conditional distribution of 3.2.Intertwining relations.We will now show that the Markov dynamic on X (see Definition 3.3), when started from an appropriate random initial state, induces an autonomous Markov dynamic on the N -th row, or 'bottom edge', of X.This will be a consequence of an intertwining relation between kernels through the theory of Markov functions, which is reviewed in Appendix B for the reader's convenience.
holds as an equality between kernels from Proof.We have to prove that (3.18) Using the definitions of K N b and Λ N a,b , we obtain, after some rearrangements and cancellations, that the left-hand side of (3.18) equals , with the usual convention y N = 0.By interchanging the order of integration, we see that the latter display equals , where J : P 3 d → C is defined by By the properties of positive definite matrices (see § 1.1), we have that w −1 − s −1 ∈ P d if and only if s − w ∈ P d ; moreover, for w ∈ P d , the latter condition is stronger than s ∈ P d .
We then make the change of variables s ′ := s − w, which preserves the Lebesgue measure on the 'independent' entries of the symmetric matrix s, so that Therefore, we have After the further, this time µ-preserving, change of variables s ′′ := T u −1 +v −1 (s ′ ), we obtain where the gamma function is well defined since by hypothesis ℜ(a + b) > (d − 1)/2.After a few cancellations, we then see that the left-hand side of (3.18) equals .
It now follows from the definitions that this equals the right-hand side of (3.18), thus concluding the proof.
A simple inductive argument shows that the intertwining (3.17) can be extended to an intertwining that involves the Π-kernel (3.10) and the Σ-kernel (2.6).From now on, we fix As usual, we also use the notation a + λ := (a + λ 1 , . . ., a + λ N ).
Corollary 3.7.The intertwining relation holds as an equality between kernels from P N d to T N d .
Proof.Taking into account (2.6), it is immediate to see that (3.19) is equivalent to (3.20) for all z ∈ P N d and test function f : , where we set To prove (3.20), we proceed by induction.For N = 1, (3.20) amounts to the identity ψ 1 λ (z) for z ∈ P d and f : P d → R. Using (2.10) and (2.15), one can easily verify that the latter is true, as both sides equal |z| a Let now N ≥ 2 and λ = (λ 1 , . . ., λ N −1 ).Assume by induction that for any y ∈ P N −1 d and any test function g : in (3.21) and integrating both sides with respect to the measure K N λ N (z; •), we obtain Using (2.9) and (3.10) for the left-hand side and (3.16) for the right-hand side, and interchanging the integration order, we then have where the latter two equalities follow from (3.18) and (2.9), respectively.The identification concludes the proof of (3.20).
Recall now that the Σ-kernels generate Whittaker functions of matrix arguments, in the sense of (2.7).By integrating the intertwining relation (3.19) and using (3.12), we immediately deduce that Whittaker functions are eigenfunctions of the integral P -operators: Corollary 3.8.We have We note that this complements the interpretation of the Whittaker functions ψ N λ , given in [O'C21], as eigenfunctions of a differential operator, namely the Hamiltonian of a quantisation in P N d of the N -particle non-Abelian Toda chain.For x ∈ T N d and z, z ∈ P N d , we now define Corollary 3.9.The intertwining relation holds as an equality between kernels from From a probabilistic point of view, (3.25) states that, for any fixed z ∈ P N d , the two following update rules are equivalent: (i) starting the process X from a (random) initial configuration dictated by the intertwining kernel Σ(z; •) and letting it evolve according to the dynamic Π; and (ii) running the dynamic P on the bottom edge (started at z) and then updating the whole triangular array according to the intertwining kernel Σ.The main result of this section is a precise account of this interpretation.
Theorem 3.10.Let X = (X(n)) n≥0 be the Markov process on T N d as in Definition 3.3.Assume that, for an arbitrary z ∈ P N d , the initial state X(0) of X is distributed according to the measure Σ N β (z; •).Then, the N -th row X N = (X N (n)) n≥0 is a time-inhomogeneous Markov process (in its own filtration) on the state space P N d , with initial state z and time-n transition kernel P N α(n),β .Moreover, for any bounded measurable function f : for all n ≥ 1 .
These follow from Corollary 3.9.
Remark 3.11.By letting N vary, it is immediate to deduce from Theorem 3.10 that every row of X evolves as a Markov process in its own filtration, under an appropriate (random) initial configuration on the previous rows.Therefore, the focus on the N -th row should only be seen as a convenient choice.

Feynman-Kac interpretation.
Here we provide a Feynman-Kac type interpretation of Whittaker functions based on the eigenfunction equation (3.22).Our result should be compared to the one obtained in [O'C21, Prop.9] in the continuous setting of Brownian particles.
random walk on P d with initial state Y i (0) = y i and inverse Wishart increments with parameter λ i .
Recalling (3.3), Y is then a time-homogeneous Markov process starting at y with transition kernel For z, z ∈ P N d , define the sub-Markov kernel where V is the 'killing potential' (3.28) Denote by P y and E y the probability and expectation, respectively, with respect to the law of Y with initial state y.Theorem 3.13.For all y ∈ P N d , we have The main purpose of this subsection is to prove (3.29).In a nutshell, using a fairly standard martingale argument, we will show that the expectation in (3.29) is the unique solution to an eigenproblem; the latter is also, essentially, solved by Whittaker functions.
Lemma 3.14.Fix an integer ℓ ≥ 0. For any y ∈ P N d , we have Remark 3.15.In particular Lemma 3.14 with ℓ = 1 implies that the infinite series inside the expectation in (3.29) converges P y -a.s.
Proof of Lemma 3.14.Since it suffices to show that, for each 1 Let us record the following properties, which hold for any a, b ∈ P d : • tr[ab] ≤ tr[a] tr[b] (submultiplicativity of the trace); • tr[a] ≤ d λ max (a); Here, λ max and λ min denote the maximum and minimum eigenvalue, respectively.Using these facts, we have, for 1 P y -a.s., where ψ is the digamma function.These are the maximum (respectively, minimum) Lyapunov exponent of a GL d -invariant random walk with inverse Wishart increments of parameter λ i+1 (respectively, λ i ).We then obtain lim sup since the digamma function is strictly increasing and, by Definition 3.12, Then, for all y ∈ P N d , Proof.Consider the process Y as in Definition 3.12, with initial state y ∈ P N d and transition kernel Θ N λ .Denote by (F(k)) k≥0 its natural filtration.It follows from the eigenfunction equation that Therefore, the process M = (M (k)) k≥0 defined by is an (F(k)) k≥0 -martingale.By the boundedness property, M is uniformly bounded and, thus, a uniformly integrable martingale.Therefore, M converges P y -a.s. and in 1-norm to a certain limit M (∞) and, for all k ≥ 0, we have . By Lemma 3.14 (with ℓ = 0), we have lim k→∞ V (Y (k); Y (k)) = 0, P y -a.s.The boundary condition then implies lim k→∞ u(Y (k)) = 1, P y -a.s., whence We conclude that, for any y ∈ P N d , Proof of Theorem 3.13.It was proven in [O'C21, proof of Prop.9] that the function is bounded and satisfies lim By Lemma 3.16, it then remains to prove that Θ N λ v = v.It follows from the definition (2.15) of the kernel P N a that for z ∈ P N d .Using the eigenfunction equation (3.22), we see that the right-hand side above equals v(z), as desired.
Proof.We may assume that N = 2, so that Y = (Y 1 , Y 2 ) starts at y = (y 1 , y 2 ).Using Theorem 3.13 and the definition of Whittaker functions, we compute the Laplace transform of the left-hand side of (3.31) as for s ∈ R, where we used the change of variables z = y 1 .The last integral equals E e −s tr [aZ] , where Z is inverse Wishart of parameter λ 2 − λ 1 .We conclude that the two sides of (3.31) have the same Laplace transform and, hence, the same law.
Remark 3.18.Up to some technical details, identity (3.31) may be also deduced from the Dufresne type identity for a random walk on P d proved in [ABO23].Let (R(n)) n≥0 be a GL d -invariant random walk on P d whose initial state R(0) is an inverse Wishart matrix with parameter λ 2 and whose increments are Beta type II matrices with parameters λ 1 and λ 2 (see [ABO23] for more details).It is then natural to expect that the eigenvalue processes of the two processes (Y 1 (n , have the same law; this is certainly true at least in the case d = 1.By summing the traces of these two processes over all n ≥ 0, [ABO23, Theorem 4.10] would then immediately provide a proof of (3.31) that does not rely upon the Feynman-Kac formula (3.29).See [O'C21, Lemma 8] for an analogous argument in the Brownian setting.

Fixed-time laws and matrix Whittaker measures
In this section, we first prove a Whittaker integral identity that allows us to introduce matrix Whittaker measures.We then obtain an asymptotic formula for a Whittaker function whose arguments go to zero or infinity in norm.Using the latter result, we next show that, for a certain singular initial state, matrix Whittaker measures appear naturally as the fixed-time laws of the bottom edge of the triangular process X introduced in § 3.1.Finally, under the same singular initial state, we study the fixed-time law of the right edge and of the left edge of X.

Matrix Whittaker measures. Whittaker functions of matrix arguments satisfy a remarkable integral identity:
The case N = n of (4.1) was noted in [O'C21, Prop.10]; however, the identity did not play a key role in that article, and the details of the proof were not provided therein.Below we provide a proof of the general case n ≥ N that involves the generalised Whittaker functions introduced in § 2.2.
In the scalar d = 1 setting, (4.1) goes back to [OSZ14, Corollary 3.5].For d = 1 and N = n, it is equivalent to an identity that was originally found in the number theoretic literature [Bum89;Sta02].
Theorem 4.1 can be also seen as an analogue, in the context of Whittaker functions, of the celebrated Cauchy-Littlewood identity for Schur functions.In the literature on symmetric functions, the latter is usually proved using either the determinantal structure of Schur functions (see [Mac79, I-(4.3)]) or the Robinson-Schensted-Knuth correspondence, a combinatorial bijection (see [Sta99, Theorem 7.12.1]).None of these tools is available, so far, in our context.To prove (4.1), we will rather proceed inductively, using the recursive definition of Whittaker functions and the eigenfunction equation (3.22).For the reader's convenience, we also include in Appendix A a proof of the classical Cauchy-Littlewood identity that similarly relies on an eigenfunction equation for Schur functions (which can be seen as a version of the so-called Pieri rule).
Proof of Theorem 4.1.We will prove (4.1) by induction on n.For a fixed integer n ≥ 1, let S(n) be the statement that (4.1) holds for all N such that n ≥ N ≥ 1 and for any choice of λ and ̺ satisfying the assumptions of the theorem.
For n = N = 1 we have where we have used the definitions of Whittaker functions and gamma functions and the µ-preserving change of variables z = T s −1 (z).This proves the base case S(1).Suppose now by induction that S(n − 1) holds for some fixed n ≥ 2. To prove that S(n) holds, let us first prove that (4.1) is valid for all N such that n > N ≥ 1.It follows from (2.17), Fubini's theorem, (2.16) and (3.22), that Since n − 1 ≥ N , using the assumption S(n − 1) in the latter integral we obtain (4.1).
Recall that we have already proved (4.1) for all N such that n > N ≥ 1. Applying this, for N = n − 1, to the latter integral, we conclude that (4.1) holds also for N = n.
We call matrix Whittaker measure with parameters λ and ̺ the measure on P N d that is absolutely continuous with respect to µ ⊗N (dz) with density where I d is the d × d identity matrix.According to the usual convention, we also denote by W N,n λ,̺ (dz) the measure itself.
By Theorem 4.1, (4.2) defines a probability distribution on P N d .This extends the definition of matrix Whittaker measures given in [O'C21, § 7.4], which corresponds to the case n = N : 4.2.Asymptotics of Whittaker functions.For any real k > 0, let and let r i (k) := (r i 1 (k), . . ., r i i (k)).Our ultimate goal is to obtain the k → ∞ leading order approximation of the Whittaker function ψ N λ (r N (k)).We rely on some results (Theorem 5.13 and Prop.5.14) that we will prove, in a more general setting, in § 5.With this purpose in mind, we use the graphical representations of the set of height-N triangular arrays T N d and of the energy function Φ N , both involved in the definition of the Whittaker function (2.4) (see Fig. 1).Given N ≥ 2, we set and consider the finite graph G = (V , E), where E consists of all (directed) edges (i, j) → (i + 1, j) and (i + 1, j + 1) → (i, j), for 1 ≤ j ≤ i ≤ N − 1.Then, T N d may be identified as the set P V d of arrays x = (x v ) v∈V , where each x v ∈ P d .Let also We may thus identify z ∈ P N d with z ∈ P Γ d , so that the set T N d (z) of all height-N triangular arrays whose N -th row equals z coincides with the set P V d (z), according to the notation (5.4).Furthermore, the energy function (2.3) can be equivalently rewritten as All the results of § 5 hold for the above 'triangular graph' structure, since: ) is an acyclic finite directed graph; • Γ is a proper subset of V containing the only source (N, N ) and sink (N, 1) of G; • the energy function Φ N is of the form (5.9).
We first prove a property of the critical points of Φ N that, in the scalar d = 1 setting, was observed in i be the determinant of the product of the i-th row of x.Then, Proof.The critical point equations of the energy function Φ N are with the convention x i−1 i = (x i−1 0 ) −1 = 0 for all 1 ≤ i < N (these correspond to (5.15) in the case of the triangular graph G).Taking determinants of both sides, we obtain Taking the product over j in the latter, many terms cancel out, yielding i j=1 By definition of p 1 , . . ., p N , the latter can be written as with the convention p 0 := 1.Finally, it is straightforward to see that equations (4.6) are equivalent to (4.5).
Let now As the components of I N d are scalar matrices, Theorem 5.13 implies: Throughout this subsection, m will always denote the above minimiser.Theorem 4.6.For any λ ∈ C N , we have where |H(m)| > 0 is the Hessian determinant of Φ N at m.
The case d = 1, N = 2 of this asymptotic result is classical; the case d = 1 and general N can be found in [O'C12, eq.20].Finally, the case d > 1, N = 2 may be inferred from the Laplace approximation of Bessel functions of matrix arguments studied in [BW03] (see also [Gra18,Appendix B] and [O'C21, Section 2.6]).
An important feature of (4.8) is that the leading order asymptotics does not depend on the parameter λ.This was already remarked in [BH95] in the special case d = 1 and N = 3, for which the full asymptotic expansion was obtained.
Proof of Theorem 4.6.By (2.4), we have Recalling (4.4), let us change variables by setting One can then easily verify, using also the invariance property of the measure µ, that (4.10) Applying Prop.5.14 with g := ∆ N λ , we obtain since the number of vertices of G that do not belong to Γ is N (N − 1)/2.The claim then follows from Corollary 4.5 (which, in particular, implies that ∆ N λ (m) = 1).
Recall now the definition (2.6) of the Σ-kernel.
Proof.As the leading order asymptotics of the Whittaker function ψ N λ (r N (k)) does not depend on λ by Theorem 4.6, we have Therefore, it suffices to prove (4.12) for ̺ = λ.Note that, using (4.10) and the fact that λ ∈ R N , the measure µ N k defined by is a probability measure on T N d (I N d ).By definition of Σ λ , we then have where in the integral we performed the change of variables (4.9).Since f is bounded and continuous, the functions {f k } k>0 are uniformly bounded and continuous; moreover, by assumption, they converge as k → ∞ to f ∞ uniformly on any compact subsets of T N d (I N d ).Therefore, by Lemma C.1, it is now enough to show that µ N k converges weakly as k → ∞ to the Dirac measure δ m , i.e. that lim for every bounded and continuous function g : This claim, in turn, follows readily from Prop.5.14, since, without loss of generality, one can assume g(m) = 0. 4.3.Fixed-time law of the 'bottom edge' process.Let us now go back to the Markov process X on T N d from Definition 3.3.Recall that, under the hypotheses of Theorem 3.10, the N -th row X N of the process X has an autonomous Markov evolution with time-n transition kernel P N α(n),β (cf.(3.24)).The transition kernel of X N from time 0 to time n is then given by the composition 2 for all 1 ≤ ℓ ≤ n and 1 ≤ i ≤ N .Iterating the eigenfunction equation (3.22) n times, we obtain the following eigenfunction equation for U N,n α,β : Consider now the initial state X N (0) = r N (k) (cf.(4.4)), which becomes singular in the limit as k → ∞.We will show that the measure U N,n α,β (r N (k); •) converges, as k → ∞, to the matrix Whittaker measure with parameters (α(1), . . ., α(n)) and β.An intuition about this fact is provided by (4.14).It follows from Theorem 4.6 that the ratio of Whittaker functions on the right-hand side of (4.14), evaluated at r N (k), converges to 1 as k → ∞.It is then easy to see that, if the convergence to matrix Whittaker measures holds as claimed above, then (4.14) reduces to the Whittaker integral identity proved in § 4.1.
Proof.We will prove that This statement is stronger than (4.15), as the supremum in (4.15) is clearly bounded from above by the integral in (4.16).
The latter expression vanishes as k → ∞ by the induction hypothesis, thus proving the induction step.It remains to prove the base case, i.e. (4.16) for n = N .Recall that the measures U N,N α,β (r N (k); •), for any k > 0, and W N,N α,β have the same finite total mass, since they are all probability distributions, and are absolutely continuous with respect to µ ⊗N .By Scheffé's theorem (see e.g.[Bil95, Theorem 16.12]), it then suffices to show the convergence of the densities: (4.17) Fix z ∈ P N d once for all.Using (4.13), we write the term e − tr[ab −1 ] in the integral.Relabelling the z i j as in (4.21)-(4.22)yields two triangular arrays x (coloured in red) and y (coloured in blue).

Define now
for z 0 , z N ∈ P N d , with the usual conventions z i N +1 := 0 for all i = 0, . . ., N .Using the definition (3.24) of the P -kernels, we then have Comparing (4.19) with (4.3), we are reduced to show that lim k→∞ Let us relabel the variables in the integral (4.18) by setting This relabelling yields two triangular arrays x, y ∈ T N d .See Fig. 3 for a graphical representation of the variables z i j and the corresponding arrays x and y.Recalling the definition (2.6) of the Σ-kernel, we have where α(1 : N ) := (−α(N ), . . ., −α(1)) and the function f : Here, each term e corresponds, graphically, to a bold arrow in Fig. 3.We now wish to apply Corollary 4.7.Notice first that f is a continuous function of x; moreover, it is bounded below by 0 and above by ψ N α(1:N ) (z) (cf.(2.4)).The associated functions f k defined in (4.11) are By dominated convergence and by the definition (2.4) of Whittaker function, we have the pointwise convergence Notice that {f k } k>0 is a collection of continuous functions, increasing with k, that converges pointwise to a continuous limit; hence, by Dini's theorem (see e.g.[Rud76, Theorem 7.13]), the convergence is uniform on compacts.Then, the assumptions of Corollary 4.7 are satisfied and we have where m is the unique global minimiser of Φ N on T N d (I N d ) (cf.Corollary 4.4).Since m 1 1 = I d by Corollary 4.5, we have This yields the desired limit (4.20).
4.4.Fixed-time laws of the 'right edge' and 'left edge' processes.Throughout this subsection, it will be convenient to work with the space of d × d positive semidefinite matrices, i.e. d × d real symmetric matrices with nonnegative eigenvalues; such a space is the closure of P d under the standard Euclidean topology, and we thus denote it by P d .
It is clear from the definition given in § 3.1 that the 'right edge' X 1 = (X 1 1 , . . ., X N 1 ) of X is a Markov process in its own filtration.Furthermore, as mentioned before, X 1 equals the system Z = (Z 1 , . . ., Z N ) of random particles in P d with one-sided interactions defined in (1.8)-(1.9),where the random weight V i (n) equals W i 1 (n), an inverse Wishart random matrix with parameter α(n) + β i .If the initial state Z(0) of this process is in P N d (respectively, P N d ), then clearly Z evolves as a process in P N d (respectively, P N d ).Analogously, the 'left edge' of X is a Markov process in its own filtration.Its 'inverse' L = (L 1 , . . ., L N ) := ((X 1 1 ) −1 , . . ., (X N N ) −1 ) is given by (3.5)-(3.6),where ).As the next lemma shows, the singular initial state of the bottom edge of X considered in § 4.3 induces (through Theorem 3.10) the initial state (1.10) on the right edge X 1 , which resembles the step or 'narrow wedge' initial configuration in systems of interacting particles/random walks.A similar statement holds for the left edge.
We now wish to apply Corollary 4.7.Since g is bounded and continuous, f also is.The associated functions f k defined in (4.11) are These functions converge as k → ∞ to f ∞ (x) := g(x 1 1 , 0 d , . . ., 0 d ) uniformly on compacts, since g is continuous on P As a consequence of Theorem 4.8 and Lemma 4.9, we obtain: Corollary 4.10.As above, let Z = (Z 1 , . . ., Z N ) and L = (L 1 , . . ., L N ) be the right edge process and the (inverse) left edge process, respectively, with initial states Then, for n ≥ N , Z N (n) and L N (n) are distributed as the first marginal and the N -th marginal, respectively, of the matrix Whittaker measure with parameters (α(1), . . ., α(n)) and β.
Proof.Again, we only prove the claim about the right edge, as the proof of the claim about the left edge is completely analogous.
Let the process X be as in Definition 3.3, with initial state X(0) distributed according to Σ β (r N (k); •).It is clear from the definition that X N 1 (n) can be written as a continuous, deterministic function of the right edge initial state X 1 (0) and of the collection of random matrices (W 1 1 (ℓ), . . ., W N 1 (ℓ)) 1≤ℓ≤n .Therefore, by Lemma 4.9 and the continuous mapping theorem [Bil99, Theorem 2.7], X N 1 (n) converges in law as k → ∞ to Z N (n).On the other hand, by Theorem 4.8, for n ≥ N , X N 1 (n) converges in law as k → ∞ to the first marginal of a matrix Whittaker measure with parameters (α(1), . . ., α(n)) and β.
Remark 4.11.The following generalisation of Corollary 4.10 is immediate: under the same hypotheses, for every 1 ≤ i ≤ N and n ≥ i, Z i (n) is distributed as the first marginal of the matrix Whittaker measure with parameters (α(1), . . ., α(n)) and (β 1 , . . ., β i ).This is due to the fact that, by definition, for any fixed i ≥ 1, the process (Z 1 , . . ., Z i ) has both an initial configuration (I d , 0 d , . . ., 0 d ) and a Markov evolution that do not depend on the choice of N ≥ i. Analogously, for every 1 ≤ i ≤ N and n ≥ i, L i (n) is distributed as the i-th marginal of the same matrix Whittaker measure.

Minimisation of energy functions and Laplace approximations
In this section, we study minimisation problems for certain energy functions of matrix arguments on directed graphs.As a consequence, we obtain Laplace approximations for integrals of exponentials of these energy functions.For our purposes, the most important application of such results consists in certain asymptotics of Whittaker functions of matrix arguments; see § 4.2.However, the results of this section may be of independent interest.For instance, the general framework we work with may be applied to obtain analogous asymptotics for orthogonal Whittaker functions, which also appeared in the study of stochastic systems -see [BZ19; BW23].5.1.Energy functions on directed graphs.Let us recall some terminology of graph theory that will be useful throughout this section.A finite directed graph G = (V, E) is a pair consisting of a nonempty finite set V of vertices and a set E ⊂ {(v, w) ∈ V 2 : v = w} of edges.Note that edges connecting a vertex to itself are not allowed, nor are multiple edges.The direction of an edge (v, w) connecting v to w is given by the ordering of the pair.For the sake of notational convenience, we also write v → w when (v, w) ∈ E, and v → w when (v, w) / ∈ E. A vertex v is called a sink if it has no outcoming edges (i.e. if v → w for all w ∈ V ) and a source if it has no incoming edges (i.e. if w → v for all w ∈ V ).For any v, w ∈ V and 0 ≤ l < ∞, we call path of length l in G from v to w any sequence (v 0 , v 1 , . . ., v l ) such that v 0 = v, v l = w, and v i−1 → v i for all 1 ≤ i ≤ l.A cycle is any path (v 0 , v 1 , . . ., v l ) such that v 0 = v l and any other two vertices are distinct.We say that G is acyclic if it has no cycles.From now on, throughout the whole section, G = (V, E) will always be an acyclic finite directed graph.
Lemma 5.1.For all v ∈ V , there exists a path in G from v to a sink; moreover, there exists a path in G from a source to v.
Proof.We will prove the existence of the first path only, as the existence of the second path follows from a similar argument.We construct the path algorithmically.Set v 0 := v.For all i = 0, 1, 2, . . ., we proceed as follows: if v i is a sink, then we stop the algorithm; otherwise, we pick v i+1 to be any vertex such that v i → v i+1 .If the algorithm never terminates, then there exist two distinct indices i, j with v i = v j , since G is finite; this implies that G has a cycle, against the hypotheses.Therefore, the procedure must stop in a finite number l of steps, thus yielding a path (v 0 , v 1 , . . ., v l ) from v 0 = v to a sink v l .
For any integer d ≥ 1, let Sym d , Diag d , and Scal d be the sets of d × d real symmetric matrices, real diagonal matrices, and real scalar matrices (i.e.multiples of the d×d identity matrix I d ), respectively.We will write Sym V d for the set of arrays x = (x v ) v∈V , where each x v ∈ Sym d .We will use the notations Diag V d and Scal V d in a similar way.Let us define the 'energy functions' (5.1) where e a denotes the usual exponential of the matrix a.The Golden-Thompson inequality (see e. (5.3)However, the two energy functions are identical only for d = 1.Notice that, by Lemma 5.1, G has at least one sink and one source, possibly coinciding.Throughout, we also assume that there exists at least one vertex of G that is neither a source nor a sink.We can thus fix a subset Γ ⊂ V that contains all the sinks and sources and such that Γ c , the complement of Γ in V , is nonempty.For any set S and any fixed array z = (z v ) v∈Γ ∈ S Γ , let (5.4) Our first result concerns the asymptotic behaviour of the energy functions on Sym Proof.By inequality (5.3), it suffices to prove the claim for χ d .As all norms on a finitedimensional space are equivalent, we may arbitrarily take where ̺(a) denotes the spectral radius of a symmetric matrix a (i.e. the largest absolute value of its eigenvalues).As the spectral radius is a norm on Sym d , it can be easily verified that (5.5) defines a norm on Sym V d .We will show that, for any sequence (x For the sake of notational simplicity, we will drop the superscript of x (n) and leave the dependence on n implicit.
By contradiction, assume that there exists a positive constant C such that, along a subsequence, χ d (x) ≤ C. Since x → ∞, there exists w ∈ Γ c such that, along a further subsequence, ̺(x w ) → ∞.This implies that, passing to a final subsequence, either λ max (x w ) → ∞ or λ max (−x w ) → ∞, where λ max (a) denotes the maximum eigenvalue of a symmetric matrix a.As w ∈ Γ c , it is neither a source nor a sink.By Lemma 5.1, there exists a path (v 0 , v 1 , . . ., v l ) of length l ≥ 1 in G from v 0 = w to a sink v l ∈ Γ.Since G has no cycles, we have v i = v j for all i = j; therefore, all directed edges where we used the bounds tr[e y ] ≥ λ max (e y ) = e λ max (y) for y ∈ Sym (5.7) By iterating (5.7) several times and using (5.6), we obtain By considering a path (u 0 , u 1 , . . ., u m ) of length m ≥ 1 from a source u 0 ∈ Γ to u m = w (which again exists by Lemma 5.1) and using similar bounds, we also have This contradicts the fact that x v l = z v l and x u 0 = z u 0 are both fixed for all x ∈ Sym V d (z), since v l , u 0 ∈ Γ.
Remark 5.3.Above we have assumed that G = (V, E) is acyclic and that Γ is a subset of V containing all the sinks and sources of G.We stress that both hypotheses are necessary for Prop.5.2 to hold.As a counterexample, let G be the cycle graph with n vertices and let Γ = ∅.If a ∈ Sym d and x = (x v ) v∈V is the array with x v = a for all v, then is constant in a; however, for the norm • defined in (5.5), if ̺(a) → ∞, then x → ∞.

Minima of energy functions.
We now study the minima of the functions (5.1)-(5.2) on the set Sym V d (z), where z ∈ Sym Γ d .In words, we wish to minimise the energy functions subject to the constraint that some of the entries of the input array (precisely, those indexed by the vertices of the subset Γ) are fixed.
We start with the simplest case d = 1, in which Sym 1 = Diag 1 = Scal 1 = R and the two energy functions coincide: We denote by ∂ v the partial derivative of a function on R V with respect to the variable x v .
Lemma 5.4.Let z ∈ R Γ .The Hessian matrix of ϕ 1 on R V (z) is positive definite everywhere.In particular, ϕ 1 is strictly convex on R V (z).
Proof.On R V (z) the variables indexed by Γ are fixed to the assigned values z, hence we can consider ϕ 1 and its Hessian as functions of (x v ) v∈Γ c .For v, w ∈ Γ c , we have Thus, the quadratic form of the Hessian of Setting α v := 0 for all v ∈ Γ, it is easy to see that the latter expression equals v,w∈V : v→w Therefore, the Hessian is positive semidefinite everywhere.To prove that it is in fact positive definite, we will show that, if the quadratic form of the Hessian vanishes at α, then α = 0.If the above expression vanishes, then α v = α w for all v, w ∈ V such that v → w.Let v ∈ Γ c .By Lemma 5.1, there exists a path from v to a sink s ∈ Γ.The value α w is then the same for all the vertices w along such a path.We then have α v = α s = 0, since s ∈ Γ.As v ∈ Γ c was arbitrary, it follows that α = (α v ) v∈Γ c = 0.
Proof.By Lemma 5.4, ϕ 1 is a strictly convex function over the convex set R V (z); therefore, it has at most one minimiser.It remains to show the existence of a minimiser.Since ϕ 1 is a continuous function, it admits at least one minimiser on every closed ball B r := {x ∈ R V (z) : x ≤ r}.By Prop.5.2, for r large enough, the minimiser on B r is also a (global) minimiser on R V (z).
The case d > 1 is much more challenging, and we are able to deal with it only under rather strong assumptions on the fixed array z.Nonetheless, this is sufficient for our ultimate purposes.
We will be using the fact that the relation between the eigenvalues and the diagonal entries of a symmetric matrix is completely characterised by the majorisation relation.Let us briefly explain this statement, referring to [HJ13, § 4.3] for proofs and details.For any α = (α 1 , . . ., α d ) ∈ R d , let us denote by α ↓ = (α ↓ 1 , . . ., α ↓ d ) its nonincreasing rearrangement, i.e. the permutation of the coordinates of α such that α (5.8) Theorem 5.6 ([HJ13, Theorem 4.3.45]).Let x ∈ Sym d .Let λ = (λ 1 , . . .λ d ) be the vector of the (real) eigenvalues of x, taken in any order.Let δ i := x(i, i) for 1 ≤ i ≤ d, so that δ = (δ 1 , . . ., δ d ) is the vector of the diagonal entries of x.Then we have λ ≻ δ, and the equality λ ↓ = δ ↓ holds if and only if x is a diagonal matrix.
We now briefly introduce the concept of Schur convexity and state the criterion that is useful for our purposes, referring e.g. to [MOA11, Ch.I.3] for more details.A function In particular, for all α, β such that α ↓ = β ↓ , we have α ≻ β ≻ α, hence H(α) = H(β); in other words, every Schur-convex function is a symmetric function.Additionally, H is called strictly Schur-convex if H(α) > H(β) for all α, β ∈ R d such that α ≻ β and α ↓ = β ↓ .Theorem 5.7 ([MOA11, Ch.I.3, § C]).Let h : R → R and As a consequence of the results just stated, we obtain: Proof.Let λ = (λ 1 , . . ., λ d ) be the vector of the eigenvalues of x, taken in any order.Let δ = (δ 1 , . . ., δ d ) be the vector of (common) diagonal entries of x and y, i.e. δ i = x(i, i) = y(i, i) for all 1 ≤ i ≤ d.Since y is diagonal, notice that the δ i 's are also its eigenvalues.Therefore, the claimed inequality tr[e x ] ≥ tr[e y ] reads as H(λ) ≥ H(δ), where The function H is strictly Schur-convex by Theorem 5.7, since the exponential function is strictly convex.Since λ ≻ δ by Theorem 5.6, we then have H(λ) ≥ H(δ), as required.Moreover, assume that H(λ) = H(δ).Then, by strict Schur-convexity of H, we have λ ↓ = δ ↓ .Again by Theorem 5.6, we conclude that x is diagonal, which in turn implies x = y.
From the latter proposition we deduce the existence and uniqueness of a minimiser of χ d on Sym V d (z), under the assumption that all the 'fixed' entries z are diagonal matrices.
Proof.The claim will immediately follow from the two following facts: (i) for any x ∈ Sym V d (z), there exists y ∈ Diag V d (z) such that χ d (x) ≥ χ d (y), with equality if and only if x = y; (ii) there exists m ∈ Diag V d (z) (as in the statement of the theorem) such that χ d (x) ≥ χ d (m) for any x ∈ Diag V d (z), with equality if and only if x = m.Proof of (i).Fix any x ∈ Sym V d (z).Define y = (y v ) v∈V so that, for all v ∈ V , y v is the diagonal matrix with the same diagonal entries as x, i.e. y v (i, i) = x v (i, i) for 1 ≤ i ≤ d.Since each z v (for v ∈ Γ) is diagonal by hypothesis, we have y ∈ Diag V d (z).For any v, w ∈ V , the matrices x v − x w ∈ Sym d and y v − y w ∈ Diag d have the same diagonal entries, hence tr[e x v −x w ] ≥ tr[e y v −y w ] by Prop.5.8; summing over v → w, we obtain that χ d (x) ≥ χ d (y).Assume now that χ d (x) = χ d (y).Then, tr[e x v −x w ] = tr[e y v −y w ] whenever v → w.Again by Prop.5.8, we then have x v − y v = x w − y w for all v → w.For any v ∈ V , by Lemma 5.1 there exists a path (v 0 , v 1 , . . ., v l ) in G from v 0 = v to a sink v l .Since all sinks are in Γ by assumption (see § 5.1) and both x and y are in Sym V d (z), we have By Prop.5.5, for all i, χ 1 has a unique minimiser m(i, i) on R V (z(i, i)).Therefore, we have and the inequality is strict whenever x = m.
In the case where the 'fixed' entries z are scalar matrices, the inequality (5.3) immediately implies the existence and uniqueness of a minimiser of ϕ d .
Proof.Since Scal Γ d ⊆ Diag Γ d , it follows from Theorem 5.9 that χ d has a unique minimiser m on Sym V d (z), which is of the form specified above.Since the scalar matrices m v and m w commute for any v, w ∈ V , we have e m v e −m w = e m v −m w , hence ϕ d (m) = χ d (m).By (5.3), we then have where the second inequality is strict if x = m.It follows that m is also the unique minimiser of ϕ d on Sym V d (z).(5.9) Take now z = (z v ) v∈Γ such that each z v is a positive multiple of I d , or equivalently log z ∈ Scal Γ d .By Corollary 5.10, Φ d has a unique minimiser m on P V d (z), where log m is the unique minimiser of ϕ d on Sym V d (log z).This implies that, on P V d (z), the Hessian of Φ d at m is positive semidefinite.We now aim to prove the stronger statement that the Hessian of Φ d at m is positive definite.
As in the previous subsections, we first work with d = 1.Recall that, by Lemma 5.4, ϕ 1 is strictly convex on R V (log z), for z ∈ P Γ 1 .The analogous statement does not hold for Φ 1 on P V 1 (z); however, the following is still true: Lemma 5.11.For z ∈ P Γ 1 , the Hessian of Φ 1 on P V 1 (z) is positive definite at any critical point.
Proof.We prove the claim by simply expressing the derivatives of Φ 1 in terms of the derivatives of ϕ 1 .For v ∈ Γ c , the first partial derivative of Φ 1 w.r.t.x v is for any x ∈ P V 1 (z) .
(5.10) Therefore, for v, w ∈ Γ c , Assume now that x is a critical point of Φ 1 on P V 1 (z), i.e. that both sides of (5.10) vanish for all v ∈ Γ c .Then, we have As the Hessian of ϕ 1 on R V (log z) is positive definite everywhere by Lemma 5.4, it follows that the Hessian of Φ 1 on P V 1 (z) is positive definite at x.
For ease of notation, given any v ∈ V and 1 ≤ i ≤ j ≤ d, we will denote by ∂ v;i,j the partial derivative of a function of x ∈ P V d with respect to the real variable x v (i, j).It follows from the definition (5.9) and from the formulas (5.13)-(5.14) that (5.15) We will now compute the second derivatives at any critical point x = (x v ) v∈V ∈ P V d (z).Using (5.12) and (5.15), we have x −1 v (i, k) u∈V : v→u x −1 u (j, ℓ) + u∈V : v→u x −1 u (i, k)x −1 v (j, ℓ) + k ↔ ℓ for v ∈ Γ c , 1 ≤ i ≤ j ≤ d, and 1 ≤ k ≤ ℓ ≤ d.Recall now that the acyclic structure of the underlying graph guarantees that, if v → w, then w → v. Therefore, for v, w ∈ Γ c such that v → w, 1 ≤ i ≤ j ≤ d, and 1 ≤ k ≤ ℓ ≤ d, we have x −1 w (i, k)x −1 w (j, ℓ) + x −1 w (i, ℓ)x −1 w (j, k) .
On the other hand, the second derivative w.r.t.x v (i, j) and x w (k, ℓ) vanishes for all v, w ∈ Γ c such that v → w and w → v.
Here, g ξ (v, w) is an explicit function of ξ, v and w; we stress that it is the same function for all d ≥ 1.It follows from (5.15) that, since x = (ξ v I d ) v∈V is a critical point of Φ d on P V d (z), ξ is a critical point of Φ 1 on P V 1 (ζ).Therefore, the matrix (with 'row index' v and 'column index' w), which is the Hessian matrix of Φ 1 on P V 1 (ζ) at ξ, is positive definite by Lemma 5.11.On the other hand, the matrix f ((i, j), (k, ℓ)) (with 'row index' (i, j) and 'column index' (k, ℓ)) is clearly positive definite as a diagonal matrix with positive diagonal entries.Therefore, the Hessian of Φ d on P V d (z) at x is positive definite, as it can be written as a Kronecker product of two positive definite matrices.
As any minimiser is a critical point, the main result of this section follows immediately from Corollary 5.10 and Lemma 5.12.We start by stating the Laplace approximation integral formula in the multivariate context, which can be found e.g. in [ES00].
Theorem 5.15 ([ES00, Theorem 4.14]).Let A be an open subset of the p-dimensional space R p .Let h : A → C and ̺ : A → R be functions such that (i) A |h(x)| e −k̺(x) dx < ∞ for some k > 0.
(5.17) (iii) h is continuous in a neighbourhood of x 0 and h(x 0 ) = 0.
(iv) ̺ is twice continuously differentiable on A and its Hessian matrix H(x 0 ) at x 0 is positive definite (in particular, its determinant |H(x 0 )| is positive).The extra product in the definition of h is the density of the measure v∈Γ c µ(dx v ) with respect to the Lebesgue measure on P V d (z).Hypothesis (ii) of Theorem 5.15 is satisfied due to Theorem 5.13 and Prop.5.2.Hypotheses (i) and (iii) are matched by the assumptions of Prop.5.14.Finally, hypothesis (iv) also holds because of Theorem 5.13.The asymptotic formula (5.16) then follows from (5.18).
Positive definite matrices.Let P d be the set of all d × d positive definite matrices, i.e. d × d real symmetric matrices with positive eigenvalues.Throughout this article, for x ∈ P d , we denote by |x| the determinant of x and by tr[x] its trace.

Figure 2 .
Figure 2. Graphical representation of a trapezoidal array x ∈ T N,n d as in (2.11), for N = 3 and n = 5.The additional variable s appears in the definition (2.13) of the energy function Φ N,n s (x), in which every summand tr[ab −1 ] corresponds to an arrow pointing from a to b in the figure.

Φ
a graphical representation of the array (2.11) and of the energy function Φ N,n s .For z ∈ P N d , let T N,n d (z) ⊂ T N,n d be the set of all trapezoidal arrays x with n-th row x n = z.For n ≥ N , λ ∈ C n , s ∈ P d and z ∈ P N d , we define ψ N,n λ;s (z) := T if s = 0 d and n = N , ψ N,N λ;0 = ψ N λ corresponds to the Whittaker function defined in (2.4).The absolute convergence of the integral in (2.14), for all λ ∈ C n , can be shown by adapting the proof of [O'C21, Prop.6-(i)].
P N a (z; d z) := P N a (z; z)µ ⊗N (d z) and P N a (z; d z) := P N a ( z; z)µ ⊗N (d z) .
according to the notation of § 2.1.One can show (an analogous computation is made in the proof of Prop.3.6) that, if ℜ(a) > (d−1)/2, , y, z; d z) = Γ d (a) N for any y, y ∈ P N −1 d and z ∈ P N d .(3.11) Using (3.2) and (3.11), we see that, if ℜ is a time-inhomogeneous Markov process with state space T N d and time-n transition kernel Π N α(n)+β .
Let N ≥ 2 and a, b ∈ C. Recalling the definitions (2.8) and (2.15) of the kernels K N b and P N a , respectively, and denoting by δ the Dirac delta kernel on P N d , let us set K N b (z; dy d z) := δ(z; d z) K N b ( z; dy) , (3.15) Λ N a,b (y, z; d y d z) := P N −1 a (y; d y) Q N a+b (y, y, z; d z) , (3.16) for z, z ∈ P N d and y, y ∈ P N −1 d .We then have the following intertwining relation.Proposition 3.6.Let N ≥ 2 and a, b ∈ C such that ℜ(a + b) > (d − 1)/2.Then, .26) Proof.The statement is an application of Theorem B.1, where the state spaces are S = T N d and T = P N d , and the function ϕ : T N d → P N d is the projection ϕ(x) := x N onto the N -th row of x, so that X N (n) = ϕ(X(n)).Hypothesis (i) of Theorem B.1, i.e. the fact that Σ β (z; ϕ −1 {z}) = 1 for any z ∈ P N d , holds because, by definition, the measure Σ β (z; •) is supported on the set T N d (z) of height-N triangular arrays with N -th row equal to z.On the other hand, by Prop.3.5, the time-n transition kernel of X is Π N α(n)+β .Therefore, in this case, hypothesis (ii) of Theorem B.1 reads as the set of intertwining relations

Corollary 4. 4 .
The function Φ N on T N d (I N d ) has a unique global minimiser, at which the Hessian is positive definite.Moreover, each component m i j of the minimiser m = (m i j ) 1≤j≤i≤N is a positive scalar matrix.

Corollary 4. 5 .
We have m 1 1 = I d and m i 1 • • • m i i = 1 for all i = 1, . . ., N .(4.7) Proof.Since m ∈ T N d (I N d ), we have m N j = I d for all j = 1, . . ., N , hence m N 1 • • • m N N = 1.On the other hand, as a minimiser, m is a critical point of Φ N on T N d (I N d ), hence (4.7) follows from Lemma 4.3.Furthermore, since m 1 1 is a multiple of I d with determinant 1, we have m 1 1 = I d .

Corollary 4. 7 .
Let f : T N d → R be a bounded and continuous function and let f k (x) := f (r i j (k)x i j ) 1≤j≤i≤N for k > 0 and x ∈ T N d .(4.11) Assume that f k k→∞ −−−→ f ∞ uniformly on any compact subsets of T N d (I N d ).Then, for any λ, ̺ ∈ R N ,
g. [Bha97]) states that tr[e a e b ] ≥ tr[e a+b ] if a and b are symmetric matrices.It follows that ϕ d (x) ≥ χ d (x) for all x ∈ Sym V d .

Proposition 5. 8 .
Suppose that x ∈ Sym d and y ∈ Diag d have the same diagonal entries.Then tr[e x ] ≥ tr[e y ], and the equality holds if and only if x = y.

Corollary 5. 10 .
Let z = (z v ) v∈Γ ∈ Scal Γ d , so that z v = ζ v I d for all v ∈ Γ and for a certain ζ = (ζ v ) v∈V ∈ R Γ .Then, the function ϕ d admits a unique minimiser on Sym V d (z).Such a minimiser is of the form

5. 3 .
Energy functions in logarithmic variables.It is a well-known fact that the functions Sym d → P d , a → e a and P d → Sym d , a → log a , namely the matrix exponential and the matrix logarithm, are both bijections on the stated domains and inverse to each other.From now on, for any set S, we will use the following compact notations: log x := (log x v ) v∈S ∈ Sym S d for x = (x v ) v∈S ∈ P S d , and e x := (e x v ) v∈S ∈ P S d for x = (x v ) v∈S ∈ Sym S d .Let us consider the analogue of ϕ d 'in logarithmic variables', that is the energy function Φ d (x) := ϕ d (log x) for x ∈ P V d .More explicitly, recalling (5.1), we define Φ d : P V d → R , Φ d (x) := v,w∈V : v→w tr[x v x −1 w ] for all x = (x v ) v∈V ∈ P V d .
v ∈ Γ c and 1 ≤ i ≤ j ≤ d.The critical point equations of Φ d on P V d (z) are then x −1 v u∈V : u→v x u x −1 v = w∈V : v→w x −1 w for all v ∈ Γ c .

Theorem 5. 13 .e
Let z ∈ (Scal + d ) Γ .Then, the function Φ d on P V d (z) has a unique (global) minimiser m, at which the Hessian is positive definite.Moreover, we have m ∈ (Scal + d ) V (z).5.4.Laplace approximation.We will now use Theorem 5.13 to study the asymptotic behaviour of integrals of exponentials of Φ d , via Laplace's approximation method.Recall the definition (1.12) of the measure µ on P d .Proposition 5.14.Let z ∈ (Scal + d ) Γ and let m be the unique global minimiser of Φ d on P V d (z) (see Theorem 5.13).Let g : P V d (z) → C be a continuous function in a neighbourhood of m, with g(m) = 0, and such thatx)| e −kΦ d (x) < ∞ for some k > 0 .v )   g(x) e −kΦ d (x) −kΦ d (m) ,(5.16)where |H(m)| > 0 is the Hessian determinant of Φ d at m and Γ c is the number of vertices in Γ c .
Proof of Prop.5.14.We will apply Theorem 5.15 withA = P V d (z) , h(x) = g(x) v∈Γ c |x v | − d(d+1) 2 , ̺ = Φ d , x 0 = m .The set P V d (z) can be clearly viewed as an open subset of R p , where p = Γ c d(d + 1)/2 is the number of 'free' real variables in A and d is the dimension of each matrix in the array.
Markov dynamics on triangular arrays.Let O d be the real orthogonal group in dimension d.Recall that a random matrix Y in P d is said to be O d -invariant (or orthogonally invariant) if k ⊤ Y k has the same distribution of Y , for every k ∈ O d .n≥1 be a family of independent and O d -invariant random matrices in P d .The GL d -invariant random walk on P d with initial state r ∈ P d and increments .24) It follows from (2.7) and (3.22) that the above kernels are normalised; therefore, they are Markov kernels when the parameters a, λ 1 , . . ., λ N are real.Notice that (3.24) may be seen as a Doob h-transform of the P -kernel (2.15).It is now immediate to deduce a renormalised version of (3.19): 1is a Wishart random matrix with parameter α(n) + β i .If the initial state L(0) of this process is in P N d (respectively, P N d