Rahman polynomials

Ryu Sasaki

Department of Physics and Astronomy, Tokyo University of Science, Noda
278-8510, Japan

###### Abstract

The Rahman polynomials are constructed explicitly as the left eigenvectors of
a certain multi-dimensional discrete time Markov chain operator
$K_{n}({\boldsymbol{x}},{\boldsymbol{y}},N)$ which is a convolution of an
$n+1$-nomial distribution $W_{n}({\boldsymbol{x}},N)$ and an $n$-tuple of
binomial distributions $\prod_{i}W_{1}(x,N)$ as
$K_{n}({\boldsymbol{x}},{\boldsymbol{y}},N)=\sum_{\boldsymbol{z}}W_{n}({\boldsymbol{x}}-{\boldsymbol{z}},N-\sum_{i}z_{i})\prod_{i}W_{1}(z_{i},y_{i})$.
The Markov chain was introduced and discussed by Hoare, Rahman and Grünbaum as
a multivariable version of the known soluble single variable one. The
anticipated solubility of the model gave Rahman polynomials the prospect of
the first multivariate hypergeometric function of Aomoto-Gelfand type
connected with solvable dynamics. The promise is now realised. The $n^{2}$
system parameters $\\{u_{i\,j}\\}$ of the Rahman polynomials are completely
determined. These $u_{i\,j}$’s are irrational functions of the original system
parameters, the probabilities of the multinomial and binomial distributions.

##  1 Introduction

The Rahman polynomials, a typical example of multivariate orthogonal
polynomials, [3, 21, 22, 18, 19, 20, 24, 25], are constructed explicitly. The
method of construction is very different from the existing theories and models
of multivariate orthogonal polynomials, see for example, [14, 15, 5, 42, 44,
6, 33, 34, 22, 18, 23, 19, 24, 34, 20, 25, 31, 4, 43, 9, 10, 11, 12, 41, 16,
17, 13, 32]. The Rahman polynomials and two multivariate orthogonal
polynomials constructed recently, the Krawtchouk [27, 39] and Meixner
polynomials [40] are a certain type of Aomoto-Gelfand [1, 8] terminating
hypergeometric functions, see (2.19). Like other orthogonal polynomials of a
discrete variable [2, 28, 30, 36], the knowledge of $n$ degree one polynomials
is enough to determine the system parameters $\\{u_{i\,j}\\}$ (3.9)
completely. The entire polynomials are expressed by $u_{i\,j}$’s polynomially,
see formula (2.19). A certain dynamical setting, a quite essential one,
determines the degree one polynomials. For the Rahman polynomials, it is the
discrete time Markov chain and for the multivariate Krawtchouk and Meixner
polynomials, it is the birth and death processes [7, 38, 39, 40]. These
dynamical settings determine the orthogonal weight functions at the same time.
But the polynomials are determined independently of the orthogonality weights.

The construction of the Rahman polynomials is straightforward. The problem
setting was done completely by the works of Hoare, Rahman and Grünbaum [22,
18, 19, 20]. I only have to calculate the degree one left eigenpolynomials of
the Markov chain operator
$K_{n}({\boldsymbol{x}},{\boldsymbol{y}};N,{\boldsymbol{\alpha}},{\boldsymbol{\beta}})$
and determine the system parameters $\\{u_{i\,j}\\}$ as a function of $2n$
system parameters ${\boldsymbol{\alpha}}$ and ${\boldsymbol{\beta}}$. The
remaining jobs are to determine the entire spectrum
$\mathcal{E}({\boldsymbol{m}})=\prod_{i=1}^{n}\lambda_{i}^{m_{i}}$ (3.15), in
which $\lambda_{i}$’s are the eigenvalues of the degree one polynomials. The
final task is to demonstrate that the entire polynomials are the left
eigenvectors of the Markov chain operator with the multiplicative spectrum.
The generating function by Mizukawa [34] plays the central role.

This paper is organised as follows. The problem setting of discrete time
Markov chain due to Hoare, Rahman and Grünbaum is briefly recapitulated in
section 2. The important role played by left eigenvectors are stated in
Theorem 2.3. The main players, Aomoto-Gelfand [1, 8] terminating
hypergeometric functions and the generating function are introduced in Theorem
2.4. The Rahman polynomials are constructed in section 3. Starting with degree
one left eigenpolynomials in §3.1, the spectrum is determined in §3.2 and the
proof of the left eigenvector equations for
$P_{\boldsymbol{m}}({\boldsymbol{x}};{\boldsymbol{u}})$ is provided in §3.3.
The final section is for some comments. This paper is prepared plainly so that
non-experts can understand.

##  2 Problem setting

Let us start with a brief review of the general framework of discrete time
Markov chains of $n$-variables. For simplicity of presentation, positive
transition matrices on finite sets only are discussed.

###  2.1 Discrete time Markov chains of $n$-variables

The main actor is the transition matrix
$\mathcal{K}_{n}({\boldsymbol{x}},{\boldsymbol{y}})$ which governs the
evolution of the probability distribution function
$\mathcal{P}({\boldsymbol{x}};\ell)$ over $\Xi$, a finite set,

$\mathcal{P}({\boldsymbol{x}};\ell)\longrightarrow\mathcal{P}({\boldsymbol{x}};\ell+1)=\sum_{{\boldsymbol{y}}\in\Xi}\mathcal{K}_{n}({\boldsymbol{x}},{\boldsymbol{y}})\mathcal{P}({\boldsymbol{y}};\ell),\quad\mathcal{P}({\boldsymbol{x}};\ell)\geq
0,\quad\sum_{{\boldsymbol{x}}\in\Xi}\mathcal{P}({\boldsymbol{x}};\ell)=1,\
\ell\in\mathbb{N}_{0},$

at time step $\ell$ to $\ell+1$. That is,
$\mathcal{K}_{n}({\boldsymbol{x}},{\boldsymbol{y}})$ specifies the one step
transition probability from point ${\boldsymbol{y}}$ to ${\boldsymbol{x}}$
which belong to a finite set $\Xi$, a subset of the non-negative integer
lattice in $n$-dimensions $\mathbb{N}_{0}^{n}$,

$\displaystyle\mathcal{K}_{n}({\boldsymbol{x}},{\boldsymbol{y}})>0,\quad{\boldsymbol{x}}=(x_{1},\ldots,x_{n}),\quad{\boldsymbol{y}}=(y_{1},\ldots,y_{n}),\quad\qquad{\bf
x},{\bf y}\in\Xi\subseteq\mathbb{N}_{0}^{n},$

satisfying

$\displaystyle\text{\bf conservation of
probability:}\quad\sum_{{\boldsymbol{x}}\in\Xi}\mathcal{K}_{n}({\boldsymbol{x}},{\boldsymbol{y}})=1,$
(2.1)
$\displaystyle\qquad\qquad\Longrightarrow\sum_{{\boldsymbol{x}}\in\Xi}\mathcal{P}({\boldsymbol{x}};\ell+1)=\sum_{{\boldsymbol{x}}\in\Xi}\mathcal{P}({\boldsymbol{x}};\ell)=1.$

Let us impose the condition that $\mathcal{P}({\boldsymbol{x}};\ell)$ does not
explode and approaches to a stationary distribution $\pi_{S}$ at large $\ell$;

$\lim_{\ell\to\infty}\mathcal{P}({\boldsymbol{x}};\ell)=\pi_{S}({\boldsymbol{x}}).$

The existence of the stationary distribution is guaranteed when
$\mathcal{K}_{n}({\boldsymbol{x}},{\boldsymbol{y}})$ satisfies the following
reversibility or detailed balance condition with a certain positive reversible
distribution $\pi_{R}({\boldsymbol{x}})$,

reversible:
$\displaystyle\mathcal{K}_{n}({\boldsymbol{x}},{\boldsymbol{y}})\pi_{R}({\boldsymbol{y}})=\mathcal{K}_{n}({\boldsymbol{y}},{\boldsymbol{x}})\pi_{R}({\boldsymbol{x}}),\quad\pi_{R}({\boldsymbol{x}})>0,\quad{\boldsymbol{x}},{\boldsymbol{y}}\in\Xi,$
(2.2)

Deviding the above reversibility condition by
$\sqrt{\pi_{R}({\boldsymbol{x}})\pi_{R}({\boldsymbol{y}})}$, one obtains a
real symmetric matrix $\mathcal{T}_{n}$;

$\displaystyle\mathcal{T}_{n}({\boldsymbol{x}},{\boldsymbol{y}})$
$\displaystyle\stackrel{{\scriptstyle\text{def}}}{{=}}\frac{1}{\sqrt{\pi_{R}({\boldsymbol{x}})}}\,\mathcal{K}_{n}({\boldsymbol{x}},{\boldsymbol{y}})\sqrt{\pi_{R}({\boldsymbol{y}})}=\frac{1}{\sqrt{\pi_{R}({\boldsymbol{y}})}}\,\mathcal{K}_{n}({\boldsymbol{y}},{\boldsymbol{x}})\sqrt{\pi_{R}({\boldsymbol{x}})}$
(2.3)
$\displaystyle=\mathcal{T}_{n}({\boldsymbol{y}},{\boldsymbol{x}}),\qquad{\boldsymbol{x}},{\boldsymbol{y}}\in\Xi.$
(2.4)

In other words, the transition matrix
$\mathcal{K}_{n}({\boldsymbol{x}},{\boldsymbol{y}})$ is related to a real
positive symmetric matrix $\mathcal{T}_{n}({\boldsymbol{x}},{\boldsymbol{y}})$
by a similarity transformation in terms of $\sqrt{\pi_{R}({\boldsymbol{x}})}$.
This means the following

###### Proposition 2.1

The eigenvalues of $\mathcal{K}_{n}$ are real and the range of the spectrum is

$-1<\text{The eigenvalues of
}\mathcal{K}({\boldsymbol{x}},{\boldsymbol{y}})\leq 1,$ (2.5)

and the eigenvector of the maximal eigenvalue 1 is the the reversible
distribution

$\displaystyle\sum_{{\boldsymbol{y}}\in\Xi}\mathcal{K}_{n}({\boldsymbol{x}},{\boldsymbol{y}})\pi_{R}({\boldsymbol{y}})=\pi_{R}({\boldsymbol{x}})\sum_{{\boldsymbol{y}}\in\Xi}\mathcal{K}_{n}({\boldsymbol{y}},{\boldsymbol{x}})=\pi_{R}({\boldsymbol{x}}).$
(2.6)

The maximal eigenvalue is simple. These are due to Perron’s theorem on
positive matrices. As $\ell$ becomes large, all the components of
$\mathcal{P}({\boldsymbol{x}};\ell)$ belonging to the eigenvalues other than 1
diminish and

$\lim_{\ell\to\infty}\mathcal{P}({\boldsymbol{x}};\ell)=\pi_{R}({\boldsymbol{x}}).$

In other words, the reversible distribution is stationary.

The above maximal eigenvector relation (2.6) also means the following

###### Remark 2.2

A constant is the left eigenvector of
$\mathcal{K}_{n}({\boldsymbol{x}},{\boldsymbol{y}})$ with the maximal
eigenvalue 1;

$v_{\boldsymbol{0}}({\boldsymbol{x}})\equiv
1,\quad\sum_{{\boldsymbol{x}}\in\Xi}\mathcal{K}_{n}({\boldsymbol{x}},{\boldsymbol{y}})v_{\boldsymbol{0}}({\boldsymbol{x}})=\sum_{{\boldsymbol{x}}\in\Xi}\mathcal{K}_{n}({\boldsymbol{x}},{\boldsymbol{y}})=1=v_{\boldsymbol{0}}({\boldsymbol{y}}).$

By generalising the Remark 2.2, one arrives at the following

###### Theorem 2.3

The complete set of left eigenvectors
$\\{v_{\boldsymbol{m}}({\boldsymbol{x}})\\}$ of
$\mathcal{K}_{n}({\boldsymbol{x}},{\boldsymbol{y}})$,

$\sum_{{\boldsymbol{x}}\in\Xi}\mathcal{K}_{n}({\boldsymbol{x}},{\boldsymbol{y}})v_{\boldsymbol{m}}({\boldsymbol{x}})=\kappa({\boldsymbol{m}})v_{\boldsymbol{m}}({\boldsymbol{x}}),\quad-1<\kappa({\boldsymbol{m}})\leq
1=\kappa({\boldsymbol{0}}),\quad\forall{\boldsymbol{m}}\in\Xi^{d},$ (2.7)

provides the complete sets of eigenvectors of
$\mathcal{K}_{n}({\boldsymbol{x}},{\boldsymbol{y}})$ and
$\mathcal{T}_{n}({\boldsymbol{x}},{\boldsymbol{y}})$;

$\displaystyle\sum_{{\boldsymbol{y}}\in\Xi}\mathcal{K}_{n}({\boldsymbol{x}},{\boldsymbol{y}})v_{\boldsymbol{m}}({\boldsymbol{y}})\pi_{R}({\boldsymbol{y}})=\kappa({\boldsymbol{m}})v_{\boldsymbol{m}}({\boldsymbol{x}})\pi_{R}({\boldsymbol{x}}),\quad\forall{\boldsymbol{m}}\in\Xi^{d},$
(2.8)
$\displaystyle\sum_{{\boldsymbol{y}}\in\Xi}\mathcal{T}_{n}({\boldsymbol{x}},{\boldsymbol{y}})v_{\boldsymbol{m}}({\boldsymbol{y}})\sqrt{\pi_{R}({\boldsymbol{y}})}=\kappa({\boldsymbol{m}})v_{\boldsymbol{m}}({\boldsymbol{x}})\sqrt{\pi_{R}({\boldsymbol{x}})},\quad\forall{\boldsymbol{m}}\in\Xi^{d}.$
(2.9)

Here $\Xi^{d}$ is the dual set of $\Xi$ and $\\#(\Xi^{d})=\\#(\Xi)$. The left
eigenvectors are orthogonal to each other with the orthogonality weight
function $\pi_{R}({\boldsymbol{x}})$;

$\sum_{{\boldsymbol{x}}\in\Xi}v_{\boldsymbol{m}}({\boldsymbol{x}})v_{{\boldsymbol{m}}^{\prime}}({\boldsymbol{x}})\pi_{R}({\boldsymbol{x}})=0,\quad{\boldsymbol{m}}\neq{\boldsymbol{m}}^{\prime},$
(2.10)

since $\mathcal{K}_{n}({\boldsymbol{x}},{\boldsymbol{y}})$ has some parameters
and for generic values of the parameters, the eigenvalues are not degenerate.

This setup of Markov chains in one dimensions has been discussed by many
authors [3, 21, 37], in connection with “cumulative Bernoulli trials.” The
multivariate version was initiated by Hoare, Rahaman and Grünbaum [22, 18, 19,
20] in a quest for a new multivariate version of the Krawtchouk polynomials,
called Rahman polynomials. In particular, the explicit setup in the next
subsection was due to Grünbaum and Rahman [20], and this paper is cite as I.

###  2.2 Quest for Rahman polynomials

The purpose of this paper is to present the explicit forms of the Rahman
polynomials by solving an $n$-variable Markov chain with the transition matrix
(operator) $K$ composed of the convolution of an $n+1$-nomial distribution
$W_{n}$ and $n$ binomial distributions $W_{1}$;

$\displaystyle
K_{n}({\boldsymbol{x}},{\boldsymbol{y}};N,{\boldsymbol{\alpha}},{\boldsymbol{\beta}})\stackrel{{\scriptstyle\text{def}}}{{=}}\sum_{{\boldsymbol{z}}\in\mathcal{X}}W_{n}({\boldsymbol{x}}-{\boldsymbol{z}};N-|z|,{\boldsymbol{\beta}})\prod_{i=1}^{n}W_{1}(z_{i};y_{i},\alpha_{i})>0,$
(2.11)

in which

$\displaystyle{\boldsymbol{x}}=(x_{1},\ldots,x_{n}),\quad{\boldsymbol{y}}=(y_{1},\ldots,y_{n}),\quad{\boldsymbol{z}}=(z_{1},\ldots,z_{n}),\quad|x|\stackrel{{\scriptstyle\text{def}}}{{=}}\sum_{i=1}^{n}x_{i},\quad
N\in\mathbb{N},$ $\displaystyle\qquad{\bf x},{\bf y},{\bf
z}\in\mathcal{X}\stackrel{{\scriptstyle\text{def}}}{{=}}\\{\boldsymbol{x}\in\mathbb{N}_{0}^{n}\
|\,|x|\leq N\\},\qquad\quad\\#\bigl{(}\mathcal{X}\bigr{)}=\binom{N+n}{n},$
$\displaystyle{\boldsymbol{\alpha}}=(\alpha_{1},\ldots,\alpha_{n}),\
{\boldsymbol{\beta}}=(\beta_{1},\ldots,\beta_{n}),\
0<\alpha_{i},\beta_{i}<1,\quad i=1,\ldots,n,\quad 0<|\beta|<1,$

and the binomial $W_{1}$ and multinomial $W_{n}$ distributions are,

$\displaystyle
W_{1}(x;y,\alpha)\stackrel{{\scriptstyle\text{def}}}{{=}}\binom{y}{x}\alpha^{x}(1-\alpha)^{y-x}>0,\quad\sum_{x=0}^{y}W_{1}(x;y,\alpha)=1,$
(2.12) $\displaystyle
W_{n}({\boldsymbol{x}};N,{\boldsymbol{\beta}})\stackrel{{\scriptstyle\text{def}}}{{=}}\frac{N!\cdot(1-|\beta|)^{N-|x|}}{x_{1}!\cdots
x_{n}!(N-|x|)!}\cdot\prod_{i=1}^{n}\beta_{i}^{x_{i}}=\binom{N}{\boldsymbol{x}}\beta_{0}^{x_{0}}\boldsymbol{\beta}^{\boldsymbol{x}}>0,$
(2.13)
$\displaystyle\sum_{{\boldsymbol{x}}\in\mathcal{X}}W_{n}({\boldsymbol{x}};N,{\boldsymbol{\beta}})=1,\quad
x_{0}\stackrel{{\scriptstyle\text{def}}}{{=}}N-|x|,\
\beta_{0}\stackrel{{\scriptstyle\text{def}}}{{=}}1-|\beta|,\
\binom{N}{\boldsymbol{x}}\stackrel{{\scriptstyle\text{def}}}{{=}}\frac{N!}{x_{1}!\cdots
x_{n}!x_{0}!}.$

To be more explicit, the summation in (2.11) means $0\leq
z_{i}\leq{\rm{min}}[x_{i},y_{i}]$, $i=1,\ldots,n$. The conservation of
probability (2.1),
$\sum_{{\boldsymbol{x}}\in\mathcal{X}}K_{n}({\boldsymbol{x}},{\boldsymbol{y}};N,{\boldsymbol{\alpha}},{\boldsymbol{\beta}})=1$,
holds as the double sum
$\sum_{{\boldsymbol{x}},{\boldsymbol{z}}\in\mathcal{X}}$ can be changed.
Summing over $\boldsymbol{x}$ for fixed $\boldsymbol{z}$ gives
$\sum_{{\boldsymbol{x}}\in{\mathcal{X}}}W_{n}({\boldsymbol{x}}-{\boldsymbol{z}};N-|z|,{\boldsymbol{\beta}})=1$.
The remaining $\boldsymbol{z}$ sum over $W_{1}$’s gives trivially 1.

The reversible distribution was identified in I [20]. It is the same
multinomial distribution as $W_{n}({\boldsymbol{x}};N,{\boldsymbol{\beta}})$
but with a different set of probability parameters ${\boldsymbol{\eta}}$,

$\displaystyle
K_{n}({\boldsymbol{x}},{\boldsymbol{y}};N,{\boldsymbol{\alpha}},{\boldsymbol{\beta}})W_{n}({\boldsymbol{y}};N,{\boldsymbol{\eta}})=K_{n}({\boldsymbol{x}},{\boldsymbol{y}};N,{\boldsymbol{\alpha}},{\boldsymbol{\beta}})W_{n}({\boldsymbol{x}};N,{\boldsymbol{\eta}}),$
(2.14)
$\displaystyle\eta_{i}=\frac{\beta_{i}}{1-\alpha_{i}}\frac{1}{D_{n}},\quad
1=1,\ldots,n,\qquad
D_{n}\stackrel{{\scriptstyle\text{def}}}{{=}}1+\sum_{k=1}^{n}\frac{\alpha_{k}\beta_{k}}{1-\alpha_{k}},$
(2.15)
$\displaystyle\frac{1-\alpha_{1}}{\beta_{1}}\eta_{1}=\frac{1-\alpha_{2}}{\beta_{2}}\eta_{2}=\cdots=\frac{1-\alpha_{n}}{\beta_{n}}\eta_{n}=\frac{1-|\eta|}{1-|\beta|}=1-\sum_{k=1}^{n}\alpha_{k}\eta_{k}=D_{n}^{-1}.$
(I.1.6)

I will not repeat the derivation of these complicated formulas (I.1.6).
Instead I show a path to (2.15) and (I.1.6) by using an equivalent expression
of the multinomial distribution (2.13), which derives from the multinomial
expansion

$\left(1+\sum_{i=1}^{n}q_{i}\right)^{N}=\sum_{{\boldsymbol{x}}\in\mathcal{X}}\binom{N}{{\boldsymbol{x}}}\prod_{i=1}^{n}q_{i}^{xi},\quad
q_{i}\in\mathbb{C},\quad i=1,\ldots,n.$

This leads to another form of the multinomial distribution

$\displaystyle\bar{W}_{n}({\boldsymbol{x}};N,{\boldsymbol{q}})\stackrel{{\scriptstyle\text{def}}}{{=}}\binom{N}{{\boldsymbol{x}}}\prod_{i=1}^{n}q_{i}^{x_{i}}(1+|q|)^{-N}$
$\displaystyle\Leftrightarrow
W_{n}({\boldsymbol{x}};N,{\boldsymbol{\beta}})=\binom{N}{{\boldsymbol{x}}}\prod_{i=1}^{n}\beta_{i}^{x_{i}}(1-|\beta|)^{N-|x|},$
$\displaystyle 0<q_{i}=\frac{\beta_{i}}{1-|\beta|}$
$\displaystyle\Leftrightarrow\beta_{i}=\frac{q_{i}}{1+|q|},\quad
i=1,\ldots,n.$ (2.16)

This type of multinomial distributions was useful for the derivation of the
multivariate Krawtchouk and Meixner polynomials based on the Birth and Death
formalism [39, 40]. In terms of these, the left hand side of the reversibility
condition (2.14) reads

$\displaystyle\sum_{{\boldsymbol{z}}\in\mathcal{x}}\bar{W}_{n}({\boldsymbol{x}}-{\boldsymbol{z}};N-|z|,{\boldsymbol{q}})\prod_{i=1}^{n}\bar{W}_{1}(z_{i};y_{i},p_{i})\bar{W}_{n}({\boldsymbol{y}};N,{\boldsymbol{r}})$
$\displaystyle=\prod_{i=1}^{n}\left[q_{i}^{x_{i}}\left(\frac{r_{i}}{1+p_{i}}\right)^{y_{i}}\right]\cdot\sum_{\boldsymbol{z}}\frac{N!(N-|z|)!\left((1+|q|)(1+|r|)\right)^{-N}}{(N-|x|)!(N-|y|)!\prod_{i=1}^{n}z_{i}!(x_{i}-z_{i})!(y_{i}-z_{i})!}.$

The reversibility requires that this expression is
${\boldsymbol{x}}\leftrightarrow{\boldsymbol{y}}$ symmetric. It is necessary
and sufficient that the unknown parameters $\\{r_{i}\\}$ should be

$r_{i}=q_{i}(1+p_{i}),\quad i=1,\ldots,n,$ (2.17)

which correspond to the known parameters $\\{\eta_{i}\\}$ in (2.15) by the
following translation rule together with (2.16),

$r_{i}=\frac{\eta_{i}}{1-|\eta|}\Leftrightarrow\eta_{i}=\frac{r_{i}}{1+|r|},\quad
p_{i}=\frac{\alpha_{i}}{1-\alpha_{i}}\Leftrightarrow\alpha_{i}=\frac{p_{i}}{1+p_{i}}\quad
i=1,\ldots,n.$ (2.18)

This might give some idea how the detailed balance condition works.

The determination of the reversible distribution (2.15) leads to the following

###### Theorem 2.4

The multivariate polynomials orthogonal with respect to the reversible
distribution (2.15) are $(n+1,2n+2)$ type terminating hypergeometric function
of Aomoto-Gelfand [1, 8, 34],

$\displaystyle
P_{\boldsymbol{m}}(\boldsymbol{x};{\boldsymbol{u}})\stackrel{{\scriptstyle\text{def}}}{{=}}\sum_{\begin{subarray}{c}\sum_{i,j}c_{ij}\leq
N\\\ (c_{ij})\in
M_{n}({\mathbb{N}_{0}})\end{subarray}}\frac{\prod\limits_{i=1}^{n}(-x_{i})_{\sum\limits_{j=1}^{n}c_{ij}}\prod\limits_{j=1}^{n}(-m_{j})_{\sum\limits_{i=1}^{n}c_{ij}}}{(-N)_{\sum_{i,j}c_{ij}}}\;\frac{\prod(u_{ij})^{c_{ij}}}{\prod
c_{ij}!},$ (2.19)
$\displaystyle\sum_{\boldsymbol{x}\in\mathcal{X}}W_{n}({\boldsymbol{x}};N,{\boldsymbol{\eta}})P_{\boldsymbol{m}}(\boldsymbol{x};{\boldsymbol{u}})P_{\boldsymbol{m}^{\prime}}(\boldsymbol{x};{\boldsymbol{u}})=0,\quad{\boldsymbol{m}}\neq{\boldsymbol{m}}^{\prime}.$
(2.20)

in which the $n^{2}$ unspecified parameters $u_{i\,j}$’s are to be determined
by the condition that they are the left eigenvectors of the transition matrix
$K_{n}({\boldsymbol{x}},{\boldsymbol{y}};N,{\boldsymbol{\alpha}},{\boldsymbol{\beta}})$
(2.11). Their explicit expressions and the so-called orthogonality conditions
will be given in Theorem 3.2 (3.10), (3.11). In the formula, $(a)_{n}$ is the
shifted factorial defined for $a\in\mathbb{C}$ and nonnegative integer $n$,
$(a)_{0}=1$, $(a)_{n}=\prod_{k=0}^{n-1}(a+k)$, $n\geq 1$ and
$M_{n}({\mathbb{N}}_{0})$ is the set of square matrices of degree $n$ with
nonnegative integer elements. The generating function of the above
hypergeometric function (2.19) is well-known [34],

$\displaystyle G(\boldsymbol{x};{\boldsymbol{u}},{\boldsymbol{t}})$
$\displaystyle\stackrel{{\scriptstyle\text{def}}}{{=}}\prod_{i=0}^{n}\left(\sum_{j=0}^{n}b_{ij}t_{j}\right)^{x_{i}}=\sum_{\boldsymbol{m}\in\mathcal{X}}\binom{N}{\boldsymbol{m}}P_{\boldsymbol{m}}(\boldsymbol{x};{\boldsymbol{u}})\boldsymbol{t}^{\boldsymbol{m}},$
(2.21)
$\displaystyle\boldsymbol{t}^{\boldsymbol{m}}\stackrel{{\scriptstyle\text{def}}}{{=}}\prod_{j=1}^{n}t_{j}^{m_{j}},\quad
t_{0}\stackrel{{\scriptstyle\text{def}}}{{=}}1.$ (2.22)

The parameters $\\{b_{i\,j}\\}$ are related to $\\{u_{i\,j}\\}$

$\displaystyle b_{0\,j}=b_{i\,0}=1\ {\rm for}\ 0\leq i,j\leq n\ {\rm and}\
b_{i\,j}=1-u_{i\,j}\ {\rm for}\ i,j=1,\ldots,n.$ (2.23)

By requiring that $P_{\boldsymbol{m}}(\boldsymbol{x};{\boldsymbol{u}})$ should
obey difference equations of Birth and Death type [28, 38], another type of
multivariate Krawtchouk polynomials are obtained [39].

##  3 Rahman polynomials

In this section I show the explicit forms of
$P_{\boldsymbol{m}}(\boldsymbol{x};{\boldsymbol{u}})$ (2.19), that is, the
expression of ${\boldsymbol{u}}=\\{u_{i\,j}\\}$ in Theorem 3.2. Based on that,
the general eigenvalues $\mathcal{E}({\boldsymbol{m}})$ are determined in
Proposition 3.7. Finally, the left eigenvalue equations

$\sum_{{\boldsymbol{x}}\in\mathcal{X}}K_{n}({\boldsymbol{x}},{\boldsymbol{y}};N,{\boldsymbol{\alpha}},{\boldsymbol{\beta}})P_{\boldsymbol{m}}(\boldsymbol{x};{\boldsymbol{u}})=\mathcal{E}({\boldsymbol{m}})P_{\boldsymbol{m}}(\boldsymbol{y};{\boldsymbol{u}}),\quad\forall{\boldsymbol{m}}\in\mathcal{X},$

are proven in Theorem 3.8.

To begin with, the structure of the hypergeometric function
$P_{\boldsymbol{m}}(\boldsymbol{x};{\boldsymbol{u}})$ (2.19) means the
following

###### Proposition 3.1

The knowledge of $n$ degree one polynomials
$P_{|m|=1}(\boldsymbol{x};{\boldsymbol{u}})$ is enough to determine
${\boldsymbol{u}}=\\{u_{i\,j}\\}$ completely.

###  3.1 Degree one left eigenpolynomials

In order to solve the degree one left eigenvalue problem of
$K_{n}({\boldsymbol{x}},{\boldsymbol{y}};N,{\boldsymbol{\alpha}},{\boldsymbol{\beta}})$
(2.11), the following summation formula is useful,

$\sum_{{\boldsymbol{x}}\in\mathcal{X}}x_{i}K_{n}({\boldsymbol{x}},{\boldsymbol{y}};N,{\boldsymbol{\alpha}},{\boldsymbol{\beta}})=\beta_{i}N-\beta_{i}\sum_{j=1}^{n}\alpha_{j}y_{j}+\alpha_{i}y_{i},\quad
i=1,\ldots,n.$ (3.1)

These are obtained by combining the recurrence relations of $W_{1}$ and
$W_{n}$,

$\displaystyle x\,W_{1}(x;t,\alpha)$ $\displaystyle=\alpha
t\,W_{1}(x-1;t-1,\alpha),\quad s,t\in\mathbb{N},$ $\displaystyle
x_{i}\,W_{n}({\boldsymbol{x}};s,{\boldsymbol{\beta}})$
$\displaystyle=\beta_{i}s\,W_{n}({\boldsymbol{x}}-{\boldsymbol{e}}_{i};s-1,{\boldsymbol{\beta}}),\quad
i=1,\ldots,n,$ (3.2)

in which ${\boldsymbol{e}}_{i}$ is the unit vector in $i$-direction. Similar
sum formulas were used in [37] to solve many generalisations of the original
single variable problem. It is interesting to note that the above formulas are
easily generalised

$\displaystyle(-x)_{k}\,W_{1}(x;t,\alpha)$
$\displaystyle=\alpha^{k}(-t)_{k}\,W_{1}(x-k;t-k,\alpha),$
$\displaystyle(-x_{i})_{k}\,W_{n}({\boldsymbol{x}};s,{\boldsymbol{\beta}})$
$\displaystyle=\beta_{i}^{k}(-s)_{k}\,W_{n}({\boldsymbol{x}}-k{\boldsymbol{e}}_{i};s-k,{\boldsymbol{\beta}}),\quad
k=1,\ldots.$

But I do not use them in this paper.

By plugging in the following form

$P_{|m|=1}({\boldsymbol{x}})=1+\frac{1}{N}\sum_{i=1}^{n}a_{i}x_{i},$ (3.3)

into the left eigenvalue equation

$\sum_{{\boldsymbol{x}}\in\mathcal{X}}K_{n}({\boldsymbol{x}},{\boldsymbol{y}};N,{\boldsymbol{\alpha}},{\boldsymbol{\beta}})P_{|m|=1}({\boldsymbol{x}})=\lambda
P_{|m|=1}({\boldsymbol{y}}),$

and using the sum formula (3.1), one obtains

$\displaystyle
1+\frac{1}{N}\sum_{i=1}^{n}a_{i}\left\\{\beta_{i}N-\beta_{i}\sum_{j=1}^{n}\alpha_{j}y_{j}+\alpha_{i}y_{i}\right\\}=\lambda\left(1+\frac{1}{N}\sum_{i=1}a_{i}y_{i}\right).$

These give rise to

$\displaystyle 1+\sum_{j=1}^{n}\beta_{j}a_{j}$ $\displaystyle=\lambda,$ (3.4)
$\displaystyle-\left(\sum_{j=1}^{n}\beta_{j}a_{j}\right)\alpha_{i}+\alpha_{i}a_{i}$
$\displaystyle=\lambda a_{i},\quad i=1,\ldots,n.$ (3.5)

For each root $\lambda_{j}$ of the characteristic equation

$Det\left(\lambda\,I_{n}-F({\boldsymbol{\alpha}},{\boldsymbol{\beta}})\right)=0,\quad
F({\boldsymbol{\alpha}},{\boldsymbol{\beta}})_{i\,j}\stackrel{{\scriptstyle\text{def}}}{{=}}-\alpha_{i}\beta_{j}+\alpha_{i}\delta_{i\,j},\quad
i,j=1,\ldots,n,$ (3.6)

the coefficients $\\{a_{i}\\}$ are determined

$a_{i\,j}=-\frac{\alpha_{i}(\lambda_{j}-1)}{\lambda_{j}-\alpha_{i}},\quad
i,j=1,\ldots,n.$

The consistency of the above formula and (3.4) implies the following sum rule

$\sum_{i=1}^{n}\frac{\alpha_{i}\beta_{i}}{\lambda_{j}-\alpha_{i}}=-1\
\Longleftarrow\sum_{i=1}^{n}\beta_{i}u_{i\,j}=1-\lambda_{j}\quad
j=1,\ldots,n,$ (3.7)

as $\lambda_{j}<1$ (2.7). Let us tentatively identify the above $j$-th
solution as $\boldsymbol{m}=\boldsymbol{e}_{j}$ left eigenpolynomial

$P_{\boldsymbol{e}_{j}}(\boldsymbol{x};{\boldsymbol{u}})=1-\frac{1}{N}\sum_{i=1}^{n}\frac{\alpha_{i}(\lambda_{j}-1)}{\lambda_{j}-\alpha_{i}}x_{i},\quad
j=1,\ldots,n.$ (3.8)

Comparison of these polynomials with the above formula (2.19) for
$P_{\boldsymbol{m}}({\boldsymbol{x}};{\boldsymbol{u}})$ leads to the following

###### Theorem 3.2

System parameters By the knowledge of the degree one left eigenvectors,

$\displaystyle\sum_{{\boldsymbol{x}}\in\mathcal{X}}K_{n}({\boldsymbol{x}},{\boldsymbol{y}};N,{\boldsymbol{\alpha}},{\boldsymbol{\beta}})P_{\boldsymbol{e}_{j}}(\boldsymbol{x};{\boldsymbol{u}})=\lambda_{j}P_{{\boldsymbol{e}}_{j}}(\boldsymbol{y};{\boldsymbol{u}}),\quad
j=1,\ldots,n,$ (3.9)

the system parameters $\\{u_{i\,j}\\}$ are completely identified,

$\displaystyle
P_{\boldsymbol{e}_{j}}(\boldsymbol{x};{\boldsymbol{u}})=1-\frac{1}{N}\sum_{i=1}^{n}u_{i\,j}x_{i}\quad\Longrightarrow
u_{i\,j}=\frac{\alpha_{i}(\lambda_{j}-1)}{\lambda_{j}-\alpha_{i}},\quad
i,j=1,\ldots,n,$ (3.10)

and each eigenvalue $\lambda_{j}$ satisfies the sum rule (3.7). As the degree
one left eigenvectors constitute the eigenvectors of the real symmetric matrix
(2.9) they are orthogonal with each other and with
$1=P_{\boldsymbol{0}}(\boldsymbol{x};{\boldsymbol{u}})$;

$\displaystyle\sum_{{\boldsymbol{x}}\in\mathcal{X}}W_{n}({\boldsymbol{x}};N,{\boldsymbol{\eta}})P_{\boldsymbol{e}_{j}}(\boldsymbol{x};{\boldsymbol{u}})=0,\
\sum_{{\boldsymbol{x}}\in\mathcal{X}}W_{n}({\boldsymbol{x}};N,{\boldsymbol{\eta}})P_{\boldsymbol{e}_{j}}(\boldsymbol{x};{\boldsymbol{u}})P_{\boldsymbol{e}_{k}}(\boldsymbol{x};{\boldsymbol{u}})=0,\
j\neq k=1,\ldots,n.$ (3.11)

It should be stressed that the system parameters $\\{u_{i\,j}\\}$ are
irrational functions of the original system parameters $\\{\alpha_{i}\\}$ and
$\\{\beta_{i}\\}$.

From the characteristic equation (3.6) one knows that the product and sum of
the eigenvalues are rational functions of the original system parameters
$\\{\alpha_{i}\\}$ and $\\{\beta_{i}\\}$.

###### Remark 3.3

$\prod_{i=1}^{n}\lambda_{i}=\prod_{i=1}^{n}\alpha_{i}\cdot(1-|\beta|),\quad\sum_{i=1}^{n}\lambda_{i}=\sum_{i=1}^{n}\alpha_{i}(1-\beta_{i}).$
(3.12)

###### Proposition 3.4

$\mathfrak{S}_{n}$ Symmetry The Rahman polynomials
$\\{P_{\boldsymbol{m}}(\boldsymbol{x})\\}$ (2.19) are invariant under the
symmetric group $\mathfrak{S}_{n}$ acting on ${\boldsymbol{m}}$, due to the
arbitrariness of the ordering of $n$ roots $\\{\lambda_{j}\\}$ of the
characteristic equation (3.6) in the parameters $u_{i\,j}$ (3.10).

####  3.1.1 Exceptional cases

For generic parameters, the eigenvalues are non-degenerate and the vanishing
of a denominator of $u_{i\,j}$, $\lambda_{j}-\alpha_{i}$ does not happen. But
as seen in the multivariate Krawtchouk and Meixner cases [39, 40], certain
exceptional parameter configurations do exist. For example, when
$\alpha_{i}=\alpha_{j}$, the characteristic matrix (3.6) with
$\lambda=\alpha_{i}$, i.e.
$\alpha_{i}I_{n}-F({\boldsymbol{\alpha}},{\boldsymbol{\beta}})$ has equal
$i$-th and $j$-th rows, $-\alpha_{i}\beta_{k}$, $k=1,\ldots,n$. Thus the
determinant vanishes and $\alpha_{i}$ becomes one of the eigenvalues, causing
some $u_{i\,k}$ singular.

###### Remark 3.5

The parameters $\\{\alpha_{i}\\}$ should be distinct so that the general form
of $P_{\boldsymbol{m}}({\boldsymbol{x}};{\boldsymbol{u}})$ (2.19) to hold.

###### Remark 3.6

The degenerates case $\alpha_{1}=\alpha_{2}$ of the bivariate Rahman
polynomials were discussed in some detain in [22].

###  3.2 Determination of the general eigenvalues

Now the whole of the Rahman polynomials are identified (2.19) and (3.10) and
their orthogonality (2.20) is guaranteed by the above (3.11) [14, 34, 39]. The
next task is to verify the left eigenvalue equations

$\sum_{{\boldsymbol{x}}\in\mathcal{X}}K_{n}({\boldsymbol{x}},{\boldsymbol{y}};N,{\boldsymbol{\alpha}},{\boldsymbol{\beta}})P_{\boldsymbol{m}}(\boldsymbol{x};{\boldsymbol{u}})=\mathcal{E}({\boldsymbol{m}})P_{\boldsymbol{m}}(\boldsymbol{y};{\boldsymbol{u}}),\quad\forall{\boldsymbol{m}}\in\mathcal{X}.$
(3.13)

For this, the knowledge of the general eigenvalues
$\mathcal{E}(\boldsymbol{m})$ is necessary. Since
$P_{\boldsymbol{m}}({\boldsymbol{0}};{\boldsymbol{u}})=1$, it is obtained by
setting ${\boldsymbol{y}}={\boldsymbol{0}}$,

$\displaystyle\mathcal{E}({\boldsymbol{m}})$
$\displaystyle=\sum_{{\boldsymbol{x}}\in\mathcal{X}}K_{n}({\boldsymbol{x}},{\boldsymbol{0}};{\boldsymbol{\alpha}},{\boldsymbol{\beta}})P_{\boldsymbol{m}}({\boldsymbol{x}};{\boldsymbol{u}})=\sum_{{\boldsymbol{x}}\in\mathcal{X}}\binom{N}{{\boldsymbol{x}}}P_{\boldsymbol{m}}({\boldsymbol{x}};{\boldsymbol{u}})\prod_{i=1}^{n}\beta_{i}^{x_{i}}\cdot(1-|\beta|)^{N-|x|}$
$\displaystyle=(1-|\beta|)^{N}\sum_{{\boldsymbol{x}}\in\mathcal{X}}\binom{N}{{\boldsymbol{x}}}P_{\boldsymbol{m}}({\boldsymbol{x}};{\boldsymbol{u}})\prod_{i=1}^{n}\left(\frac{\beta_{i}}{1-|\beta|}\right)^{x_{i}}$
$\displaystyle=(1-|\beta|)^{N}\sum_{{\boldsymbol{x}}\in\mathcal{X}}\binom{N}{{\boldsymbol{x}}}P_{\boldsymbol{x}}({\boldsymbol{m}};{\boldsymbol{u}}^{T})\prod_{i=1}^{n}\left(\frac{\beta_{i}}{1-|\beta|}\right)^{x_{i}}.$
(3.14)

Here, the duality

$P_{\boldsymbol{m}}({\boldsymbol{x}};{\boldsymbol{u}})=P_{\boldsymbol{x}}({\boldsymbol{m}};{\boldsymbol{u}}^{T})$

of the $(n+1,2n+2)$ hypergeometric function (2.19) is used. That is,
exchanging ${\boldsymbol{x}}\leftrightarrow{\boldsymbol{m}}$ is realised by
$u_{i\,j}\leftrightarrow u_{j\,i}$. The summation in r.h.s. of (3.14) has the
same structure as that of the generating function
$G({\boldsymbol{x}};{\boldsymbol{u}},{\boldsymbol{t}})$ in (2.21) except for
the exchanges ${\boldsymbol{x}}\leftrightarrow{\boldsymbol{m}}$ and
$u_{i\,j}\to u_{j\,i}$, ${\boldsymbol{t}}\to{\boldsymbol{\beta}}/(1-|\beta|)$.
This means

$\displaystyle\mathcal{E}({\boldsymbol{m}})$
$\displaystyle=(1-|\beta|)^{N}G\bigl{(}{\boldsymbol{m}},{\boldsymbol{u}}^{T},{\boldsymbol{\beta}}/(1-|\beta|)\bigr{)}$
$\displaystyle=(1-|\beta|)^{N}\prod_{i=1}^{n}\left(1+\frac{|\beta|}{1-|\beta|}\right)^{N-|m|}\left(1+\sum_{j=1}^{n}b_{i\,j}^{T}\beta_{j}/(1-|\beta|)\right)^{m_{i}}$
$\displaystyle=\prod_{i=1}^{n}\left(1-|\beta|+\sum_{j=1}^{n}b_{i\,j}^{T}\beta_{j}\right)^{m_{i}}=\prod_{i=1}^{n}\left(1-\sum_{j=1}^{n}u_{j\,i}\beta_{j}\right)^{m_{i}}=\prod_{i=1}^{n}\lambda_{i}^{m_{i}}.$

The explicit expression of $u_{i\,j}$ (3.10) and the sum rule (3.7)

$1-\sum_{j=1}^{n}u_{j\,i}\beta_{j}=1-\sum_{j=1}^{n}\frac{\alpha_{j}\beta_{j}(\lambda_{i}-1)}{\lambda_{i}-\alpha_{j}}=1-(\lambda_{i}-1)=\lambda_{i},$

leads to

###### Proposition 3.7

If $P_{\boldsymbol{m}}({\boldsymbol{x}};{\boldsymbol{u}})$ is the left
eigenvector of
$K_{n}({\boldsymbol{x}},{\boldsymbol{y}};N,{\boldsymbol{\alpha}},{\boldsymbol{\beta}})$
(3.13), it has a multiplicative spectrum

$\mathcal{E}({\boldsymbol{m}})=\prod_{i=1}^{n}\lambda_{i}^{m_{i}}.$ (3.15)

###  3.3 Left eigenpvalue equation for
$P_{\boldsymbol{m}}({\boldsymbol{x}};{\boldsymbol{u}})$

The next and the main task is to prove the left eigenvalue equations for all
the Rahman polynomials

$\sum_{{\boldsymbol{x}}\in\mathcal{X}}K_{n}({\boldsymbol{x}},{\boldsymbol{y}};N,{\boldsymbol{\alpha}},{\boldsymbol{\beta}})P_{\boldsymbol{m}}(\boldsymbol{x};{\boldsymbol{u}})=P_{\boldsymbol{m}}(\boldsymbol{y};{\boldsymbol{u}})\cdot\prod_{i=1}^{n}\lambda_{i}^{m_{i}},\quad\forall{\boldsymbol{m}}\in\mathcal{X}.$
(3.16)

Let us multiply both sides by
$\binom{N}{\boldsymbol{m}}\prod_{i=1}^{n}t_{i}^{m_{i}}$ and take
$\boldsymbol{m}$ summation to obtain

$\displaystyle\sum_{{\boldsymbol{x}}\in\mathcal{X}}K_{n}({\boldsymbol{x}},{\boldsymbol{y}};N,{\boldsymbol{\alpha}},{\boldsymbol{\beta}})\sum_{{\boldsymbol{m}}\in\mathcal{X}}\binom{N}{\boldsymbol{m}}P_{\boldsymbol{m}}(\boldsymbol{x};{\boldsymbol{u}})\prod_{i=1}^{n}t_{i}^{m_{i}}$
$\displaystyle=\sum_{{\boldsymbol{m}}\in\mathcal{X}}\binom{N}{\boldsymbol{m}}P_{\boldsymbol{m}}(\boldsymbol{y};{\boldsymbol{u}})\cdot\prod_{i=1}^{n}(\lambda_{i}t_{i})^{m_{i}}.$

By comparing with the expansion formula (2.21) of the generating function

$G({\boldsymbol{x}};{\boldsymbol{u}},{\boldsymbol{t}})=(1+|t|)^{N-|x|}\prod_{i=1}^{n}T_{i}({\boldsymbol{t}})^{x_{i}},\quad
T_{i}({\boldsymbol{t}})\stackrel{{\scriptstyle\text{def}}}{{=}}1+\sum_{j=1}^{n}b_{i\,j}t_{j},$
(3.17)

the above equation reads

$\displaystyle\sum_{{\boldsymbol{x}}\in\mathcal{X}}K_{n}({\boldsymbol{x}},{\boldsymbol{y}};N,{\boldsymbol{\alpha}},{\boldsymbol{\beta}})G({\boldsymbol{x}};{\boldsymbol{u}},{\boldsymbol{t}})=G({\boldsymbol{y}};{\boldsymbol{u}},{\boldsymbol{\lambda}}{\boldsymbol{t}}).$
(3.18)

The main goal is achieved by proving this equation, which is rather tedious
but straightforward.

$\displaystyle\sum_{{\boldsymbol{x}}\in\mathcal{X}}K_{n}({\boldsymbol{x}},{\boldsymbol{y}};N,{\boldsymbol{\alpha}},{\boldsymbol{\beta}})G({\boldsymbol{x}};{\boldsymbol{u}},{\boldsymbol{t}})$
$\displaystyle=\sum_{{\boldsymbol{x}},{\boldsymbol{z}}\in\mathcal{X}}\binom{N-|z|}{{\boldsymbol{x}}-{\boldsymbol{z}}}\prod_{i=1}^{n}\beta_{i}^{x_{i}-z_{i}}\left((1-|\beta|)(1+|t|)\right)^{N-|x|}T_{i}({\boldsymbol{t}})^{x_{i}}\cdot\prod_{i=1}\binom{y_{i}}{z_{i}}\alpha_{i}^{z_{i}}(1-\alpha_{i})^{y_{i}-z_{i}}$
$\displaystyle=\sum_{{\boldsymbol{x}},{\boldsymbol{z}}\in\mathcal{X}}\binom{N-|z|}{{\boldsymbol{x}}-{\boldsymbol{z}}}\\!\prod_{i=1}^{n}\left(\beta_{i}T_{i}({\boldsymbol{t}})\right)^{x_{i}-z_{i}}\\!\left((1-|\beta|)(1+|t|)\right)^{N-|x|}\cdot\prod_{i=1}\binom{y_{i}}{z_{i}}\\!\left(\alpha_{i}T_{i}({\boldsymbol{t}})\right)^{z_{i}}(1-\alpha_{i})^{y_{i}-z_{i}}$
$\displaystyle=\sum_{{\boldsymbol{z}}\in\mathcal{X}}\Bigl{(}(1-|\beta|)(1+|t|)+\sum_{i=1}^{n}\beta_{i}T_{i}({\boldsymbol{t}})\Bigr{)}^{N-|z|}\cdot\prod_{i=1}\binom{y_{i}}{z_{i}}\left(\alpha_{i}T_{i}({\boldsymbol{t}})\right)^{z_{i}}(1-\alpha_{i})^{y_{i}-z_{i}}$
$\displaystyle=(1+|\lambda
t|)^{N}\prod_{i=1}^{n}\sum_{z_{i}=1}^{y_{i}}\binom{y_{i}}{z_{i}}\left(\frac{\alpha_{i}T_{i}({\boldsymbol{t}})}{1+|\lambda
t|}\right)^{z_{i}}(1-\alpha_{i})^{y_{i}-z_{i}}.$ (3.19)

Here the simplification of the first term is used.

$\displaystyle\sum_{i=1}^{n}\beta_{i}T_{i}({\boldsymbol{t}})$
$\displaystyle=\sum_{i=1}^{n}\beta_{i}\left(1+\sum_{j=1}^{n}(1-u_{i\,j})t_{j}\right)=|\beta|(1+|t|)-\sum_{i\,j}\beta_{i}u_{i\,j}t_{j}$
$\displaystyle=|\beta|(1+|t|)-|t|+|\lambda t|,$
$\displaystyle\qquad\Longrightarrow(1-|\beta|)(1+|t|)+\sum_{i=1}^{n}\beta_{i}T_{i}({\boldsymbol{t}})=1+|\lambda
t|.$

At the last equality, the sum rule (3.7) is used. Now, the $z_{i}$ summation
gives

$\displaystyle\sum_{z_{i}=1}^{y_{i}}\binom{y_{i}}{z_{i}}\left(\frac{\alpha_{i}T_{i}({\boldsymbol{t}})}{1+|\lambda
t|}\right)^{z_{i}}(1-\alpha_{i})^{y_{i}-z_{i}}$
$\displaystyle=\left(1-\alpha_{i}+\frac{\alpha_{i}T_{i}({\boldsymbol{t}})}{1+|\lambda
t|}\right)^{y_{i}}$ $\displaystyle=(1+|\lambda
t|)^{-y_{i}}\bigl{(}(1-\alpha_{i})(1+|\lambda
t|)+\alpha_{i}T_{i}({\boldsymbol{t}})\bigr{)}.$

A straightforward calculation based on the explicit form of $u_{i\,j}$ (3.10)
shows

$(1-\alpha_{i})(1+|\lambda
t|)+\alpha_{i}T_{i}({\boldsymbol{t}})=T_{i}({\boldsymbol{\lambda}}{\boldsymbol{t}}).$
(3.20)

Finally the desired result is obtained

$\displaystyle\sum_{{\boldsymbol{x}}\in\mathcal{X}}K_{n}({\boldsymbol{x}},{\boldsymbol{y}};N,{\boldsymbol{\alpha}},{\boldsymbol{\beta}})G({\boldsymbol{x}};{\boldsymbol{u}},{\boldsymbol{t}})=(1+|\lambda
t|)^{N-|y|}\prod_{i=1}^{n}T_{i}({\boldsymbol{\lambda}}{\boldsymbol{t}})^{y_{i}}=G({\boldsymbol{y}};{\boldsymbol{u}},{\boldsymbol{\lambda}}{\boldsymbol{t}}).$
(3.21)

###### Theorem 3.8

The Rahman polynomials
$\\{P_{\boldsymbol{m}}({\boldsymbol{x}};{\boldsymbol{u}})\\}$ (2.19) are the
left eigenvectors of the Markov chain operator
$K_{n}({\boldsymbol{x}},{\boldsymbol{y}};N,{\boldsymbol{\alpha}},{\boldsymbol{\beta}})$
with the multiplicative eigenvalues $\prod_{i=1}^{n}\lambda_{i}^{mi}$ (3.15),
in which $\lambda_{i}$’s are the eigenvalues of the degree one polynomials,
i.e. the roots of the $\\\ \times n$ characteristic equation (3.6),

$\sum_{{\boldsymbol{x}}\in\mathcal{X}}K_{n}({\boldsymbol{x}},{\boldsymbol{y}};N,{\boldsymbol{\alpha}},{\boldsymbol{\beta}})P_{\boldsymbol{m}}(\boldsymbol{x};{\boldsymbol{u}})=P_{\boldsymbol{m}}(\boldsymbol{y};{\boldsymbol{u}})\cdot\prod_{i=1}^{n}\lambda_{i}^{m_{i}},\quad\forall{\boldsymbol{m}}\in\mathcal{X}.$
(3.16)

Likewise
$\\{P_{\boldsymbol{m}}({\boldsymbol{x}};{\boldsymbol{u}})W_{n}({\boldsymbol{x}};N,{\boldsymbol{\eta}})\\}$
are the eigenvectors of
$K_{n}({\boldsymbol{x}},{\boldsymbol{y}};N,{\boldsymbol{\alpha}},{\boldsymbol{\beta}})$
with the same eigenvalue $\prod_{i=1}^{n}\lambda_{i}^{mi}$ (3.15),

$\sum_{{\boldsymbol{y}}\in\mathcal{X}}K_{n}({\boldsymbol{x}},{\boldsymbol{y}};N,{\boldsymbol{\alpha}},{\boldsymbol{\beta}})P_{\boldsymbol{m}}(\boldsymbol{y};{\boldsymbol{u}})W_{n}({\boldsymbol{y}};N,{\boldsymbol{\eta}})=\prod_{i=1}^{n}\lambda_{i}^{m_{i}}\cdot
P_{\boldsymbol{m}}(\boldsymbol{x};{\boldsymbol{u}})W_{n}({\boldsymbol{x}};N,{\boldsymbol{\eta}}),\quad\forall{\boldsymbol{m}}\in\mathcal{X},$
(3.22)

in which $W_{n}({\boldsymbol{x}};N,{\boldsymbol{\eta}})$ is the reversible
distribution for
$K_{n}({\boldsymbol{x}},{\boldsymbol{y}};N,{\boldsymbol{\alpha}},{\boldsymbol{\beta}})$
(2.14).

##  4 Comments

The most important point of this paper is that the Rahman polynomials are
obtained as the left eigenvectors of the Markov chain operator
$K_{n}({\boldsymbol{x}},{\boldsymbol{y}};N,{\boldsymbol{\alpha}},{\boldsymbol{\beta}})$.
Therefore, they are completely determined by the initial system parameters
$\\{{\boldsymbol{\alpha}}\\}$ and $\\{{\boldsymbol{\beta}}\\}$ alone as shown
in Theorems 3.2 3.8. In contrast, in some literature, the eigenvalue problem
was addressed [21](3.17), [22](1.5), [18]§6, [20](1.7). This naturally
involves the reversible distribution
$W_{n}({\boldsymbol{x}};N,{\boldsymbol{\eta}})$

$\sum_{{\boldsymbol{y}}\in\mathcal{X}}K_{n}({\boldsymbol{x}},{\boldsymbol{y}};N,{\boldsymbol{\alpha}},{\boldsymbol{\beta}})P_{\boldsymbol{m}}(\boldsymbol{y};{\boldsymbol{u}})W_{n}({\boldsymbol{y}};N,{\boldsymbol{\eta}})=\prod_{i=1}^{n}\lambda_{i}^{m_{i}}\cdot
P_{\boldsymbol{m}}(\boldsymbol{x};{\boldsymbol{u}})W_{n}({\boldsymbol{x}};N,{\boldsymbol{\eta}}),\quad\forall{\boldsymbol{m}}\in\mathcal{X},$

and the parameters ${\boldsymbol{\eta}}$ come into the play. Efforts were made
to construct Rahman polynomials from the fact that they had the orthogonality
measure $W_{n}({\boldsymbol{x}};N,{\boldsymbol{\eta}})$ and
${\boldsymbol{\eta}}$’s came into the polynomial parameters. Of course some
combinations of the original parameters ${\boldsymbol{\alpha}}$ and
${\boldsymbol{\beta}}$ can be expressed by ${\boldsymbol{\eta}}$’s based on
the relation I.(1.6). Therefore, it might be possible to mingle some
${\boldsymbol{\eta}}$’s into the expression of the system parameters
$u_{i\,j}$. But it is totally unnecessary.

To say symbolically, “The parameters ${\boldsymbol{\eta}}$’s are chosen by
$K_{n}({\boldsymbol{x}},{\boldsymbol{y}};N,{\boldsymbol{\alpha}},{\boldsymbol{\beta}})$
to provide the orthogonality weight for the Rahman polynomials. The
orthogonality weights never determine multivariate orthogonal polynomials, as
shown explicitly for the cases of the Krawtchouk and Meixner polynomials [39,
40].”

In some works, the importance of the left eigenvalue problems was well
recognised [3](2.2), [21](3.23), [20]§6, but unfortunately the idea was not
well developed. A superficial reason might be that they adopted notation like
$K_{n}({\boldsymbol{j}},{\boldsymbol{i}};N,{\boldsymbol{\alpha}},{\boldsymbol{\beta}})$
[20](1.7) which was rather distant from polynomials in ${\boldsymbol{x}}$,
$P_{\boldsymbol{m}}({\boldsymbol{x}};{\boldsymbol{u}})$ as compared with mine,
$K_{n}({\boldsymbol{x}},{\boldsymbol{y}};N,{\boldsymbol{\alpha}},{\boldsymbol{\beta}})$.
A similar situation existed in the treatments of birth and death processes.
They used $\lambda_{n}$ and $\mu_{n}$ for the birth/death rates [29] instead
of $B(x)$ and $D(x)$ [38] and $B_{j}({\boldsymbol{x}})$ and
$D_{j}({\boldsymbol{x}})$ for the multivariate cases [39, 40].

It should be emphasised that the eigenvalues were correctly identified [22,
20] in terms of $u_{i\,j}$, as the eigenvalues are common in the eigenvalue
problem and the left eigenvalue problem. For the bivariate case [22],

$\lambda_{m,n}=(1-\beta_{1}t-\beta_{2}u)^{m}(1-\beta_{1}v-\beta_{2}w)^{n},$
([22].4.23)

Which is the same as the sum rule (3.7) and
$\lambda_{m,n}=\lambda_{1}^{m}\lambda_{2}^{n}$, as
$\left(\begin{array}[]{cc}t&v\\\
u&w\end{array}\right)=\left(\begin{array}[]{cc}u_{1\,1}&u_{1\,2}\\\
u_{2\,1}&u_{2\,2}\end{array}\right)$. For the general multivariate case
[20]§7, $\lambda_{\boldsymbol{m}}=\prod_{i=1}^{n}(1-\omega_{i})^{m_{i}}$, in
which $\omega_{i}=\sum_{j=1}\beta_{j}u_{i\,j}$. This is again the same as my
multiplicative spectrum (3.15) by using the sum rule (3.7). Their $u_{i\,j}$
is my $u_{j\,i}$. The explicit expressions of $u_{i\,j}$ in terms of the
original system parameters ${\boldsymbol{\alpha}}$ and ${\boldsymbol{\beta}}$
were not available in their works.

## References

  * [1] Aomoto K., Kita M., Hypergeometric funtions, Springer, Berlin, 1994 (in Japanese).
  * [2] Andrews G.E., Askey R. and Roy R., Special Functions, Encyclopedia of mathematics and its applications, Cambridge Univ. Press, Cambridge, (1999).
  * [3] Cooper R. D., Hoare M. R. and Rahman M., “Stochastic Processes and Special Functions: On the Probabilistic Origin of Some Positive Kernels Associated with Classical Orthogonal Polynomials,” J. Math. Anal. Appl. 61 (1977) 262-291.
  * [4] Diaconis P. and Griffiths R. C., “An introduction to multivariate Krawtchouk polynomials and their applications,” J. Stat. Planning and Inference 154 (2014) 39-53, arXiv:1309.0112[math.PR].
  * [5] Dunkl C. F. “A difference equation and Hahn polynomials in two variables,” Pacific J. Math. 92 (1981) 57–71.
  * [6] Dunkl C. F. and Xu Y., Orthogonal polynomials of several variables, Encyclopedia of mathematics and its applications, Cambridge Univ. Press, Cambridge, (2001).
  * [7] Feller W., An Introduction to Probability Theory and its Applications, I, (2nd ed.), Wiley, New York, (1966).
  * [8] Gelfand I. M., “General theory of hypergeometric functions,” Sov. Math. Dokl. 33 (1986) 573-577.
  * [9] Genest V. X., Vinet L. and Zhedanov A., “The multivariate Krawtchouk polynomials as matrix elements of the rotation group representations on oscillator states,” J. Phys. A46 (2013) 0505203, arXiv:1306.4256[math-ph].
  * [10] Genest V. X., Miki H., and Vinet L., “The multivariate Meixner polynomials as matrix elements of $SO(d,1)$ representations on oscillator states,” J. Phys. A47 (2014) 045207, arXiv:1310.6953[math-ph].
  * [11] Genest V. X., Miki H., Vinet L. and Zhedanov A., “Spin lattices, state transfer and bivariate Krawtchouk polynomials,” Canadian J. Phys. 93 (2014) 979–984, arXiv:1410.4703[math-ph].
  * [12] Genest V. X., Miki H., Vinet L. and Zhedanov A., “The multivariate Charlier polynomials as matrix elements of the Euclidean group representation on oscillator states,” J. Phys. A47 (2014) 215204 (16pp), arXiv:1401.7901[math-ph].
  * [13] Genest V. X., Post S., and Vinet L., “An algebraic interpretation of the multivariate q-Krawtchouk polynomials,” Ramanujan J. 43 (2017) 415–445. arXiv:1508.07770[math.CA].
  * [14] Griffiths R. C., “Orthogonal polynomials on the multinomial distribution,” Austral. J. Statist. 13 (1971) 27–35.
  * [15] Griffiths R. C., “Orthogonal Polynomials on the Negative Multinomial Distribution,” J. Multivariate Anal. 5 (1975) 271–277.
  * [16] Griffiths R., “Lancaster distributions and Markov chains with multivariate Poisson-Charlier, Meixner and Hermite-Chebycheff polynomial eigenfunctions,” J. Appr. Theor. 207 (2016) 139–164, arXiv:1412.3931v3[math.PR].
  * [17] Griffiths R. C., “Multivariate Krawtchouk Polynomials and Composition Birth and Death Processes,” Symmetry 8 (2016) 33 19pp, arXiv:1603.00196[math.PR].
  * [18] Grünbaum F. A., “The Rahman polynomials are bispectral,” SIGMA 3 (2007) 065, 11pp, arXiv:0705.0468[math.CA].
  * [19] Grünbaum F. A. and Rahman M., “On a family of 2-variable orthogonal Krawtchouk polynomials,” SIGMA 6 (2010) 090, 12 pages, arXiv:1007.4327[math.CA].
  * [20] Grünbaum F. A. and Rahman M., “A System of multivariable Krawtchouk polynomials and a probabilistic application,” SIGMA 7 (2011) 119, 17pp, arXiv:1106.1835[math.PR].
  * [21] Hoare M. R. and Rahman M., “ Cumulative Bernoulli trials and Krawtchouk processes,” Stochastic Process. Appl. 16 (1983) 113–139.
  * [22] Hoare M. R. and Rahman M., “A probabilistic origin for a new class of bivariate polynomials,” SIGMA 4 (2008) 089, 18 pages, arXiv:0812.3879[math.CA].
  * [23] Iliev P. and Xu Y., “Discrete orthogonal polynomials and difference equations of several variables,” Adv. Math. 212 (2007) 1–36, arXiv:math.CA/0508039.
  * [24] Iliev P. and Terwilliger P., “The Rahman polynomials and the Lie algebra $sl_{3}(C)$,” Trans. Amer. Math. Soc. 364 (2012) 4225–4238, arXiv:1006.5062[math.RT].
  * [25] Iliev P., “A Lie theoretic interpretation of multivariate hypergeometric polynomials,” Compositio Math. 148 (2012) 991-1002, arXiv:1101.1683[math.RT].
  * [26] Iliev P., “Meixner polynomials in several variables satisfying bispectral difference equations,” Adv. in Appl. Math. 49 (2012) 15–23, arXiv:1112.5589[math.CA]
  * [27] Iliev P., “Gaudin model for the multinomial distribution,” arXiv:2303.08206[math-ph].
  * [28] Ismail M. E. H., Classical and quantum orthogonal polynomials in one variable, Encyclopedia of mathematics and its applications, Cambridge, (2005).
  * [29] Karlin S. and McGregor J. L., “The differential equations of birth-and-death processes,” Trans. Amer. Math. Soc. 85 (1957) 489–546; “Linear growth, birth-and-death processes,” J. Math. Mech. 7 (1958) 643–662; “Ehrenfest urn models,” J. Appl. Prob. 19 (1965) 477–487.
  * [30] Koekoek R., Lesky P. A. and Swarttouw R. F., Hypergeometric orthogonal polynomials and their $q$-analogues, Springer Monographs in Mathematics, Springer-Verlag, Berlin, (2010).
  * [31] Miki H., Tsujimoto S., Vinet L. and Zhedanov A., “An Algebraic Model for the Multiple Meixner Polynomials of the First Kind,” J. Phys. A45 (2012) 325205, arXiv:1203.0357[math-ph].
  * [32] Miki H., Tsujimoto S., Vinet L. and Zhedanov A., “Perfect state transfer in two dimensions and the bivariate dual-Hahn polynomials,” Prog. Theor. Exper. Phys. 2022 (2022) 053A01, arXiv:2009.06195[math-ph].
  * [33] Mizukawa H., “Zonal spherical functions on the complex reflection groups and $(m+1,n+1)$-hypergeometric functions,” Adv. Math. 184 (2004) 1–17.
  * [34] Mizukawa H., “Orthogonal relations for multivariate Krawtchouk polynomials,” SIGMA 7 (2011) 017 5pp, arXiv:1009.1203[math.CO].
  * [35] Mizukawa H. and Tanaka H., $(n+1,m+1)$-hypergeometric functions associated to character algebras, Proc. Amer. Math. Soc. 132 (2004) 2613–2618.
  * [36] Odake S. and Sasaki R., “Orthogonal Polynomials from Hermitian Matrices,” J. Math. Phys. 49 (2008) 053503 (43 pp), arXiv:0712.4106[math.CA].
  * [37] Odake S. and Sasaki R., “Markov chains generated by convolutions of orthogonality measures,” J. Phys. A: Math. Gen. 55 (2022) 275201 (42pp), arXiv:2106.04082[math.PR]
  * [38] Sasaki R., “Exactly Solvable Birth and Death Processes,” J. Math. Phys. 50 (2009) 103509 (18 pp), arXiv:0903.3097[math-ph].
  * [39] Sasaki R., “Multivariate Kawtchouk polynomials as Birth and Death polynomials,” arXiv:2305.08581v2[math.CA].
  * [40] Sasaki R., “Multivariate Meixner polynomials as Birth and Death polynomials,” arXiv:2310.04968[math.CA].
  * [41] Shibukawa, G., “Multivariate Meixner, Charlier and Krawtchouk Polynomials According to Analysis on Symmetric Cones,” J. Lie Theory 26 (2016) 439–477, arXiv:1404.7491[math.CA].
  * [42] Tratnik M.V., “Some multivariable orthogonal polynomials of the Askey tableau-discrete families,” J. Math. Phys. 32 (1991), 2337–2342.
  * [43] Xu Y., “Hahn, Jacobi, and Krawtchouk polynomials of several variables,” Journal of Approximation Theory 195 (2015) 19–42, arXiv:1309.1510[math.CA].
  * [44] Zhedanov A., “9j-symbols of the oscillator algebra and Krawtchouk polynomials in two variables,” J. Phys. A: Math. Gen. 30 (1997) 8337-8353.