Zhelobenko - Stern formulas and $B_{n}$ Toda wave functions

A.Galiullin∘, S. Khoroshkin∘∗, M.Lyachko∘
∘National Research University Higher School of Economics, Myasnitskaya 20,
Moscow, 101000, Russia;
∗Institute for Information Transmission Problems RAS (Kharkevich Institute),
Bolshoy Karetny per. 19, Moscow, 127994, Russia

###### Abstract

Using Zhelobenko - Stern formulas for the action of the generators of
orthogonal Lie algebra in corresponding Gelfand - Tsetlin basis, we derive
Mellin - Barnes presentations for the wave functions of $B_{n}$ Toda lattice.
They are in accordance with Iorgov - Shadura formulas.

##  1 Introduction

In the paper [GKL] Gerasimov, Kharchev and Lebedev applied the famous formulas
[GT] for the action of generators of general Lie algebra $gl(n)$ in Gelfand -
Tsetlin basis of irreducible finite - dimensional representations of general
linear group $GL(n,\mathbb{C})$ to obtain Mellin - Barnes presentation of the
wave functions of open $A_{n}$ Toda chain. Using Gelfand -Tsetlin formulas,
they constructed an infinite - dimensional representation of Lie algebra
$gl(n,\mathbb{C})$ in the space of meromorphic functions on $n(n-1)/2$
variables, found there two dual Whittaker vectors and realized, according
Kostant theory [K], the Toda wave function as certain matrix element in this
representation. The same formulas were earlier established by Kharchev and
Lebedev in the technique of Yang-Baxter formalism [KL].

Besides, Toda wave functions admits another presentation, known by the name
Gauss - Givental by means of integrals over spatial variables. It was found
first in [Giv]. Gauss - Givental presentation was then derived in [GLO] for
wave functions of Toda systems related to $B_{n}$, $C_{n}$, $D_{n}$ root
systems.

Mellin transform of these formulas was computed in [KK]. However, the formulas
presented in [KK] differ from that of [GKL] and are not satisfactory by
several reasons. In particular, one cannot find in this presentation Sklyanin
measure and thus the results of [KK] cannot be used to establish the
completeness and orthogonality of the wave function and develop the
corresponding integral transform.

The goal of this paper is to try to fill this gap at least partially for
$B_{n}$ Toda system using representation theoretical tools similar to that of
[GKL]. The only known result in this direction is the paper [IS] of Iorgov and
Shadura where they constructed $B_{n}$ wave function by its decomposition over
related $A_{n}$ Toda wave function. As well as in [KL], this work was done in
a framework of the Yang - Baxter formalism.

Our starting point is an analog of Gelfand - Tsetlin formulas for orthogonal
groups published without a proof by Zhelobenko and Stern in [ZS]. These
formulas look much more complicated compared to [GT] and we did not find
numerous applications of them in the literature. However, after their check we
constructed infinite dimensional representation of the orthogonal Lie algebra
and found there two dual Whittaker vectors. With their help we constructed the
integrals, presenting $B_{n}$ wave functions in which we see all expected
ingredients of Sklyanin measure. The resulting formula can be presented as an
iterative procedure in two ways.

Firstly, it is an iterative procedure over the rank of orthogonal group and
this is probably the most interesting result of this paper. Each step can be
interpreted as an action of the raising integral operator, which kernel is
itself an integral over intermediate variables. Such type of structure we also
observe in Gauss - Givental representation [GLO, (1.74)]. Second, we can
consider two successive iterative integrals combining them in other parity.
Then intermediate step becomes precisely a degeneration of $B_{n}$ Gustafson
integral and can be explicitly evaluated. In this way we arrive at Iorgov -
Shadura formula.

Note two subtle points of our construction. First, Zhelobenko-Stern formulas
are written for the generators of orthogonal Lie algebras in their orthogonal
realization, while Whittakker vectors refer to simple root generators. An
existence of Whittakker vectors in a factorized form was not evident from the
beginning. By the same reasons action of the Cartan subalgebra in
corresponding infinite dimensional representation cannot be written, contrary
to $gl(n)$, in terms of multiplications by linear functions. Fortunately, it
is so for the action on Whittaker vectors.

Despite the fact that Zhelobenko -Stern formulas are written uniformly for all
orthogonal Lie algebras, we managed to find Whittaker vectors only for Lie
algebras $so(2n-1)$. More precisely, we managed firstly to found factorized
expressions for ’degenerate’ Whittaker vectors for $so(2n)$, for which one of
the simple generators act by zero, so that they are essentially Whittaker
vectors for embedded $gl(n)$ Lie algebra. Restricting these vectors to
$so(2n-1)$ we get ’nondegenerate’ Whittaker vectors for this Lie algebra.

##  2 Gelfand-Tsetlin type representation

###  2.1 Zhelobenko-Stern formulas

It is well known that each irreducible representation of the orthogonal groups
$SO(2n+1)$ and $SO(2n)$ is parametrized by its signature, given by ordered
sequences of integers, respectively

$\begin{split}p_{1}\geq p_{2}\geq\ldots p_{n-1}\geq p_{n}\geq 0\\\ p_{1}\geq
p_{2}\geq\ldots p_{n-1}\geq|p_{n}|\end{split}$ (2.1)

and the restriction of irreducible representation of $SO(2n+1)$ to $SO(2n)$
has simple spectrum described by all signatures $q_{1},\ldots q_{n}$,
satisfying interleaving inequalities

$p_{1}\geq q_{1}\geq p_{2}\geq\ldots q_{n-1}\geq p_{n}\geq q_{n}\geq-p_{n}.$
(2.2)

Analogously, the restriction of irreducible representation of $SO(2n)$ to
$SO(2n-1)$ has simple spectrum described by all signatures $q_{1},\ldots
q_{n}$, satisfying interleaving inequalities

$p_{1}\geq q_{1}\geq p_{2}\geq\ldots q_{n-1}\geq p_{n}.$ (2.3)

This enables one to construct an orthogonal basis of irreducible
representation of the orthogonal group $SO(n)$ parametrised by Gelfand-Tsetlin
tableaux

$\mathbf{p}=\left(\begin{array}[]{cccc}p_{n-1,1}&p_{n-1,2}&\ldots&p_{n-1,[\frac{n}{2}]}\\\
\vdots&\vdots&\rotatebox{45.0}{\ldots}\\\ p_{3,1}&p_{3,2}\\\ p_{2,1}\\\
p_{1,1}\end{array}\right)=\left(\begin{array}[]{c}\mathbf{p}_{n-1}\\\
\vdots\\\ \mathbf{p}_{3}\\\ \mathbf{p}_{2}\\\
\mathbf{p}_{1}\end{array}\right)$ (2.4)

The upper row $\mathbf{p}_{n-1}$ indicates the signature of the irreducible
representation of $SO(n)$ and is fixed for all tableaux parametrising its
basic vectors, the second row indicates the signature of the restriction to
$SO(n-1)$, etc., the integer $p_{11}$ indicates the irreducible of $SO(2)$.
Integers $p_{ij}$ should satisfy the row by row interleaving inequalities
(2.1)–(2.2).

It is natural to shift the signatures by the corresponding half sum of
positive roots of the related root system, that is, we set

$\begin{split}m_{2n,k}&=p_{2n,k}+(n-k)+\frac{1}{2},\\\
m_{2n-1,k}&=p_{2n-1,k}+(n-k)\end{split}$ (2.5)

Zhelobenko and Stern [ZS, Chapter II, Section 5.8] presented without a proof a
precise expression for the matrix elements of generators of the Lie algebra
$so(n)$ in the corresponding orthogonal basis.

The Lie algebra $so(n)$ is generated, as a vector space, by elements

$I_{ij}=e_{ij}-e_{ji},\qquad i>j$

As a Lie algebra, it is generated by elements

$I_{i+1,i},\qquad i=1,\ldots,n-1$

with defining relations

$\begin{split}&\left[I_{i+1,i},\left[I_{i+2,i+1},I_{i+1,i}\right]\right]=I_{i+2,i+1}\quad
i=1,\ldots,n-2\\\
&\left[I_{i+2,i+1},\left[I_{i+1,i},I_{i+2,i+1}\right]\right]=I_{i+1,i}\quad
i=1,\ldots,n-2\\\ &\left[I_{i+1,i},I_{j+1,j}\right]=0\quad|i-j|>1\end{split}$
(2.6)

After a renormalization, eliminating square roots in the coefficients and
correcting misprints, their formulas look like

$\begin{split}I_{2k+1,2k}&=-\sum_{j=1}^{k}\frac{\prod\limits_{r=1}^{k-1}(m_{2k-2,r}+m_{2k-1,j}+\frac{1}{2})\prod\limits_{r=1}^{k}(m_{2k-1,j}-m_{2k,r}+\frac{1}{2})}{2\prod\limits_{r\neq
j}(m_{2k-1,j}-m_{2k-1,r})(m_{2k-1,j}+m_{2k-1,r}+1)}e^{\partial_{m_{2k-1,j}}}-\\\
&-\sum_{j=1}^{k}\frac{\prod\limits_{r=1}^{k-1}(m_{2k-1,j}-m_{2k-2,r}-\frac{1}{2})\prod\limits_{r=1}^{k}(m_{2k,r}+m_{2k-1,j}-\frac{1}{2})}{2\prod\limits_{r\neq
j}(m_{2k-1,j}-m_{2k-1,r})(m_{2k-1,j}+m_{2k-1,r}-1)}e^{-\partial_{m_{2k-1,j}}}\end{split}$
(2.7)
$\begin{split}I_{2k+2,2k+1}&=\sum_{j=1}^{k}\frac{\prod\limits_{r=1}^{k+1}((m_{2k,j}+\frac{1}{2})^{2}-m_{2k+1,r}^{2})}{2m_{2k,j}(m_{2k,j}+\frac{1}{2})\prod\limits_{r\neq
j}(m^{2}_{2k,j}-m^{2}_{2k,r})}e^{\partial_{m_{2k,j}}}+\\\
&+\sum_{j=1}^{k}\frac{\prod\limits_{r=1}^{k}((m_{2k,j}-\frac{1}{2})^{2}-m_{2k-1,r}^{2})}{2m_{2k,j}(m_{2k,j}-\frac{1}{2})\prod\limits_{r\neq
j}(m^{2}_{2k,j}-m^{2}_{2k,r})}e^{-\partial_{m_{2k,j}}}+\\\
&+i\frac{\prod\limits_{r=1}^{k}m_{2k-1,r}\prod\limits_{r=1}^{k+1}m_{2k+1,r}}{\prod\limits_{r=1}^{k}(m_{2k,r}+\frac{1}{2})(m_{2k,r}-\frac{1}{2})}\end{split}$
(2.8)

Here the operators $e^{\pm\partial m_{kj}}$ are operators of shifts of the
entries of Gelfand -Tsetlin tableau: the operator $e^{\pm\partial m_{kj}}$
changes $m_{kj}$ by $m_{kj}\pm 1$ (and respectively $p_{kj}$ by $p_{kj}\pm
1$). We can extend the RHS of relations (2.7) and (2.8) to arbitrary complex
parameters $m_{ij}$ and regard them as operators acting in the space of
rational functions on $m_{ij}$.

###  2.2 Representation in meromorphic functions

Following [GKL], we renormalize the variables

$m_{ij}=\frac{\nu_{ij}}{ic}$ (2.9)

in order to have an additional scaling variable in the representation. Then we
have

$\begin{split}I_{2k+1,2k}&=-\frac{1}{ic}\sum_{j=1}^{k}\frac{\prod\limits_{r=1}^{k-1}(\nu_{2k-2,r}+\nu_{2k-1,j}+\frac{{ic}}{2})\prod\limits_{r=1}^{k}(\nu_{2k-1,j}-\nu_{2k,r}+\frac{{ic}}{2})}{2\prod\limits_{r\neq
j}(\nu_{2k-1,j}-\nu_{2k-1,r})(\nu_{2k-1,j}+\nu_{2k-1,r}+{ic})}e^{{ic}\partial_{\nu_{2k-1,j}}}-\\\
&-\frac{1}{{ic}}\sum_{j=1}^{k}\frac{\prod\limits_{r=1}^{k-1}(\nu_{2k-1,j}-\nu_{2k-2,r}-\frac{ic}{2})\prod\limits_{r=1}^{k}(\nu_{2k,r}+\nu_{2k-1,j}-\frac{{ic}}{2})}{2\prod\limits_{r\neq
j}(\nu_{2k-1,j}-\nu_{2k-1,r})(\nu_{2k-1,j}+\nu_{2k-1,r}-{ic})}e^{-{ic}\partial_{\nu_{2k-1,j}}}\end{split}$
(2.10)
$\begin{split}iI_{2k+2,2k+1}&=\frac{1}{c}\sum_{j=1}^{k}\frac{\prod\limits_{r=1}^{k+1}((\nu_{2k,j}+\frac{{ic}}{2})^{2}-\nu_{2k+1,r}^{2})}{2\nu_{2k,j}(\nu_{2k,j}+\frac{{ic}}{2})\prod\limits_{r\neq
j}(\nu^{2}_{2k,j}-\nu^{2}_{2k,r})}e^{{ic}\partial_{\nu_{2k,j}}}+\\\
&+\frac{1}{c}\sum_{j=1}^{k}\frac{\prod\limits_{r=1}^{k}((\nu_{2k,j}-\frac{{ic}}{2})^{2}-\nu_{2k-1,r}^{2})}{2\nu_{2k,j}(\nu_{2k,j}-\frac{{ic}}{2})\prod\limits_{r\neq
j}(\nu^{2}_{2k,j}-\nu^{2}_{2k,r})}e^{-{ic}\partial_{\nu_{2k,j}}}+\\\
&-\frac{1}{ic}\frac{\prod\limits_{r=1}^{k}\nu_{2k-1,r}\prod\limits_{r=1}^{k+1}\nu_{2k+1,r}}{\prod\limits_{r=1}^{k}(\nu_{2k,r}+\frac{{ic}}{2})(\nu_{2k,r}-\frac{{ic}}{2})}\end{split}$
(2.11)

######  Proposition 1.

The operators (2.7) and (2.8) satisfy the defining relations (2.6) of the
generators of orthogonal Lie algebra $so(n)$

Surely, this statement follows from its validity in finite dimensional
representations, since the relations are then satisfied on sufficiently many
integer points. However, since the proof of the formulas is missing in [ZS],
we checked the defining relation (2.6) directly.

######  Proposition 2.

The center of $SO(2n+1)$ acts by multiplication on symmetric polynomials in
$\nu_{2n,k}^{2}$. The center of $SO(2n)$ acts by multiplication of polynomials
on $\nu^{2}_{2n-1,k}$, symmetric with respect to the permutations of the
variables, and by powers of the monomial
$\nu_{2n-1,1}\nu_{2n-,2}\cdots\nu_{2n-1,n}$.

This follows from Harish-Chandra isomorphism, see e.g. [D, Section 7.4].

Define the following automorphism of the space of meromorphic functions on
$\nu_{kl}$ and $c$:

$\tau(\nu_{2k+1,j})=\nu_{2k+1,j},\qquad\tau(\nu_{2k,j})=-\nu_{2k,j},\qquad\tau(c)=-c$
(2.12)

######  Lemma 1.

We have the relations

$\tau I_{k,k+1}=-I_{k,k+1}\tau$ (2.13)

##  3 Whittaker vectors

###  3.1 Two chains of groups

Zhelobenko - Stern construction of the Gelfand - Tsetlin basis for orthogonal
groups uses the chain of embeddings

$\operatorname{{\mathbf{i}}}_{n}:SO(N)\hookrightarrow SO(N+1)$ (3.1)

where the compact group $SO(N)$ is embedded into the compact group $SO(N+1)$
as the stabilizer of the vector $e_{N+1}$ so that the generators $I_{ij}$,
$i,j\leq N$ of the Lie algebra $so(N)$ are identified with the corresponding
generators of the Lie algebra $so(N+1)$.

However, for the construction of Whittaker vectors in the related infinite
dimensional representations of $so(N)$ in meromorphic functions we pass to
another, noncompact real form $SO(N,J)$ of the group $SO(N,\mathbb{C})$ and
use the chain of the corresponding Lie algebras compatible with the natural
chain of Lie group $SO(N)$. Here

$J=\begin{pmatrix}0&0&...&0&1\\\ 0&0&...&1&0\\\ &&...&&\\\ 0&1&...&0&0\\\
1&0&...&0&0\end{pmatrix}$

The Lie algebra $so(N)$ is generated by elements $I_{ij}$, with the relation
$I_{ji}=-I_{ij}$, so that the elements $I_{ij}$ with $i>j$ are chosen as a
linear basis of Lie algebra $so(N)$. The Lie algebra $so(N,J)$, acting in the
space with the basis $f_{1},\ldots,f_{N}$ is generated by the elements

$F_{ij}=f_{ij}-f_{\hat{j}\hat{i}},\qquad\text{where}\qquad
f_{ij}(f_{k})=\delta_{jk}f_{i}\qquad\text{and}\qquad\hat{i}=N+1-i$

with the relation $F_{\hat{j},\hat{i}}=-F_{i,j}$ so that the elements $F_{ij}$
with $i+j\leq N$ can be chosen as a linear basis of the Lie algebra $so(N,J)$.
The elements $F_{ij}$ for $i<j$ form a positive nilpotent subalgebra, elements
$F_{kk}$ form a Cartan subalgebra.

The chain of embedding

$\operatorname{{\mathbf{j}}}_{n}:SO(N,J)\hookrightarrow SO(N+1,J)$ (3.2)

is different. The group $SO(2n,J)$ is the stabilizer of the vector $f_{n+1}$
in the group $SO(2n+1,J)$, while the group $SO(2n-1,J)$ is the stabilizer of
the element $f_{n-1}+(-1)^{n-1}f_{n}$ in $SO(2n,J)$. Let us describe the maps

$\operatorname{{\mathbf{s}}}_{N}:so(N)\to so(N,J)$ (3.3)

of complex Lie algebras, which intertwine the embeddings (3.1) and (3.2). On
the level of bases of the vector space $\mathbb{C}^{2n}$ the map
$\operatorname{{\mathbf{s}}}_{2n}$ corresponds to the transformation of
initial orthogonal basis $e_{1},\ldots,e_{2n}$ of $\mathbb{C}^{2n}$ to the
defining basis $f_{1},\ldots,f_{2n}$ of the form $J$,

$(f_{1},\ldots,f_{n},f_{n+1},\ldots,f_{2n}),\qquad(f_{i},f_{\hat{j}})=\delta_{ij}$

by the relation

$\begin{split}f_{j}=i^{j}\cdot\frac{ie_{2j-1}+e_{2j}}{\sqrt{2}},\qquad
f_{\hat{j}}=i^{-j}\cdot\frac{-ie_{2j-1}+e_{2j}}{\sqrt{2}},\qquad
j=1,,\ldots,n\end{split}$ (3.4)

For the group $SO(2n+1)$ we transform the initial orthogonal basis
$e_{1},\ldots,e_{2n+1}$ of $\mathbb{C}^{2n+1}$ to the defining basis
$f_{1},\ldots,f_{2n+1}$ of the form $J$,

$(f_{1},\ldots,f_{n},f_{n+1},\ldots,f_{2n+1}),\qquad(f_{i},f_{\hat{j}})=\delta_{ij}$

by the relation

$\begin{split}f_{j}=i^{j}\cdot\frac{ie_{2j-1}+e_{2j}}{\sqrt{2}},\qquad
f_{\hat{j}}=i^{-j}\cdot\frac{-ie_{2j-1}+e_{2j}}{\sqrt{2}},\qquad
j=1,,\ldots,n,\qquad f_{n+1}=e_{n+1}\end{split}$

Correspondingly, the transformation formula from Lie algebra elements $I_{ij}$
to $F_{ij}$ are given by conjugation of the matrix $(I_{ij})$ be means of the
corresponding transition matrix. In particular, we have the following
expressions for the generators of Lie algebra $so(2n,J)$:

$\displaystyle F_{j,j+1}=$
$\displaystyle\frac{1}{2}(I_{2j+1,2j}-I_{2j+2,2j-1})+\frac{i}{2}(I_{2j+2,2j}+I_{2j+1,2j-1}),$
$\displaystyle F_{j+1,j}=$
$\displaystyle\frac{1}{2}(I_{2j+2,2j-1}-I_{2j+1,2j})+\frac{i}{2}(I_{2j+2,2j}+I_{2j+1,2j-1},)$
$\displaystyle F_{j,2n-j}=$
$\displaystyle(-1)^{j}\Big{(}\frac{1}{2}(-I_{2j+1,2j}-I_{2j+2,2j-1})+\frac{i}{2}(I_{2j+2,2j}-I_{2j+1,2j-1})\Big{)},$
(3.5) $\displaystyle F_{2n-j,j}=$
$\displaystyle(-1)^{j}\Big{(}\frac{1}{2}(I_{2j+2,2j-1}+I_{2j+1,2j})+\frac{i}{2}(I_{2j+2,2j}-I_{2j+1,2j-1})\Big{)},$
$\displaystyle F_{j,j}=-iI_{2j,2j-1}$ (3.6)

Here $j=1,\ldots,n-1$. Besides, instead of the use of the last simple root
generator $F_{n-1,n}^{2n-1}$ of the Lie algebra $so(2n-1,J)$ it is convenient
to use its image

$\operatorname{{\mathbf{j}}}_{2n-1}(F_{n-1,n}^{2n-1})=\frac{1}{\sqrt{2}}\big{(}F_{n-1,n}^{2n}-(-1)^{n-1}F_{n-1,n+1}^{2n}\big{)}$
$\operatorname{{\mathbf{j}}}_{2n-1}(F_{n,n-1}^{2n-1})=\frac{1}{\sqrt{2}}\big{(}F_{n,n-1}^{2n}-(-1)^{n-1}F_{n+1,n-1}^{2n}\big{)}$

in the Lie algebra $so(2n)$.

Note that for the Lie algebra $so(n)$ the automorphism $\tau$ represents the
longest element of the Weyl group,

$\tau F_{j,j+1}=F_{j+1,j}\tau.$ (3.7)

###  3.2 Right and left Whittaker vectors

For further convenience we denote by $\bm{\nu}_{n}$ the tuples of variables

$\bm{\nu}_{2k}=\\{\nu_{2k,1},\ldots\nu_{2k,k}\\},\qquad\bm{\nu}_{2k-1}=\\{\nu_{2k-1,1},\ldots\nu_{2k-1,k}\\}$
(3.8)

and by $\hat{\bm{\nu}}_{n}$ the Gelfand - Tsetlin array

$\widehat{\bm{\nu}}_{n}=\begin{pmatrix}\bm{\nu}_{n}\\\ \bm{\nu}_{n-1}\\\
\vdots\\\ \bm{\nu}_{2}\\\ \bm{\nu}_{1}\end{pmatrix}$ (3.9)

For any two sets

$\mathbf{X}=(x_{1},\ldots,x_{n}),\qquad\mathbf{Y}=(y_{1},\ldots,y_{m})$

set

$s(\mathbf{X};\mathbf{Y})=\prod\limits_{x\in\mathbf{X},\,y\in\mathbf{Y}}(ic)^{\frac{x+y}{ic}}\cdot\Gamma\left(\frac{x+y}{ic}+\frac{1}{2}\right)$
(3.10)

With their help for any $k>0$ define the meromorphic functions $W_{k}^{\pm}$
and $V_{k}^{\pm}$ by the relation

$\begin{split}W_{k}^{\pm}&=e^{\mp\frac{\pi}{2c}\sum_{j}\nu_{2k,j}}s(\pm\bm{\nu}_{2k-1},\bm{\nu}_{2k})s(-\bm{\nu}_{2k},\pm\bm{\nu}_{2k+1}),\\\
V_{k}^{\pm}&=e^{\mp\frac{\pi}{2c}\sum_{j}\nu_{2k,j}}s(\pm\bm{\nu}_{2k-1},\bm{\nu}_{2k}).\\\
\end{split}$ (3.11)

Set

$\begin{split}w_{n}&=e^{-\frac{\pi n}{2c}\delta_{2n}}W_{1}^{+}W_{2}^{-}\cdots
W_{n-1}^{(-)^{n}}V_{n}^{(-)^{n+1}},\\\
w_{n}^{\prime}&=e^{\frac{\pi}{c}\sum\limits_{k=1}^{n}\delta_{2k-1}}\tau(w_{n})=e^{-\frac{\pi
n}{2c}\delta_{2n}}W_{1}^{-}W_{2}^{+}\cdots
W_{n-1}^{(-)^{n-1}}V_{n}^{(-)^{n}}\end{split}$ (3.12)

Here

$\delta_{k}=\sum_{j}\gamma_{kj},\qquad\gamma_{kj}\in{\bm{\gamma}}_{k}.$

######  Theorem 1.

The functions $w_{n}$ and $w^{\prime}_{n}$ are left and right $SO(2n+1)$
Whittaker vectors:

$\displaystyle F_{k,k+1}w_{n}=\frac{(-1)^{k+1}}{ic}w_{n},\qquad k<n,$
$\displaystyle F_{n,n+1}w_{n}=\frac{(-1)^{n+1}}{\sqrt{2}ic}w_{n},$ (3.13)
$\displaystyle
F_{k+1,k}w^{\prime}_{n}=\frac{(-1)^{k+1}}{ic}w^{\prime}_{n},\qquad k<n,$
$\displaystyle
F_{n+1,n}w^{\prime}_{n}=\frac{(-1)^{n+1}}{\sqrt{2}ic}w^{\prime}_{n}.$ (3.14)

Proof is given in Appendix A

Precise formulas for Whittaker vector and Whittaker function look better after
the following change of variables

$\gamma_{2k-1,j}=(-1)^{k+1}\nu_{2k-1,j},\qquad\gamma_{2k,j}=\nu_{2k,j}.$
(3.15)

Then the Whittaker vectors can be written as

$\begin{split}&w_{n}=e^{-\frac{\pi
n}{2c}\delta_{2n}}e^{\sum\limits^{n}_{k=1}(-1)^{k}\delta_{2k}}\prod\limits^{n}_{k=1}s({\bm{\gamma}}_{2k-1},{\bm{\gamma}}_{2k})\prod\limits^{n-1}_{k=1}s(-{\bm{\gamma}}_{2k},-{\bm{\gamma}}_{2k+1})\\\
&w^{\prime}_{n}=e^{-\frac{\pi
n}{2c}\delta_{2n}}e^{\frac{\pi}{2c}\sum\limits^{n}_{k=1}(-1)^{k+1}\delta_{2k}}\prod\limits^{n}_{k=1}s(-{\bm{\gamma}}_{2k-1},{\bm{\gamma}}_{2k})\prod\limits^{n-1}_{k=1}s(-{\bm{\gamma}}_{2k},{\bm{\gamma}}_{2k+1})\end{split}$
(3.16)

###  3.3 Action of Cartan subalgebra

In Gelfand - Tsetlin representation of $gl_{n}$ see [GKL], the Cartan
subalgebra acts by multiplication by linear functions on the variables
$\gamma_{ij}$. It is not so for Gelfand-Tsetlin representations of $so(n)$,
which we study here. However, the Cartan subalgebra of $so(n)$ acts in a
similar way on Whittaker vectors. The Cartan subalgebra of $so(n)$ is
generated by the elements $F_{kk}=-iI_{2k,2k-1}$.

######  Proposition 3.

$\displaystyle F_{kk}w_{n}$
$\displaystyle=\frac{(-1)^{k-1}}{ic}\left(\sum\limits_{j=1}^{k}\nu_{2k-1,j}+\sum\limits_{j=1}^{k-1}\nu_{2k-3,j}+(-1)^{k}(k-1)ic\right)w_{n}$
(3.17) $\displaystyle F_{kk}w^{\prime}_{n}$
$\displaystyle=\frac{(-1)^{k-1}}{ic}\left(\sum\limits_{j=1}^{k}\nu_{2k-1,j}+\sum\limits_{j=1}^{k-1}\nu_{2k-3,j}-(-1)^{k}i(k-1)c\right)w^{\prime}_{n}$
(3.18)

In the variables $\gamma_{kj}$ the relations (3.17) look as follows

$\displaystyle F_{kk}w_{n}$
$\displaystyle=\frac{1}{ic}\left(\sum\limits_{j=1}^{k}\gamma_{2k-1,j}-\sum\limits_{j=1}^{k-1}\gamma_{2k-3,j}-(k-1)ic\right)w_{n}$
(3.19) $\displaystyle F_{kk}w^{\prime}_{n}$
$\displaystyle=\frac{1}{ic}\left(\sum\limits_{j=1}^{k}\gamma_{2k-1,j}-\sum\limits_{j=1}^{k-1}\gamma_{2k-3,j}+(k-1)ic\right)w^{\prime}_{n}$
(3.20)

Proof is given in Appendix B

##  4 Whittaker wave function

We construct Toda wave function as a matrix coefficient (4.4) of elements of
Cartan subgroup with respect to two Whittaker vectors.

###  4.1 Invariant pairing

The invariant pairing looks the same in both $\nu$ and $\gamma$ variables.

Define the functions $\tilde{\mu}({\bm{\gamma}}_{2k})$ and
$\tilde{\mu}({\bm{\gamma}}_{2k+1})$ by the relations

$\begin{split}\tilde{\mu}({\bm{\gamma}}_{2k})=e^{\frac{\pi}{c}\sum_{j}\gamma_{2k,j}}&\prod_{r}\Gamma^{-1}\Big{(}\frac{2\gamma_{2k,r}}{ic}\Big{)}\Gamma^{-1}\Big{(}\frac{-2\gamma_{2k,r}}{ic}\Big{)}\times\\\
&\prod_{r\not=s}\Gamma^{-1}\Big{(}\frac{\gamma_{2k,r}-\gamma_{2k,s}}{ic}\Big{)}\prod_{r<s}\Gamma^{-1}\Big{(}\frac{\gamma_{2k,r}+\gamma_{2k,s}}{ic}\Big{)}\Gamma^{-1}\Big{(}\frac{-\gamma_{2k,r}-\gamma_{2k,s}}{ic}\Big{)},\\\
\tilde{\mu}({\bm{\gamma}}_{2k+1})=&\prod_{r\not=s}\Gamma^{-1}\Big{(}\frac{\gamma_{2k+1,r}-\gamma_{2k+1,s}}{ic}\Big{)}\\\
&\prod_{r<s}\Gamma^{-1}\Big{(}1+\\!\frac{\gamma_{2k+1,r}+\gamma_{2k+1,s}}{ic}\Big{)}\Gamma^{-1}\Big{(}1-\\!\frac{\gamma_{2k+1,r}+\gamma_{2k+1,s}}{ic}\Big{)},\end{split}$
(4.1)

that is

$\begin{split}\tilde{\mu}({\bm{\gamma}}_{2k})=e^{\frac{\pi}{c}\delta_{2k}}&\prod_{r<s}\left|\Gamma\Big{(}\frac{\gamma_{2k,r}-\gamma_{2k,s}}{ic}\Big{)}\right|^{-2}\left|\Gamma\Big{(}\frac{\gamma_{2k,r}+\gamma_{2k,s}}{ic}\big{)}\right|^{-2}\prod_{r}\left|\Gamma\Big{(}\frac{2\gamma_{2k,r}}{ic}\Big{)}\right|^{-2},\\\
\tilde{\mu}({\bm{\gamma}}_{2k+1})=&\prod_{r<s}\left|\Gamma\Big{(}\frac{\gamma_{2k+1,r}-\gamma_{2k+1,s}}{ic}\Big{)}\right|^{-2}\left|\Gamma\Big{(}1+\\!\frac{\gamma_{2k+1,r}+\gamma_{2k+1,s}}{ic}\Big{)}\right|^{-2}\,\end{split}$
(4.2)

and define a scalar product on functions in $M$ as

$(f,g)_{2n}=\int_{\mathbb{R}^{n^{2}}}\bar{f}(\widehat{{\bm{\gamma}}}_{2n-1})g(\widehat{{\bm{\gamma}}}_{2n-1})\tilde{\mu}(\widehat{{\bm{\gamma}}}_{2n-1})d\widehat{{\bm{\gamma}}}_{2n-1}$
(4.3)

where

$\tilde{\mu}(\widehat{{\bm{\gamma}}}_{n})=\prod_{k=1}^{n}\tilde{\mu}({\bm{\gamma}}_{k})$

Then

######  Lemma 2.

The operators $I_{2k+1,2k}$ and $iI_{2k+2,2k+1}$ are skew symmetric with
respect to the pairing (4.3). In particular, operators $F_{kl}$ are skew
symmetric with respect to the pairing (4.3).

###  4.2 Integral formula

Due to (3.16) the product $\bar{w^{\prime}}_{n}w_{n}$ looks as

$\begin{split}\bar{w^{\prime}}_{n}w_{n}=&e^{-\frac{\pi
n}{c}\delta_{2n}}\prod\limits^{n-1}_{k=1}s({\bm{\gamma}}_{2k-1},{\bm{\gamma}}_{2k})\bar{s}(-{\bm{\gamma}}_{2k-1},{\bm{\gamma}}_{2k})s(-{\bm{\gamma}}_{2k},-{\bm{\gamma}}_{2k+1})\bar{s}(-{\bm{\gamma}}_{2k},{\bm{\gamma}}_{2k+1})\times\\\
&s({\bm{\gamma}}_{2n-1},{\bm{\gamma}}_{2n})\bar{s}(-{\bm{\gamma}}_{2n-1},{\bm{\gamma}}_{2n})\end{split}$

or, in terms of Gamma functions

$\begin{split}&w_{n}\cdot\overline{w_{n}^{\prime}}=e^{-\frac{\pi}{c}\sum\limits^{n-1}_{k=1}\delta_{2k}}\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}c^{\frac{1}{ic}\left(\sum\limits^{n-1}_{k=1}2k(\delta_{2k-1}-\delta_{2k+1})+2n\delta_{2n-1}\right)}\\\
&\times\left(\prod^{n-1}_{k=1}\Gamma\left(\frac{\pm{\bm{\gamma}}_{2k}+{\bm{\gamma}}_{2k-1}}{ic}+\frac{1}{2}\right)\Gamma\left(\frac{\pm{\bm{\gamma}}_{2k}-{\bm{\gamma}}_{2k+1}}{ic}+\frac{1}{2}\right)\right)\cdot\Gamma\left(\frac{\pm{\bm{\gamma}}_{2n}+{\bm{\gamma}}_{2n-1}}{ic}+\frac{1}{2}\right)\end{split}$

Thus the Whittaker function

$\Psi_{{\bm{\gamma}}_{2n}}=(w_{n}^{\prime},\,e^{-\sum_{k=1}^{n}x_{k}F_{kk}}w_{n})_{2n}$
(4.4)

is given by the integral

$\begin{split}\Psi_{{\bm{\gamma}}_{2n}}(\mathbf{x}_{n})=&\int_{\mathbb{R}^{n^{2}}}\prod_{k=1}^{2n-1}{\mu}({\bm{\gamma}}_{k})d{\bm{\gamma}}_{k}\,\cdot
c^{\frac{2}{ic}\sum_{k=1}^{n}\delta_{2k-1}}e^{\frac{1}{ic}\sum_{k=1}^{n}(\delta_{2k-3}-\delta_{2k-1}+(k-1)ic)x_{k}}\times\\\
&\prod_{k=1}^{n-1}\prod_{i,j=1}^{k}\prod_{l=1}^{k+1}\Gamma\left(\frac{\pm\gamma_{2k,i}+\gamma_{2k-1,j}}{ic}+\frac{1}{2}\right)\Gamma\left(\frac{\pm\gamma_{2k,i}-\gamma_{2k+1,l}}{ic}+\frac{1}{2}\right)\times\\\
&\prod_{i,j=1}^{n}\Gamma\left(\frac{\pm\gamma_{2n,i}+\gamma_{2n-1,j}}{ic}+\frac{1}{2}\right)\end{split}$
(4.5)

Here

$\begin{split}{\mu}({\bm{\gamma}}_{2k})=e^{-\frac{\pi}{c}\delta_{2k}}\tilde{\mu}({\bm{\gamma}}_{2k})&=\prod_{r<s}\left|\Gamma\Big{(}\frac{\gamma_{2k,r}-\gamma_{2k,s}}{ic}\Big{)}\right|^{-2}\left|\Gamma\Big{(}\frac{\gamma_{2k,r}+\gamma_{2k,s}}{ic}\big{)}\right|^{-2}\prod_{r}\left|\Gamma\Big{(}\frac{2\gamma_{2k,r}}{ic}\Big{)}\right|^{-2},\\\
{\mu}({\bm{\gamma}}_{2k+1})=\tilde{\mu}({\bm{\gamma}}_{2k+1})&=\prod_{r<s}\left|\Gamma\Big{(}\frac{\gamma_{2k+1,r}-\gamma_{2k+1,s}}{ic}\Big{)}\right|^{-2}\left|\Gamma\Big{(}1+\\!\frac{\gamma_{2k+1,r}+\gamma_{2k+1,s}}{ic}\Big{)}\right|^{-2}\,\end{split}$
(4.6)

The measure functions ${\mu}({\bm{\gamma}}_{j})$ do not contain exponential
factors.

The convergence of the integral (4.5) can be proved by the arguments given in
[IS, Appendix A]. Namely, let us complete the sequence

${\bm{\gamma}}_{1}=\\{\gamma_{11}\\},{\bm{\gamma}}_{2}=\\{\gamma_{21}\\},{\bm{\gamma}}_{3}=\\{\gamma_{31},\gamma_{32}\\},\ldots{\bm{\gamma}}_{2n}=\\{\gamma_{2n,1},\ldots\gamma_{2n,n}\\},$

to the sequence

${\bm{\gamma}}^{\prime}_{1}=\\{\gamma^{\prime}_{1,-1},\gamma^{\prime}_{11}\\},{\bm{\gamma}}^{\prime}_{2}=\\{\gamma^{\prime}_{2,-1},\gamma^{\prime}_{2,0},\gamma^{\prime}_{21}\\},,\ldots{\bm{\gamma}}^{\prime}_{2n}=\\{\gamma_{2n,-n},\ldots\gamma_{2n,n}\\},$

where

$\gamma^{\prime}_{n,k}=\gamma_{n,k},\
\gamma^{\prime}_{n,-k}=-\gamma_{n,k}\qquad\text{for}\qquad
k>0,\qquad\text{and}\qquad\gamma^{\prime}_{2m,0}=0$

as is customary in the representation theory of orthogonal groups. Then the
inequality [IS, 34] applied to this sequence is transformed to the bound

$\begin{split}&\sum_{k=1}^{2n}\sum_{i,j}|\pm\gamma_{2k,i}-\gamma_{2k\pm
1,j}|-\sum_{r=1}^{2n-1}\sum_{i\not=j}|\pm\gamma_{r,i}-\gamma_{r,j}|-2\sum_{k=1}^{n-1}\sum_{i}|\gamma_{2k,i}|\geq\\\
&C({\bm{\gamma}}_{2n})+\frac{2}{n}\sum_{r=1}^{2n-1}\sum_{i}|\gamma_{r,i}|,\end{split}$
(4.7)

where the constant $C({\bm{\gamma}}_{2n})$ depends on the values of
$\gamma_{2n,i}$. Due to the asymptotics of the Gamma function

$\Gamma(ix)\sim e^{-\pi|x|/{2}}$

in imaginary direction we observe that the integrand of (4.5) can be bounded
by (4.7) as

$C^{\prime}({\bm{\gamma}}_{2n})e^{-\frac{\pi-\varepsilon}{cn}\sum_{r=1}^{2n-1}\sum_{i}|\gamma_{r,i}|}$
(4.8)

for any small positive $\varepsilon>0$ and a proper positive constant
$C^{\prime}({\bm{\gamma}}_{2n})$, which implies absolute convergence of the
integral (4.5)

###  4.3 Toda equation

Denote by $H_{n}^{B}$ the Toda hamiltonian

$H_{B_{n}}=\sum_{k=1}^{n}\left(-\frac{\partial^{2}}{\partial
x_{k}^{2}}+(2n-2k+1)\frac{\partial}{\partial
x_{k}}\right)+\sum_{k=1}^{n-1}\frac{2}{c^{2}}e^{x_{k}-x_{k+1}}+\frac{1}{c^{2}}e^{x_{n}}$
(4.9)

######  Theorem 2.

The Whittaker function (4.5) is a wave function for $B_{n}$ Toda Hamiltonian:

$H_{B_{n}}\Psi_{{\bm{\gamma}}_{2n}}(\mathbf{x}_{n})=\frac{1}{c^{2}}\sum_{i=1}^{n}\gamma_{2n,i}^{2}+\frac{n(2n-1)(2n+1)}{12}$
(4.10)

Proof is a standard game with the matrix element

$(w_{n}^{\prime}\,L_{2n+1}\,e^{-\sum_{k}x_{k}F_{kk}}w_{n})_{2n}$

where

$L_{2n+1}=\sum_{i,j=1}^{2n+1}I_{ij}^{2}=\sum_{i,j=1}^{2n+1}F_{ij}F_{ji}$

is Laplace operator of $SO(2n+1)$. For this one should also know the
eigenvalue of $L$ in the representation in meromorphic functions. But it is
known from the theory of highest weight representations of $so(2n+1)$. It
gives us the eigenvalue

$\sum_{i=1}^{n}m_{2n,i}^{2}-(\rho,\rho)$

where

$\rho=(n-\frac{1}{2},\ldots,\frac{1}{2})$

is a half sum of positive roots for $so(2n+1)$. Thus the eigenvalue is equal
to

$-\frac{1}{c^{2}}\sum_{i=1}^{n}\gamma_{2n,i}^{2}-\frac{n(2n-1)(2n+1)}{12}$

$\Box$

Remark. Note that the function

$\tilde{\Psi}_{{\bm{\gamma}}_{2n}}(\mathbf{x}_{n})=e^{-(\rho,\mathbf{x}_{n})}{\Psi}_{{\bm{\gamma}}_{2n}}(\mathbf{x}_{n})=e^{-\sum\limits^{n}_{k=1}(n-k+\frac{1}{2})x_{k}}{\Psi}_{{\bm{\gamma}}_{2n}}(\mathbf{x}_{n})$

is the solution of more familiar spectral problem

$\left(-\sum_{k=1}^{n}\frac{\partial^{2}}{\partial
x_{k}^{2}}+\sum_{k=1}^{n-1}\frac{2}{c^{2}}e^{x_{k}-x_{k+1}}+\frac{1}{c^{2}}e^{x_{n}}\right)\tilde{\Psi}_{{\bm{\gamma}}_{2n}}(\mathbf{x}_{n})=\frac{1}{c^{2}}\left(\sum_{i=1}^{n}\gamma_{2n,i}^{2}\right)\tilde{\Psi}_{{\bm{\gamma}}_{2n}}(\mathbf{x}_{n})$

###  4.4 Iterative procedures

1. The integral (4.5) can be formulated as an iterative integral presentation of the Whittaker wave function,

$\begin{split}\Psi_{{\bm{\gamma}}_{2n}}(x_{1},\ldots,x_{n})=\int_{\mathbb{R}^{2n-1}}{\mu}({\bm{\gamma}}_{2n-1})d{\bm{\gamma}}_{2n-1}\,{\mu}({\bm{\gamma}}_{2n-2})d{\bm{\gamma}}_{2n-2}\\\
\prod_{i=1}^{n-1}\prod_{l=1}^{n}\Gamma\left(\frac{\pm\gamma_{2n-2,i}-\gamma_{2n-1,l}}{ic}+\frac{1}{2}\right)\prod_{i,j=1}^{n}\Gamma\left(\frac{\pm\gamma_{2n,i}+\gamma_{2n-1,j}}{ic}+\frac{1}{2}\right)\times\\\
c^{\frac{2}{ic}\delta_{2n-1}}e^{-\frac{\delta_{2n-1}}{ic}x_{n}+((n-1)+\sum_{k=1}^{n-2}k)x_{n}}\Psi_{{\bm{\gamma}}_{2n-2}}(x_{1}-x_{n},\ldots,x_{n-1}-x_{n})\end{split}$
(4.11)

or

$\Psi_{{\bm{\gamma}}_{2n}}(x_{1},\ldots,x_{n})=\Lambda_{n}(x_{n})\left(\Psi_{{\bm{\gamma}}_{2n-2}}(x_{1}-x_{n},\ldots,x_{n-1}-x_{n})\right)$
(4.12)

where $\Lambda(x_{n})$ is an integral operator

$(\Lambda_{n}(x)f)({\bm{\gamma}}_{2n})=\int_{\mathbb{R}^{n-1}}K({\bm{\gamma}}_{2n};{\bm{\gamma}}_{2n-2}|x)f({\bm{\gamma}}_{2n-2})\mu({\bm{\gamma}}_{2n-2})d{\bm{\gamma}}_{2n-2}$
(4.13)

with the kernel

$\begin{split}K({\bm{\gamma}}_{2n};{\bm{\gamma}}_{2n-2}|x)=e^{\frac{n(n-1)}{2}x}\int_{\mathbb{R}^{n}}{\mu}({\bm{\gamma}}_{2n-1})d{\bm{\gamma}}_{2n-1}c^{\frac{2}{ic}\delta_{2n-1}}e^{-\frac{\delta_{2n-1}}{ic}x}\\\
\prod_{i=1}^{n-1}\prod_{l=1}^{n}\Gamma\left(\frac{\pm\gamma_{2n-2,i}-\gamma_{2n-1,l}}{ic}+\frac{1}{2}\right)\prod_{i,j=1}^{n}\Gamma\left(\frac{\pm\gamma_{2n,i}+\gamma_{2n-1,j}}{ic}+\frac{1}{2}\right)\end{split}$
(4.14)

2. Another recurrent procedure uses the observation that the above construction of the Whittaker vectors and Whittaker functions, restricted to $SO(2n)$ produced actually $GL_{n}$ Whittaker vectors and functions. It can be seen from the relation (A.1), (A.2). In this way we arrive by using Gustafson integrals [Gust] to Iorgov - Shadura formula [IS], which expresses the $B_{n}$ Toda wave function via $A_{n}$ Toda wave function and contains in total twice less integrals.

The restriction $\Psi_{{\bm{\gamma}}_{2n-1}}(\mathbf{x}_{n})$ of the wave
function (4.5) to $SO(2n)$ is given by the integral

$\begin{split}\Psi_{{\bm{\gamma}}_{2n-1}}(\mathbf{x}_{n})=&\int_{\mathbb{R}^{n^{2}-n}}\prod_{k=1}^{2n-2}{\mu}({\bm{\gamma}}_{k})d{\bm{\gamma}}_{k}\,\cdot
c^{\frac{2}{ic}\sum_{k=1}^{n-1}\delta_{2k-1}}e^{\frac{1}{ic}\sum_{k=1}^{n}(\delta_{2k-3}-\delta_{2k-1}+(k-1)ic)x_{k}}\times\\\
&\prod_{k=1}^{n-1}\prod_{i,j=1}^{k}\prod_{l=1}^{k+1}\Gamma\left(\frac{\pm\gamma_{2k,i}+\gamma_{2k-1,j}}{ic}+\frac{1}{2}\right)\Gamma\left(\frac{\pm\gamma_{2k,i}-\gamma_{2k+1,l}}{ic}+\frac{1}{2}\right)\end{split}$
(4.15)

and the functions $\Psi_{{\bm{\gamma}}_{2n}}(\mathbf{x}_{n})$ and
$\Psi_{{\bm{\gamma}}_{2n-1}}(\mathbf{x}_{n})$ are related as follows

$\Psi_{{\bm{\gamma}}_{2n}}(\mathbf{x}_{n})=\int_{\mathbb{R}^{n}}{\mu}({\bm{\gamma}}_{2n-1})d{\bm{\gamma}}_{2n-1}c^{\frac{2}{ic}\delta_{2n-1}}\prod_{i,j=1}^{n}\Gamma\left(\frac{\pm\gamma_{2n,i}+\gamma_{2n-1,j}}{ic}+\frac{1}{2}\right)\Psi_{{\bm{\gamma}}_{2n-1}}(\mathbf{x}_{n})$
(4.16)

For each $k$ the integral over ${\bm{\gamma}}_{2k}$ in (4.15) can be
explicitly calculated by means of the degenerate $B_{n}$ Gustafson integral

$\begin{split}\frac{1}{(2\pi)^{n}}\int\limits_{\mathbb{R}^{n}}\frac{\prod\limits_{i=1}^{2n+1}\prod\limits_{j=1}^{n}\Gamma(a_{i}+iz_{j}))\Gamma(a_{i}-iz_{j}))}{\prod\limits_{1\leq
i<j\leq
n}\big{|}\Gamma(i(z_{i}+z_{j}))\Gamma(i(z_{i}-z_{j}))\big{|}^{2}\prod\limits_{j=1}^{n}\big{|}\Gamma(2iz_{j})\big{|}^{2}}d\mathbf{x}_{n}={n!2^{n}\prod\limits_{1\leq
i<j\leq 2n+1}\Gamma(a_{i}+a_{j})}\;,\end{split}$ (4.17)

where all $a_{i}$ are assumed to have a positive real part. The integral
(4.17) is a limiting case

$a_{2n+2}=\varepsilon+iL,\qquad L\to\infty,\qquad\varepsilon\to+0$

of a general $B_{n}$ Gustafson integral [Gust]

$\frac{1}{(2\pi)^{n}}\int\limits_{\mathbb{R}^{n}}\frac{\prod\limits_{i=1}^{2n+2}\prod\limits_{j=1}^{n}\Gamma(a_{i}+iz_{j}))\Gamma(a_{i}-iz_{j}))}{\prod\limits_{1\leq
i<j\leq
n}\big{|}\Gamma(i(z_{i}+z_{j}))\Gamma(i(z_{i}-z_{j}))\big{|}^{2}\prod\limits_{j=1}^{n}\big{|}\Gamma(2iz_{j})\big{|}^{2}}d\mathbf{x}_{n}=\frac{n!2^{n}\prod\limits_{1\leq
i<j\leq
2n+2}\Gamma(a_{i}+a_{j})}{\Gamma\Big{(}\sum\limits_{i=1}^{2n+2}a_{i}\Big{)}}.$

Using

$a_{1}=\frac{\gamma_{2k-1,1}}{ic}+\frac{1}{2},\ldots,a_{k}=\frac{\gamma_{2k-1,k}}{ic}+\frac{1}{2},\qquad
a_{k+1}=-\frac{\gamma_{2k+1,1}}{ic}+\frac{1}{2},\ldots
a_{2k+1}=-\frac{\gamma_{2k+1,k+1}}{ic}+\frac{1}{2},$

we get

$\begin{split}&\int\limits_{\mathbb{R}^{k}}\frac{\prod\limits^{k}_{i=1}\left(\prod\limits^{k}_{j=1}\Gamma\left(\frac{\pm\gamma_{2k,i}+\gamma_{2k-1,j}}{ic}+\frac{1}{2}\right)\prod\limits^{k+1}_{j=1}\Gamma\left(\frac{\pm\gamma_{2k,i}-\gamma_{2k+1,j}}{ic}+\frac{1}{2}\right)\right)}{\prod\limits_{r<s}\left|\Gamma\Big{(}\frac{\gamma_{2k,r}-\gamma_{2k,s}}{ic}\Big{)}\right|^{2}\left|\Gamma\Big{(}\frac{\gamma_{2k,r}+\gamma_{2k,s}}{ic}\big{)}\right|^{2}\prod\limits_{r}\left|\Gamma\Big{(}\frac{2\gamma_{2k,r}}{ic}\Big{)}\right|^{2}}d\gamma_{2k,1}\ldots
d\gamma_{2k,k}=\\\ &=c^{k}(2\pi)^{k}2^{k}\cdot
k!\cdot\prod\limits_{r<s}\Gamma\left(\frac{\gamma_{2k-1,r}+\gamma_{2k-1,s}}{ic}+1\right)\cdot\prod\limits_{i,j}\Gamma\left(\frac{\gamma_{2k-1,i}-\gamma_{2k+1,j}}{ic}+1\right)\\\
&\qquad\qquad\times\prod\limits_{r<s}\Gamma\left(-\frac{\gamma_{2k+1,r}+\gamma_{2k+1,s}}{ic}+1\right)\end{split}$
(4.18)

and

$\Psi_{{\bm{\gamma}}_{2n-1}}(\mathbf{x}_{n})=d_{n}\cdot
c^{\frac{n+1}{ic}\delta_{2n-1}}\prod_{1\leq r<s\leq
n}\Gamma\left(1-\frac{\gamma_{2n-1,r}+\gamma_{2n-1,s}}{ic}\right)\Phi_{{\bm{\gamma}}_{2n-1}}(\mathbf{x}_{n})$
(4.19)

where

$\begin{split}\Phi_{{\bm{\gamma}}_{2n-1}}(\mathbf{x}_{n})=\int\limits_{\mathbb{R}^{n(n-1)/2}}e^{\frac{1}{ic}\sum\limits^{n}_{k=1}(\delta_{2k-3}-\delta_{2k-1}+(k-1)ic)x_{k}}d\gamma_{1}\ldots
d\gamma_{2n-3}\times\\\
\frac{\prod\limits^{n-1}_{k=1}\prod\limits^{k}_{i=1}\prod\limits^{k+1}_{j=1}c^{\frac{\gamma_{2k-1,i}-\gamma_{2k+1,j}}{ic}}\Gamma\left(\frac{\gamma_{2k-1,i}-\gamma_{2k+1,j}}{ic}+1\right)}{\prod\limits^{n-1}_{k=1}\prod\limits_{r<s}\left|\Gamma\left(\frac{\gamma_{2k-1,r}-\gamma_{2k-1,s}}{ic}\right)\right|^{2}}\end{split}$
(4.20)

and

$d_{n}=\prod_{k=1}^{n-1}c^{k}(2\pi)^{k}2^{k}\cdot k!$ (4.21)

Both functions $\Psi_{\mathbf{x}_{n}}(\gamma_{2n-1})$ and
$\Phi_{{\bm{\gamma}}_{2n-1}}(\mathbf{x}_{n})$ are solutions of $GL(n)$ Toda
equations. In order to compare the final results with Iorgov - Shadura formula
[IS] we perform the change of integration variables

$\gamma_{2k-1,j}\to\gamma_{2k-1,j}-\left(n-k+\frac{1}{2}\right)ic.$

Then

$\begin{split}\Phi_{{\bm{\gamma}}_{2n-1}}(\mathbf{x}_{n})=e^{\sum\limits^{n}_{k=1}(n-k+\frac{1}{2})x_{k}}&\int\limits_{C_{n}}e^{\frac{1}{ic}\sum\limits^{n}_{k=1}(\delta_{2k-3}-\delta_{2k-1})x_{k}}d\gamma_{1}\ldots
d\gamma_{2n-3}\times\\\
&\frac{\prod\limits^{n-1}_{k=1}\prod\limits^{k}_{i=1}\prod\limits^{k+1}_{j=1}c^{\frac{\gamma_{2k-1,i}-\gamma_{2k+1,j}}{ic}}\Gamma\left(\frac{\gamma_{2k-1,i}-\gamma_{2k+1,j}}{ic}\right)}{\prod\limits^{n-1}_{k=1}\prod\limits_{r<s}\left|\Gamma\left(\frac{\gamma_{2k-1,r}-\gamma_{2k-1,s}}{ic}\right)\right|^{2}}\end{split}$
(4.22)

and

$\begin{split}\Psi_{{\bm{\gamma}}_{2n}}(\mathbf{x}_{n})=d_{n}&\int\limits_{C^{\prime}_{n}}\frac{\prod\limits^{n}_{i=1}\prod\limits^{n}_{j=1}\Gamma\left(\frac{\gamma_{2n-1,i}+\gamma_{2n,j}}{ic}\right)\Gamma\left(\frac{\gamma_{2n-1,i}-\gamma_{2n,j}}{ic}\right)}{\prod\limits_{r<s}\Gamma\left(\frac{\gamma_{2n-1,r}+\gamma_{2n-1,s}}{ic}\right)\prod\limits_{r<s}\left|\Gamma\left(\frac{\gamma_{2n-1,r}-\gamma_{2n-1,s}}{ic}\right)\right|^{2}}\times\\\
&c^{\frac{n+1}{ic}\delta_{2n-1}}\Phi_{{\bm{\gamma}}_{2n-1}}(\mathbf{x}_{n})d\gamma_{2n-1}\end{split}$
(4.23)

Here the contour $C_{n}$ in (4.22) is a deformation of
$\mathbb{R}^{n(n-1)/2}$, such that the integration over the variable
$\gamma_{2k-1}$ is performed in such a way that the singularity of
$\Gamma\left(\frac{\gamma_{2k-1,i}-\gamma_{2k+1,j}}{ic}\right)$ is under the
line of integration and the singularity of
$\Gamma\left(\frac{\gamma_{2k-3,i}-\gamma_{2k-1,j}}{ic}\right)$ is above the
line of integration over $\gamma_{2k-1}$. The contour $C^{\prime}_{n}$ in
(4.23) is a deformation of $\mathbb{R}^{n}$ where the singularities of the
nominators are under the contours of integrations over all variables.

The relations (4.22) - (4.23) are in accordance with Iorgov - Shadura formula
[IS, (26),(27)]. More precisely, in Iorgov - Shadura description the boundary
wall corresponds to the first coordinate $x_{1}$, while we work with the
boundary wall related to the last coordinate. One can observe the coincidence
of formulas after the change of variables

$x_{k}\to-x_{n+1-k}$

in their (or in ours) formulas.

##  5 Examples

$\ \mathbf{n=1}$. For one particle the system and the wave function coincide
with that of $sl(2)$ Toda system:

$\Psi_{\gamma_{21}}(x_{1})=\int_{\mathbb{R}}d\gamma_{11}c^{\frac{2}{ic}\gamma_{11}}e^{-\frac{1}{ic}\gamma_{11}x_{1}}\Gamma\left(\frac{\pm\gamma_{21}+\gamma_{11}}{ic}+\frac{1}{2}\right)$

$\ \mathbf{n=2}$. The wave function for $B_{2}$ Toda system is given by $4$
fold integral, which can be reduced to three fold by using Gustafson integral,

$\begin{split}\Psi_{\gamma_{41},\gamma_{42}}(x_{1},x_{2})=\int_{\mathbb{R}^{4}}&d\gamma_{11}\frac{d\gamma_{21}}{\Big{|}\Gamma\Big{(}\frac{2\gamma_{21}}{ic}\Big{)}\Big{|}^{2}}\frac{d\gamma_{31}d\gamma_{32}}{\Big{|}\Gamma\Big{(}\frac{\gamma_{31}-\gamma_{32}}{ic}\Big{)}\Gamma\Big{(}1+\frac{\gamma_{31}+\gamma_{32}}{ic}\Big{)}\Big{|}^{2}}\times\\\
&c^{\frac{2}{ic}(\gamma_{11}+\gamma_{31}+\gamma_{32})}e^{\frac{1}{ic}\left(-\gamma_{11}x_{1}+(\gamma_{11}-\gamma_{31}-\gamma_{32}+ic)x_{2}\right)}\times\\\
\Gamma\left(\frac{\pm\gamma_{21}+\gamma_{11}}{ic}+\frac{1}{2}\right)&\Gamma\left(\frac{\pm\gamma_{21}-\gamma_{31}}{ic}+\frac{1}{2}\right)\Gamma\left(\frac{\pm\gamma_{21}-\gamma_{32}}{ic}+\frac{1}{2}\right)\prod_{i,j=1}^{2}\Gamma\left(\frac{\pm\gamma_{4i}+\gamma_{3j}}{ic}+\frac{1}{2}\right)\end{split}$

that is

$\begin{split}&\Psi_{\gamma_{41},\gamma_{42}}(x_{1},x_{2})=\\\
&\int_{\mathbb{R}^{3}}\frac{d\gamma_{21}}{\Big{|}\Gamma\Big{(}\frac{2\gamma_{21}}{ic}\Big{)}\Big{|}^{2}}\frac{d\gamma_{31}d\gamma_{32}}{\Big{|}\Gamma\Big{(}\frac{\gamma_{31}-\gamma_{32}}{ic}\Big{)}\Gamma\Big{(}1+\frac{\gamma_{31}+\gamma_{32}}{ic}\Big{)}\Big{|}^{2}}c^{\frac{2}{ic}(\gamma_{31}+\gamma_{32})}e^{\frac{1}{ic}\left((-\gamma_{31}-\gamma_{32}+ic)x_{2}\right)}\\\
&\Gamma\left(\frac{\pm\gamma_{21}-\gamma_{31}}{ic}+\frac{1}{2}\right)\Gamma\left(\frac{\pm\gamma_{21}-\gamma_{32}}{ic}+\frac{1}{2}\right)\prod_{i,j=1}^{2}\Gamma\left(\frac{\pm\gamma_{4i}+\gamma_{3j}}{ic}+\frac{1}{2}\right)\Psi_{\gamma_{21}}(x_{1}-x_{2})\end{split}$

or

$\begin{split}\Psi_{\gamma_{41},\gamma_{42}}(x_{1},x_{2})=\int_{\mathbb{R}^{2}+\varepsilon}\frac{\prod_{i,j=1}^{2}\Gamma\left(\frac{\pm\gamma_{4i}+\gamma_{3j}}{ic}\right)}{\Gamma\left(\frac{\gamma_{31}+\gamma_{32}}{ic}\right)\left|\Gamma\left(\frac{\gamma_{31}-\gamma_{32}}{ic}\right)\right|^{2}}c^{\frac{3}{ic}(\gamma_{31}+\gamma_{32})}\Phi_{\gamma_{31},\gamma_{32}}(x_{1},x_{2})d\gamma_{31}d\gamma_{32}\end{split}$

where $\Phi_{\gamma_{31},\gamma_{32}}(x_{1},x_{2})$ is the wave function of
$GL(2)$ Toda system

$\begin{split}\Phi_{\gamma_{31},\gamma_{32}}(x_{1},x_{2})=e^{\frac{3x_{1}+x_{2}}{2}}&\int_{\mathbb{R}+\varepsilon}d\gamma_{11}e^{\frac{1}{ic}(-\gamma_{11}x_{1}+(\gamma_{11}-\gamma_{31}-\gamma_{32})x_{2})}c^{\frac{1}{ic}(2\gamma_{11}-\gamma_{31}-\gamma_{32})}\\\
&\Gamma\left(\frac{\gamma_{11}-\gamma_{31}}{ic}\right)\Gamma\left(\frac{\gamma_{11}-\gamma_{32}}{ic}\right)\end{split}$

## Acknowlegements

The authors thank N.Belousov, S.Derkachov and S.Kharchev for their interest in
the work and for fruitul discussions. The work of S. Kh. is an output of a
research project implemented as part of the Basic Research Program at the
National Research University Higher School of Economics (HSE University).

##  Appendix A Calculations of Whittaker vector

In this section we prove Theorem 1. First we note that the relation (3.14)
follows from (3.13) by using the automorphism $\tau$. Proof of the equality
(3.13) reduces to check of the following equalities, where $k=1,\ldots n-1$ in
the relations (A.2) and (A.3)

$\displaystyle(F_{k,k+1}-(-1)^{k}F_{k,2n-k})w_{n}=(I_{2k+1,2k}+iI_{2k+1,2k-1})w_{n}=\frac{(-1)^{k+1}}{ic}w_{n}$
(A.1)
$\displaystyle(F_{k,k+1}+(-1)^{k}\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}F_{k,2n-k})w_{n}=(-I_{2k+2,2k-1}+iI_{2k+2,2k})w_{n}=\frac{(-1)^{k+1}}{ic}w_{n}$
(A.2) $\displaystyle
F_{n,n+1}w_{n}=\frac{1}{\sqrt{2}}(I_{2n+1,2n}+iI_{2n+1,2n-1})w_{n}=\frac{(-1)^{n+1}}{ic\sqrt{2}}w_{n}$
(A.3)

Besides, the relation (A.3) is a particular case of (A.1). So we have to prove
the relations (A.1) and (A.2).

The proof of (A.1)–(A.2) requires certain calculations. For their
visualization we introduce some intermediate notations. First rewrite the
operators (2.10) and (2.11) as

$\begin{split}I_{2k+1,2k}=\sum\limits_{\varepsilon=\pm
1}\sum\limits_{j=1}^{k}P^{\varepsilon}_{kj},\qquad
I_{2k+2,2k+1}=\sum\limits_{j=1}^{k}Q_{kj}+\sum\limits_{j=1}^{k}R_{kj}+T_{k}\end{split}$
(A.4)

where

$\begin{split}P^{\varepsilon}_{k,j}=&-\frac{1}{ic}\frac{\prod\limits_{r=1}^{k-1}(\nu_{2k-1,j}+\varepsilon(\nu_{2k-2,r}+\frac{ic}{2}))\prod\limits_{r=1}^{k}(\nu_{2k-1,j}-\varepsilon(\nu_{2k,r}-\frac{ic}{2}))}{2\prod\limits_{r\neq
j}(\nu_{2k-1,j}-\nu_{2k-1,r})(\nu_{2k-1,j}+\nu_{2k-1,r}+\varepsilon
ic)}e^{\varepsilon ic\partial_{\nu_{2k-1,j}}},\\\
Q_{k,j}=&\frac{1}{ic}\frac{\prod\limits_{r=1}^{k+1}((\nu_{2k,j}+\frac{ic}{2})^{2}-\nu_{2k+1,r}^{2})}{2\nu_{2k,j}(\nu_{2k,j}+\frac{ic}{2})\prod\limits_{r\neq
j}(\nu_{2k,j}^{2}-\nu_{2k,r}^{2})}e^{ic\partial_{\nu_{2k,j}}}\\\
R_{k,j}=&\frac{1}{ic}\frac{\prod\limits_{r=1}^{k}((\nu_{2k,j}-\frac{ic}{2})^{2}-\nu_{2k-1,r}^{2})}{2\nu_{2k,j}(\nu_{2k,j}-\frac{ic}{2})\prod\limits_{r\neq
j}(\nu_{2k,j}^{2}-\nu_{2k,r}^{2})}e^{-ic\partial_{\nu_{2k,j}}}\\\
T_{k}=&\frac{1}{c}\frac{\prod\limits_{r=1}^{k}\nu_{2k-1,r}\prod\limits_{r=1}^{k+1}\nu_{2k+1,r}}{\prod\limits_{r=1}^{k}(\nu_{2k,r}^{2}+\frac{c^{2}}{4})}\end{split}$

Set also for $\varepsilon,\delta=\pm 1$

$J^{\varepsilon}_{k,\delta,j}=\sum\limits_{s=1}^{k}\frac{\varepsilon
ic}{\nu_{2k+\delta,j}-\varepsilon(\nu_{2k,s}+\frac{ic}{2})}Q_{k,s}+\sum\limits_{s=1}^{k}\frac{\varepsilon
ic}{\nu_{2k+\delta,j}+\varepsilon(\nu_{2k,s}-\frac{ic}{2})}R_{k,s}+\frac{\varepsilon
ic}{\nu_{2k+\delta,j}}T_{k}.$ (A.5)

In these notation we can present the following expressions for other
generators of Lie algebra $so(n)$ needed in the equations on Whittaker
vectors. They can be checked by straightforward calculations.

######  Lemma 3.

We have the following relations

$\begin{split}&I_{2k+1,2k-1}=\sum\limits_{\varepsilon=\pm
1}\sum\limits_{j=1}^{k}P^{\varepsilon}_{kj}J_{k-1,1,j}^{\varepsilon}\qquad
k=1,\ldots,n\\\ &I_{2k+2,2k}=-\sum\limits_{\varepsilon=\pm
1}\sum\limits_{j=1}^{k}P^{\varepsilon}_{kj}J_{k,-1,j}^{\varepsilon}\qquad
k=1,\ldots,n-1\\\ &I_{2k+2,2k-1}=-\sum\limits_{\varepsilon=\pm
1}\sum\limits_{j=1}^{k}P^{\varepsilon}_{kj}J_{k,-1,j}^{\varepsilon}J_{k-1,1,j}^{\varepsilon}\qquad
k=1,\ldots,n-1\\\ &\end{split}$ (A.6)

For more brevity in futher formulas we denote

$\theta_{k}=(-1)^{k},\qquad\theta_{k+1}=(-1)^{k+1}.$

The following statement is one of the important steps in the proof of Theorem
1.

######  Lemma 4.

For any $k=0,\ldots,n-1$ ($\delta=1$ when $k=0$) we have the relations

$\displaystyle J^{\theta_{k}}_{k,\delta,j}w_{n}=iw_{n}$ (A.7) $\displaystyle
J^{\theta_{k+1}}_{k,\delta,j}w_{n}=i\left(1-\frac{\prod\limits_{r=1}^{k+1}(\nu_{2k+\delta,j}+\nu_{2k+1,r})\prod\limits_{r=1}^{k}(\nu_{2k+\delta,j}+\nu_{2k-1,r})}{\nu_{2k+\delta,j}\prod\limits_{r=1}^{k}((\nu_{2k+\delta,j}+(-1)^{k}\frac{ic}{2})^{2}-\nu_{2k,r}^{2})}\right)w_{n}\qquad\delta=\pm
1$ (A.8)

Proof . Note that the operator $J^{\theta_{k}}_{k,\delta,j}$ containes only
shifts of variables $\nu_{2k,s}$, and thus only the factor
$W_{k}^{\theta_{k+1}}$ can change. Let us check the equality (A.7). We then
have

$\begin{split}&J^{\theta_{k}}_{k,\delta,j}w_{n}=\\\
\bigg{(}\sum\limits_{s=1}^{k}&\frac{\theta_{k}ic}{\nu_{2k+\delta,j}-\theta_{k}(\nu_{2k,s}+\frac{ic}{2})}Q_{k,s}+\sum\limits_{s=1}^{k}\frac{\theta_{k}ic}{\nu_{2k+\delta,j}+\theta_{k}(\nu_{2k,s}-\frac{ic}{2})}R_{k,s}+\frac{\theta_{k}ic}{\nu_{2k+\delta,j}}T_{k}\bigg{)}w_{n}\\\
=&\Bigg{(}\sum\limits_{s=1}^{k}\frac{\theta_{k}ic}{\nu_{2k+\delta,j}-\theta_{k}(\nu_{2k,s}+\frac{ic}{2})}\cdot\frac{1}{ic}\cdot\frac{\prod\limits_{r=1}^{k+1}((\nu_{2k,s}+\frac{ic}{2})^{2}-\nu_{2k+1,r}^{2})}{2\nu_{2k,s}(\nu_{2k,s}+\frac{ic}{2})\prod\limits_{r\neq
s}(\nu_{2k,s}^{2}-\nu_{2k,r}^{2})}\times\\\
i^{\theta_{k}}\times&\frac{\prod\limits_{r=1}^{k}(\nu_{2k,s}-\theta_{k}\nu_{2k-1,r}+\frac{ic}{2})}{\prod\limits_{r=1}^{k+1}(-\nu_{2k,s}-\theta_{k}\nu_{2k+1,r}-\frac{ic}{2})}+\sum\limits_{s=1}^{k}\frac{\theta_{k}ic}{\nu_{2k+\delta,j}+\theta_{k}(\nu_{2k,s}-\frac{ic}{2})}\cdot\frac{1}{ic}\times\\\
&\frac{\prod\limits_{r=1}^{k+1}((\nu_{2k,s}-\frac{ic}{2})^{2}-\nu_{2k+1,r}^{2})}{2\nu_{2k,s}(\nu_{2k,s}-\frac{ic}{2})\prod\limits_{r\neq
s}(\nu_{2k,s}^{2}-\nu_{2k,r}^{2})}\cdot
i^{-\theta_{k}}\times\frac{\prod\limits_{r=1}^{k+1}(-\nu_{2k,s}-\theta_{k}\nu_{2k+1,r}+\frac{ic}{2})}{\prod\limits_{r=1}^{k}(\nu_{2k,s}-\theta_{k}\nu_{2k-1,r}-\frac{ic}{2})}+\\\
&i\frac{\prod\limits_{r=1}^{k}(-\theta_{k}\nu_{2k-1,r})\prod\limits_{r=1}^{k+1}(-\theta_{k}\nu_{2k+1,r})}{(-\theta_{k}\nu_{2k+\delta,j})\prod\limits_{r=1}^{k}(-\frac{ic}{2}-\nu_{2k,r})(-\frac{ic}{2}+\nu_{2k,r})}\Bigg{)}w_{n}=iw_{n}\end{split}$

The last line is obtained by the use of the following well know identity:

$\sum\limits^{m}_{i=1}\frac{\prod\limits^{m-1}_{j=1}(x_{i}-y_{j})}{\prod\limits_{r\neq
i}(x_{i}-x_{r})}=1$ (A.9)

where for indeterminates $x_{i}$ we choose $2k+2$ variables
$\\{\pm\nu_{2k,s},-\frac{ic}{2},\theta_{k}\nu_{2k+\delta}-\frac{ic}{2}\\}$,
and for indeterminates $y_{i}$ choose
$\\{\theta_{k}\nu_{2k-1}-\frac{ic}{2},\theta_{k}\nu_{2k+1}-\frac{ic}{2}\\}$.

Let us check the equality (A.8). We have

$\begin{split}&J^{\theta_{k+1}}_{k,\delta,j}w_{n}=\\\
&\bigg{(}\sum\limits_{s=1}^{k}\frac{\theta_{k+1}ic}{\nu_{2k+\delta,j}-\theta_{k+1}(\nu_{2k,s}+\frac{ic}{2})}Q_{k,s}+\sum\limits_{s=1}^{k}\frac{\theta_{k+1}ic}{\nu_{2k+\delta,j}+\theta_{k+1}(\nu_{2k,s}-\frac{ic}{2})}R_{k,s}+\frac{\theta_{k+1}ic}{\nu_{2k+\delta,j}}T_{k}\bigg{)}w_{n}\end{split}$
$\begin{split}&=\Bigg{(}\sum\limits_{s=1}^{k}\frac{\theta_{k+1}ic}{\nu_{2k+\delta,j}-\theta_{k+1}(\nu_{2k,s}+\frac{ic}{2})}\cdot\frac{1}{ic}\cdot\frac{\prod\limits_{r=1}^{k+1}((\nu_{2k,s}+\frac{ic}{2})^{2}-\nu_{2k+1,r}^{2})}{2\nu_{2k,s}(\nu_{2k,s}+\frac{ic}{2})\prod\limits_{r\neq
s}(\nu_{2k,s}^{2}-\nu_{2k,r}^{2})}\times\end{split}$
$\begin{split}&i^{\theta_{k+1}}\times\frac{\prod\limits_{r=1}^{k}(\nu_{2k,s}-\theta_{k}\nu_{2k-1,r}+\frac{ic}{2})}{\prod\limits_{r=1}^{k+1}(-\nu_{2k,s}-\theta_{k}\nu_{2k+1,r}-\frac{ic}{2})}-\\\
&\sum\limits_{s=1}^{k}\frac{\theta_{k+1}ic}{\nu_{2k+\delta,j}+\theta_{k+1}(\nu_{2k,s}-\frac{ic}{2})}\cdot\frac{1}{ic}\cdot\frac{\prod\limits_{r=1}^{k+1}((\nu_{2k,s}-\frac{ic}{2})^{2}-\nu_{2k+1,r}^{2})}{2\nu_{2k,s}(\nu_{2k,s}-\frac{ic}{2})\prod\limits_{r\neq
s}(\nu_{2k,s}^{2}-\nu_{2k,r}^{2})}\times\\\
&i^{-\theta_{k+1}}\cdot\frac{\prod\limits_{r=1}^{k+1}(-\nu_{2k,s}-\theta_{k}\nu_{2k+1,r}+\frac{ic}{2})}{\prod\limits_{r=1}^{k}(\nu_{2k,s}-\theta_{k}\nu_{2k-1,r}-\frac{ic}{2})}+i\frac{\prod\limits_{r=1}^{k}(-\theta_{k}\nu_{2k-1,r})\prod\limits_{r=1}^{k+1}(-\theta_{k}\nu_{2k+1,r})}{(-\theta_{k+1}\nu_{2k+\delta})\prod\limits_{r=1}^{k}(-\frac{ic}{2}-\nu_{2k,r})(-\frac{ic}{2}+\nu_{2k,r})}\Bigg{)}w_{n}\\\
&=i\Bigg{(}1-\frac{\prod\limits_{r=1}^{k+1}(\nu_{2k+\delta,j}+\nu_{2k+1,r})\prod\limits_{r=1}^{k}(\nu_{2k+\delta,j}+\nu_{2k-1,r})}{\nu_{2k+\delta,j}\prod\limits_{r=1}^{k}((\nu_{2k+\delta,j}+\theta_{k}\frac{ic}{2})^{2}-\nu_{2k,r}^{2})}\Bigg{)}w_{n}\end{split}$

The last line of the equality is obtained after simplifications of the ratios
by use of the relation $(\ref{A8})$, where for inderteminates $x_{i}$ we
substitute $x$ variables
$\\{\pm\nu_{2k,s},-\frac{ic}{2},\theta_{k+1}\nu_{2k+\delta}-\frac{ic}{2}\\}$,
and for $y_{i}$ we substitute
$\\{\theta_{k}\nu_{2k-1}-\frac{ic}{2},\theta_{k}\nu_{2k+1}-\frac{ic}{2}\\}$.

Proof of (A.1). According to Lemmas 3 and 4, for any $k=1,\ldots,n$ we have:

$\displaystyle(I_{2k+1,2k}+iI_{2k+1,2k-1})w_{n}=\sum\limits_{\varepsilon=\pm
1}\sum\limits_{j=1}^{k}P^{\varepsilon}_{kj}(1+iJ^{\varepsilon}_{k-1,1,j})w_{n}=$
$\displaystyle\sum\limits_{j=1}^{k}P^{\theta_{k}}_{kj}\frac{\prod\limits_{r=1}^{k}(\nu_{2k-1,j}+\nu_{2k-1,r})\prod\limits_{r=1}^{k-1}(\nu_{2k-1,j}+\nu_{2k-3,r})}{\nu_{2k-1,j}\prod\limits_{r=1}^{k-1}((\nu_{2k-1,j}+\theta_{k-1}\frac{ic}{2})^{2}-\nu_{2k-2,r}^{2})}w_{n}=$
(A.10)
$\begin{split}&-\frac{1}{ic}\sum\limits_{j=1}^{k}\frac{\prod\limits_{r=1}^{k-1}(\nu_{2k-1,j}+\theta_{k}(\nu_{2k-2,r}+\frac{ic}{2}))\prod\limits_{r=1}^{k}(\nu_{2k-1,j}-\theta_{k}(\nu_{2k,r}-\frac{ic}{2}))}{2\prod\limits_{r\neq
j}(\nu_{2k-1,j}-\nu_{2k-1,r})(\nu_{2k-1,j}+\nu_{2k-1,r}+\theta_{k}ic)}e^{\theta_{k}ic\partial_{\nu_{2k-1,j}}}\times\\\
&2\frac{\prod\limits_{r\neq
j}^{k}(\nu_{2k-1,j}+\nu_{2k-1,r})\prod\limits_{r=1}^{k-1}(\nu_{2k-1,j}+\nu_{2k-3,r})}{\prod\limits_{r=1}^{k-1}((\nu_{2k-1,j}+\theta_{k-1}\frac{ic}{2})^{2}-\nu_{2k-2,r}^{2})}w_{n}=\\\
&-\frac{1}{ic}\sum\limits_{j=1}^{k}\frac{\prod\limits_{r=1}^{k-1}(\nu_{2k-1,j}+\nu_{2k-3,r}+\theta_{k}ic)}{\prod\limits_{r\neq
j}(\nu_{2k-1,j}-\nu_{2k-1,r})}\cdot\frac{\prod\limits_{r=1}^{k}\left(\nu_{2k-1,j}+\theta_{k+1}\nu_{2k,r}+\theta_{k}\frac{ic}{2}\right)}{\prod\limits_{r=1}^{k-1}\left(\nu_{2k-1,j}+\theta_{k+1}\nu_{2k-2,r}+\theta_{k}\frac{ic}{2}\right)}e^{\theta_{k}ic\partial_{\nu_{2k-1,j}}}w_{n}\end{split}$

In the product, presenting the function $w_{n}$, only the factors
$W_{k-1}^{\theta_{k}}$ and $W_{k}^{\theta_{k+1}}$ depend on the variables
$\nu_{2k-1,j}$. Thus, using (3.11), functional relations on the Euler Gamma
function and (A.9), we get

$\begin{split}&(I_{2k+1,2k}+iI_{2k+1,2k-1})w_{n}=\frac{\theta_{k+1}}{ic}\sum\limits_{j=1}^{k}\frac{\prod\limits_{r=1}^{k-1}(\nu_{2k-1,j}+\nu_{2k-3,r}+\theta_{k}ic)}{\prod\limits_{r\neq
j}(\nu_{2k-1,j}-\nu_{2k-1,r})}w_{n}=\frac{\theta_{k+1}}{ic}w_{n}\end{split}$
(A.11)

$\Box$

Proof of (A.2). Again, according to Lemmas 3 and 4, for any $k=1,\ldots,n$ we
have:

$\begin{split}&(-I_{2k+2,2k-1}+iI_{2k+2,2k})w_{n}=\sum\limits_{\varepsilon=\pm
1}\sum\limits_{j=1}^{k}P^{\varepsilon}_{kj}J^{\varepsilon}_{k,-1,j}(J^{\varepsilon}_{k-1,1,j}-i)w_{n}=\\\
&-i\sum\limits_{j=1}^{k}P^{\theta_{k}}_{kj}J^{\theta_{k}}_{k,-1,j}\frac{\prod\limits_{r=1}^{k}(\nu_{2k-1,j}+\nu_{2k-1,r})\prod\limits_{r=1}^{k-1}(\nu_{2k-1,j}+\nu_{2k-3,r})}{\nu_{2k-1,j}\prod\limits_{r=1}^{k-1}((\nu_{2k-1,j}+\theta_{k-1}\frac{ic}{2})^{2}-\nu_{2k-2,r}^{2})}w_{n}.\end{split}$
(A.12)

Since the operator $J^{\theta_{k}}_{k,-1,j}$ contains shifts only of the
variables $\nu_{2k,s}$, we can rewright the result as

$-i\sum\limits_{j=1}^{k}P^{\theta_{k}}_{kj}\frac{\prod\limits_{r=1}^{k}(\nu_{2k-1,j}+\nu_{2k-1,r})\prod\limits_{r=1}^{k-1}(\nu_{2k-1,j}+\nu_{2k-3,r})}{\nu_{2k-1,j}\prod\limits_{r=1}^{k-1}((\nu_{2k-1,j}+\theta_{k-1}\frac{ic}{2})^{2}-\nu_{2k-2,r}^{2})}J^{\theta_{k}}_{k,-1,j}w_{n}$
(A.13)

Now Lemma 4 says that

$J^{\theta_{k}}_{k,-1,j}w_{n}=iw_{n}.$

Thus $(-I_{2k+2,2k-1}+iI_{2k+2,2k})w_{n}$ equals to

$=\sum\limits_{j=1}^{k}P^{\theta_{k}}_{kj}\frac{\prod\limits_{r=1}^{k}(\nu_{2k-1,j}+\nu_{2k-1,r})\prod\limits_{r=1}^{k-1}(\nu_{2k-1,j}+\nu_{2k-3,r})}{\nu_{2k-1,j}\prod\limits_{r=1}^{k-1}((\nu_{2k-1,j}+\theta_{k-1}\frac{ic}{2})^{2}-\nu_{2k-2,r}^{2})}w_{n}=\frac{\theta_{k+1}}{ic}w_{n}$
(A.14)

The latter equality was proved during the derivation of (A.11) from (A.10).

$\Box$

This ends the proof of Theorem 1. $\Box$

Remark. Analyzing the proof of Theorem 1 we see that the derivations over
variables $\nu_{2k-1,j}$ enter the game only in the last stage of
calculations. Moreover, we can freely add to Whittaker vectors factors of the
form

$e^{i\frac{\alpha_{j}}{c}\delta_{2j-1}},\qquad\text{where}\qquad\delta_{2j-1}=\sum_{i}\nu_{2j-1,i}=(-1)^{j+1}\sum_{i}\gamma_{2j-1,i}$

where $\alpha$ is arbitrary real number. This does not affect to the
convergence of integrals and does not change the action of Cartan subalgebra.
In Toda equation we earn thus arbitrary positive constants
$c_{j}=e^{\alpha_{j}}$ at exponentials $e^{x_{j-1}-x_{j}}$ which can be
equivalently obtained by successice shifts of the variables $x_{j}$.

##  Appendix B Action of Cartan subalgebra

It is sufficient to calculate $F_{k,k}w_{n}$, and then use the automorphism
$\tau$.

$\begin{split}&-iI_{2k,2k-1}w_{n}=\\\
&-\Bigg{(}\frac{1}{c}\sum_{j=1}^{k-1}\frac{\prod\limits_{r=1}^{k}((\nu_{2k-2,j}+\frac{ic}{2})^{2}-\nu_{2k-1,r}^{2})}{2\nu_{2k-2,j}(\nu_{2k-2,j}+\frac{ic}{2})\prod\limits_{r\neq
j}(\nu_{2k-2,j}^{2}-\nu_{2k-2,r}^{2})}i^{\theta_{k-1}}\frac{\prod\limits_{r=1}^{k-1}(\nu_{2k-2,j}-\theta_{k-1}\nu_{2k-3,r}+\frac{ic}{2})}{\prod\limits_{r=1}^{k}(-\nu_{2k-2,j}-\theta_{k-1}\nu_{2k-1,r}-\frac{ic}{2})}\\\
&+\frac{1}{c}\sum_{j=1}^{k-1}\frac{\prod\limits_{r=1}^{k-1}((\nu_{2k-2,j}-\frac{ic}{2})^{2}-\nu_{2k-3,r}^{2})}{2\nu_{2k-2,j}(\nu_{2k-2,j}-\frac{ic}{2})\prod\limits_{r\neq
j}(\nu_{2k-2,j}^{2}-\nu_{2k-2,r}^{2})}i^{-\theta_{k-1}}\frac{\prod\limits_{r=1}^{k}(-\nu_{2k-2,j}-\theta_{k-1}\nu_{2k-1,r}+\frac{ic}{2})}{\prod\limits_{r=1}^{k-1}(\nu_{2k-2,j}-\theta_{k-1}\nu_{2k-3,r}-\frac{ic}{2})}-\\\
&-\frac{1}{ic}\frac{\prod\limits_{r=1}^{k-1}\nu_{2k-3,r}\prod\limits_{r=1}^{k}\nu_{2k-1,r}}{\prod\limits_{r=1}^{k-1}(\nu_{2k-2,r}^{2}+\frac{c^{2}}{4})}\Bigg{)}w_{n}=\end{split}$
(B.1)
$\begin{split}=-&\frac{1}{ic}\Bigg{(}\sum_{j=1}^{k-1}\frac{\prod\limits_{r=1}^{k}(\nu_{2k-2,j}-\theta_{k-1}\nu_{2k-1,r}+\frac{ic}{2})\prod\limits_{r=1}^{k-1}(\nu_{2k-2,j}-\theta_{k-1}\nu_{2k-3,r}+\frac{ic}{2})}{2\nu_{2k-2,j}(\nu_{2k-2,j}+\frac{ic}{2})\prod\limits_{r\neq
j}(\nu_{2k-2,j}^{2}-\nu_{2k-2,r}^{2})}+\\\
&\sum_{j=1}^{k-1}\frac{\prod\limits_{r=1}^{k-1}(-\nu_{2k-2,j}-\theta_{k-1}\nu_{2k-3,r}+\frac{ic}{2})\prod\limits_{r=1}^{k}(-\nu_{2k-2,j}-\theta_{k-1}\nu_{2k-1,r}+\frac{ic}{2})}{2\nu_{2k-2,j}(\nu_{2k-2,j}-\frac{ic}{2})\prod\limits_{r\neq
j}(\nu_{2k-2,j}^{2}-\nu_{2k-2,r}^{2})}+\\\
&\frac{\prod\limits_{r=1}^{k-1}(-\theta_{k-1}\nu_{2k-3,r})\prod\limits_{r=1}^{k}(-\theta_{k-1}\nu_{2k-1,r})}{\prod\limits_{r=1}^{k-1}(-\frac{ic}{2}-\nu_{2k-2,r})(-\frac{ic}{2}+\nu_{2k-2,r})}\Bigg{)}w_{n}=\\\
&\frac{\theta_{k-1}}{ic}\left(\sum\limits^{k}_{j=1}\nu_{2k-1,j}+\sum\limits^{k-1}_{j=1}\nu_{2k-3,j}-\theta_{k-1}(k-1)ic\right).\end{split}$
(B.2)

Here we again use the identity (A.9) where for $x_{i}$ we substitute $2k-1$
variables $\\{\pm\nu_{2k-2,s},-\frac{ic}{2}\\}$, а and for $y_{i}$ use
$\\{\theta_{k-1}\nu_{2k-1}-\frac{ic}{2},\theta_{k-1}\nu_{2k-3}-\frac{ic}{2}\\}$.
$\Box$

## References

  * [D] Dixmier, Jacques,. Algebres enveloppantes. Vol. 192. Paris: Gauthier-Villars, 1974.
  * [GT] I.M.Gelfand, M.L. Tsetlin, Finite-dimensional representations of the group of unimodular matrices, DAN SSSR, 1950, v. 51, 825-828 (in Russian)
  * [GLO] Gerasimov, A. A., Lebedev, D. R., Oblezin, S. V. (2012). New integral representations of Whittaker functions for classical Lie groups. Russian Mathematical Surveys, 67(1), 1.
  * [GKL] A. Gerasimov, S. Kharchev, D. Lebedev, Scattering Method: The Open Toda Chain and the Hyperbolic Sutherland, International Mathematics Research Notices 2004.17 (2004): 823-854.
  * [Giv] A. Givental, “Stationary phase integrals, quantum Toda lattices, flag manifolds and the mirror conjecture”, Topics in singularity theory, Amer. Math. Soc. Transl. Ser. 2, 180, Amer. Math. Soc., Providence, RI, 1997, 103–115;
  * [Gust] R.A. Gustafson, Some q-Beta and Mellin-Barnes Integrals on Compact Lie Groups and Lie Algebras, Transactions of the American Mathematical Society, Vol. 341, No. 1. (1994), pp. 69-119.
  * [IS] Iorgov, N. Z., Shadura V.N., "Wave functions of the Toda chain with boundary interaction." Theoretical and mathematical physics 142 (2005): 289-305.
  * [KK] Kharchev, S., and S. Khoroshkin. "Mellin-Barnes presentations for Whittaker wave functions." Advances in Mathematics 375 (2020): 107368.
  * [KL] Kharchev, S., Lebedev, D. (2000). Eigenfunctions of $GL(N,R)$ Toda chain: Mellin-Barnes representation. Journal of Experimental and Theoretical Physics Letters, 71, 235-238.
  * [K] B. Kostant, Quantization and representation theory, In: Representation Theory of Lie Groups, Proc. of Symp., Oxford, 1977, pp. 287-317, London Math. Soc. Lecture Notes series, 34, Cambridge, 1979
  * [ZS] D.P. Zhelobenko, A.I. Stern, Representation of Lie groups(in Russian), "Science", 1983


# PREDICTION OF CELLULAR IDENTITIES FROM TRAJECTORY AND CELL FATE INFORMATION

###### Abstract

Determining cell identities in imaging sequences is an important yet
challenging task. The conventional method for cell identification is via cell
tracking, which is complex and can be time-consuming. In this study, we
propose an innovative approach to cell identification during early C. elegans
embryogenesis using machine learning. We employed random forest, MLP, and LSTM
models, and tested cell classification accuracy on 3D time-lapse confocal
datasets spanning the first $4$ hours of embryogenesis. By leveraging a small
number of spatial-temporal features of individual cells, including cell
trajectory and cell fate information, our models achieve an accuracy of over
$90\%$, even with limited data. We also determine the most important feature
contributions and can interpret these features in the context of biological
knowledge. Our research demonstrates the success of predicting cell identities
in 4D imaging sequences directly from simple spatio-temporal features.

Index Terms—  cell identification, machine learning, C. elegans embryogenesis,
classification, time-lapse imaging

##  1 Introduction

Cell identification is a crucial goal when imaging embryogenesis [1, 2, 3].
However, identifying cells in images, even in small, stereotyped organisms
like the nematode C. elegans [4], is challenging. The standard approach to
identify cells during development is via cell tracking [1, 5, 6]. Generally,
in cell tracking, the entire image sequence is processed, with cells detected
and then traced back to the very first frame of the image sequence. By
comparing the cell tracking result with the invariant lineage or other
biological prior, cells can be identified. This is a complex approach, which
involves several steps and can be time-consuming.

Fig. 1: An example of the original embryo image (left) and the image rotated
to canonical orientation (right).

Is it possible to recognize cells directly? On image object classification
tasks [7], researchers have demonstrated impressive results that outperformed
human-level performance [8]. Can we classify cells to identify them directly
as we classify objects in images? This is difficult as the appearances of
nuclei in fluorescence data are often very similar, yet classifying nuclear
identity may nevertheless be possible. Studies have shown that cells exhibit
stereotypical position and migration patterns during early C. elegans
embryogenesis [9], and the variability of cell position is highly
deterministic [10]. Here we propose a machine learning based cell
classification approach to identify cells in C. elegans embryos, based on a
small number of spatial-temporal features of individual cell nuclei, including
trajectory and cell fate information. We exploited random forest [11],
multilayer perceptron (MLP), and long short-term memory (LSTM) [12] as our
models for cell classification and tested our models with $28$ time-lapse
3-dimensional (4D) confocal imaging sequences [10]. All three models achieve
outstanding performances of over $90\%$ accuracy. Moreover, we justified our
models by examining feature contributions and linking feature importance to a
cell’s naming convention in the worm embryo.

##  2 Method

###  2.1 Data

The data used for this study comprise $28$ time-lapse 3-dimensional (4D)
confocal imaging sequences for wild-type C. elegans embryos, with a temporal
resolution of $1.25$ min per frame [10]. In the dataset, cell nuclei have
already been detected and tracked until approximately the $350$-cell stage in
each embryo. We rotated embryo images together with $(x,y,z)$ nuclear
positions into canonical orientation, by assigning anterior-posterior (A-P),
left-right (L-R), and dorsal-ventral (D-V) axes to X, Y, Z axes respectively.
Fig. 1 shows an example of the original and the rotated image.

Fig. 2: Trajectories (colored lines) of ABar cell in $28$ embryos.

###  2.2 Data preparation for learning

We cleaned the data, processed and extracted features for each cell for model
training and testing. There are $781$ unique cells in the $28$ embryos. We
filtered out the cells that only appear in some embryos, or do not have mother
or daughter cells, leaving $346$ cells remaining. The filtered cells are
mostly very early-stage or very late-stage cells. Out of these $346$ cells, we
further excluded $12$ cells with lifespans longer than $50$ frames, leaving us
a total of $334$ unique cells for our study.

The features we used in the models are the cell trajectory $(t,x,y,z)$, start
time, lifespan, division orientation to its mother cell, and the division
orientations of its two daughter cells. We extracted and calculated these
features from the cell tracking results $(cell\ label,t,x,y,z)$ associated
with the images. For cells with lifespans shorter than $50$ frames, we pad
negative $10000$ to the end of their trajectory $(t,x,y,z)$. The start time of
each cell is its birth time, in units of frames, with respect to the end of
$4$-cell stage. The division orientation is a 3D unit vector representing a
cell’s birth location relative to its mother cell. Therefore, the feature
dimension is $211$ (the first $200$ being trajectory features) if trajectory
features are of $(t,x,y,z)$ format, or $161$ (the first $150$ being trajectory
features) if trajectory features are of $(x,y,z)$ format. The trajectories of
the ABar cell across all $28$ embryos are shown in Fig. 2.

Fig. 3: Model architecture for MLP and LSTM. The dashed line depicts whether
extra features, like start time, lifespan, and division orientations, are used
in the model.

###  2.3 Models

For our 4D cell classification task, we consider 3 models, random forest [11],
multilayer perceptron (MLP), and long short-term memory (LSTM) [12]. Our MLP
and LSTM models consist of two layers. We use a single layer MLP (with ReLU
activation and residual connection [8]) or single layer LSTM to learn patterns
from trajectory features. After learning these patterns, they are concatenated
with other features and fed into a linear layer. The linear layer then uses a
softmax function to output the predicted cell label. The model architecture
for MLP and LSTM is illustrated in Fig. 3. For MLP and random forest, the
trajectory feature input is $(t,x,y,z)$. For LSTM, we designed two models: one
accepts $(x,y,z)$ format trajectory features, which we denote as LSTM; the
other accepts $(t,x,y,z)$ format trajectory features, which we denote as
LSTMt.

###  2.4 Cross-validation

We randomly sampled $4$ embryo datasets for model testing and left the
remaining $24$ embryo samples for $6$-fold cross-validation. We use cross-
validation to find the optimal hyperparameters for the models. The cross-
validation results are listed in Table 1. Here we performed a feature ablation
study by deleting some features, and inspected model performances. For the
random forest model, we examined the number of trees from $10$ to $200$ in
cross-validation and picked a random forest with $169$ trees according to
validation accuracy. For MLP, LSTM, and LSTMt, we searched learning rate,
weight_decay ($L_{2}$ regularization), batch size, dropout rate [13], and
number of hidden layers. We use AdamW [14] as the optimization method to
employ $L_{2}$ regularization, and an exponential learning rate scheduler to
decay the learning rate by $0.999$ every epoch. All models are trained for
$3000$ epochs. The optimal weight decay is $0.1$ for MLP, $0.15$ for LSTM, and
$0.05$ for LSTMt. For MLP, an optimal dropout rate of $0.5$ is utilized. For
all three models, the optimal learning rate is $0.001$ and the optimal batch
size is $128$. All three models perform best with a single hidden layer,
probably due to the small size of our dataset.

From the feature ablation study result in Table 1, we can see models with full
features produce the best performances. One notable point is that, when the
feature corresponding to division orientation relative to the mother cell is
added, the validation accuracy exhibits a large increase compared to when
adding other features.

---
|  | Random
---
Forest
MLP | LSTM | LSTMt
Traj | 0.885 | 0.813 | 0.785 | 0.817
Traj + SF | 0.887 | 0.813 | 0.827 | 0.844
Traj + SF + LF | 0.887 | 0.816 | 0.831 | 0.845
Traj + SF + LF + DM | 0.918 | 0.870 | 0.899 | 0.901
Full | 0.921 | 0.879 | 0.907 | 0.898

Table 1: Cross-validation results. Validation accuracies above $90\%$ are
highlighted. Here Traj denotes trajectory, SF denotes start time, LF denotes
lifespan, DM denotes division orientation to mother cell, Full denotes all
features included.

##  3 Results

###  3.1 Test accuracy

---
|  | Random
---
Forest
MLP | LSTM | LSTMt
Traj | 0.897 | 0.850 | 0.836 | 0.859
Traj + SF | 0.892 | 0.852 | 0.856 | 0.884
Traj + SF + LF | 0.894 | 0.844 | 0.841 | 0.882
Traj + SF + LF + DM | 0.926 | 0.907 | 0.924 | 0.924
Full | 0.932 | 0.912 | 0.923 | 0.932

Table 2: Test results. Test accuracies above $90\%$ are highlighted. Here Traj
denotes trajectory, SF denotes start time, LF denotes lifespan, DM denotes
division orientation to mother cell, Full denotes all features included.

To evaluate the test performance of the models, we trained them on all $24$
cross-validation embryo samples using the optimal hyperparameters found in the
cross-validation step (Section 2.4). The models were then tested on $4$ test
embryo samples. Table 2 shows the test accuracy with feature ablations. It is
evident that all four models achieved remarkable accuracy when using full
features. In fact, the random forest and LSTMt models even achieved above
$93\%$ accuracy. Furthermore, by looking at different models, we can see
clearly that LSTM/LSTMt and random forest show superior performances over MLP.
This is because LSTM excels at capturing temporal dependencies, and the small
size of our dataset, which requires simple models to avoid overfitting, makes
the traditional machine learning random forest model well suited to the task.

Considering that we only have $28$ embryo image sequences to work with, and
there are $334$ unique cells to classify, the excellent cell classification
accuracy is a surprising result. The simplicity of our models together with
their superb performances imply that cells can be identified correctly with
machine learning based classification methods based solely on spatial-temporal
features of individual cells, including trajectory, start time, lifespan, and
division orientations.

###  3.2 Feature importance

Fig. 4: Top-$10$ features in random forest ranked by feature importance. The
blue error bars show standard deviations of feature importance for each
feature.

Table 2 also illustrates the feature contributions, similar to the cross-
validation result (Table 1). From the table, we can see that the highest
accuracy for each model is achieved when using all the features. A significant
increase in accuracy occurs when adding a cell’s division orientation to its
mother cell, and this is consistent across all the models. To investigate why
this happens, we ranked the feature importance in random forest [11], and the
top-$10$ features are displayed in Fig. 4. From the figure, it is obvious that
feature $202$ (the first $200$ features are the $50$-frame trajectory, feature
$200$ is the start time, feature $201$ is the lifespan), which is cell’s
division orientation to its mother cell along X-axis, is the most critical
feature. This feature indicates whether a cell is born towards the anterior or
posterior (A or P) direction. The result is reasonable as most cells in worm
embryos undergo A-P divisions (these cells obtain their names by appending ‘a’
or ‘p’ at the end of their mother cells’ names) [4]. Meanwhile, we can see
that the remaining top-$10$ features are all trajectory features, which
explains why models that use trajectory features alone, can achieve an
accuracy of over $85\%$ (see Table 2). This implies that it is possible to
design a high-performance cell classification model by only using cell
trajectories.

##  4 CONCLUSIONS AND DISCUSSION

In summary, we proposed the first machine learning based cell classification
approach to identify cells in C. elegans embryos by utilizing a cell nuclei’s
spatiotemporal features. We demonstrated the power of our approach with four
machine learning models: random forest, MLP, LSTM, and LSTMt, all of which
achieved outstanding performance of over $90\%$ accuracy. Furthermore, we
investigated the feature importance in our models and found that the cell’s
division orientation to its mother cell is the most crucial feature, which
accords with cell’s naming conventions in worm embryos. Moreover, our result
showed that even when using trajectory features alone, cells can be classified
at a high accuracy of $85\%$.

Our approach necessitates the incorporation of trajectory features, which
means we need to detect and track accurately the cell of interest, to acquire
its trajectory $(t,x,y,z)$. In contrast to cell tracking approaches, our
approach only requires tracking the cell of interest and does not need to
track and trace the cell all the way back to its root/ancestor cell at the
very first frame of the image sequence, as in cell tracking. Further, our
approach can output the cell name directly, whereas, in cell tracking, an
additional step of comparing cell tracking results with the lineage of C.
elegans is required to get the cell name. Given these advantages and the
simplicity of our approach, we believe our approach will be broadly useful in
cell identification, especially when there is a need to identify individual
cells of interest. This happens when cells are labeled a fluorescent marker,
such as in nematode neuronal or muscle cells.

##  5 Compliance with ethical standards

No ethical approval was required for this study.

##  6 CODE AVAILABILITY

The code is publicly available at: https://github.com/pmcesky/4D-Cell.

##  7 Acknowledgments

This work is supported by the Gordon and Betty Moore Foundation under grant
number GBMF8815. This work was supported by the Howard Hughes Medical
Institute. This article is subject to HHMI’s Open Access to Publications
policy. HHMI lab heads have previously granted a nonexclusive CC BY 4.0
license to the public and a sublicensable license to HHMI in their research
articles. Pursuant to those licenses, the author-accepted manuscript of this
article can be made freely available under a CC BY 4.0 license immediately
upon publication.

## References

  * [1] Zhirong Bao, John I Murray, Thomas Boyle, Siew Loon Ooi, Matthew J Sandel, and Robert H Waterston,  “Automated cell lineage tracing in caenorhabditis elegans,”  Proceedings of the National Academy of Sciences, vol. 103, no. 8, pp. 2707–2712, 2006.
  * [2] Yicong Wu, Peter Wawrzusin, Justin Senseney, Robert S Fischer, Ryan Christensen, Anthony Santella, Andrew G York, Peter W Winter, Clare M Waterman, Zhirong Bao, et al.,  “Spatially isotropic four-dimensional imaging with dual-view plane illumination microscopy,”  Nature biotechnology, vol. 31, no. 11, pp. 1032–1038, 2013.
  * [3] Mark W Moyle, Kristopher M Barnes, Manik Kuchroo, Alex Gonopolskiy, Leighton H Duncan, Titas Sengupta, Lin Shao, Min Guo, Anthony Santella, Ryan Christensen, et al.,  “Structural and developmental principles of neuropil assembly in c. elegans,”  Nature, vol. 591, no. 7848, pp. 99–104, 2021.
  * [4] John E Sulston, Einhard Schierenberg, John G White, and J Nichol Thomson,  “The embryonic cell lineage of the nematode caenorhabditis elegans,”  Developmental biology, vol. 100, no. 1, pp. 64–119, 1983.
  * [5] Thomas J Boyle, Zhirong Bao, John I Murray, Carlos L Araya, and Robert H Waterston,  “Acetree: a tool for visual analysis of caenorhabditis elegans embryogenesis,”  BMC bioinformatics, vol. 7, no. 1, pp. 1–14, 2006.
  * [6] Martin Maška, Vladimír Ulman, Pablo Delgado-Rodriguez, Estibaliz Gómez-de Mariscal, Tereza Nečasová, Fidel A Guerrero Peña, Tsang Ing Ren, Elliot M Meyerowitz, Tim Scherr, Katharina Löffler, et al.,  “The cell tracking challenge: 10 years of objective benchmarking,”  Nature Methods, pp. 1–11, 2023.
  * [7] Olga Russakovsky, Jia Deng, Hao Su, Jonathan Krause, Sanjeev Satheesh, Sean Ma, Zhiheng Huang, Andrej Karpathy, Aditya Khosla, Michael Bernstein, et al.,  “Imagenet large scale visual recognition challenge,”  International journal of computer vision, vol. 115, pp. 211–252, 2015.
  * [8] Kaiming He, Xiangyu Zhang, Shaoqing Ren, and Jian Sun,  “Deep residual learning for image recognition,”  in Proceedings of the IEEE conference on computer vision and pattern recognition, 2016, pp. 770–778.
  * [9] Guoye Guan, Ming-Kin Wong, Vincy Wing Sze Ho, Xiaomeng An, Lu-Yan Chan, Binghui Tian, Zhiyuan Li, Leihan Tang, Zhongying Zhao, and Chao Tang,  “System-level quantification and phenotyping of early embryonic morphogenesis of caenorhabditis elegans,”  bioRxiv, 2019.
  * [10] Xiaoyu Li, Zhiguang Zhao, Weina Xu, Rong Fan, Long Xiao, Xuehua Ma, and Zhuo Du,  “Systems properties and spatiotemporal regulation of cell position variability during embryogenesis,”  Cell reports, vol. 26, no. 2, pp. 313–321, 2019.
  * [11] Leo Breiman,  “Random forests,”  Machine learning, vol. 45, pp. 5–32, 2001.
  * [12] Sepp Hochreiter and Jürgen Schmidhuber,  “Long short-term memory,”  Neural computation, vol. 9, no. 8, pp. 1735–1780, 1997.
  * [13] Nitish Srivastava, Geoffrey Hinton, Alex Krizhevsky, Ilya Sutskever, and Ruslan Salakhutdinov,  “Dropout: a simple way to prevent neural networks from overfitting,”  The journal of machine learning research, vol. 15, no. 1, pp. 1929–1958, 2014.
  * [14] Ilya Loshchilov and Frank Hutter,  “Decoupled weight decay regularization,”  in International Conference on Learning Representations, 2019.