A method for computing the Perron root for primitive matrices

TL;DR

This paper introduces an efficient iterative method to compute the Perron root of primitive matrices by constructing a matrix with equal row and column sums, reducing computational load compared to traditional methods.

Contribution

The paper presents a novel iterative algorithm that computes the Perron root without explicit knowledge of the auxiliary matrix, improving efficiency over existing methods.

Findings

Algorithm converges similarly to the power method

Reduces computational load in Perron root calculation

Provides a new approach for first eigenvector computation

Abstract

Following the Perron-Frobenius theorem, the spectral radius of a primitive matrix is a simple eigenvalue. It is shown that for a primitive matrix $A$, there is a positive rank one matrix $X$ such that $B = A \circ X$, where $\circ$ denotes the Hadamard product of matrices, and such that the row (column) sums of matrix $B$ are the same and equal to the Perron root. An iterative algorithm is presented to obtain matrix $B$ without an explicit knowledge of $X$. The convergence rate of this algorithm is similar to that of the power method but it uses less computational load. A byproduct of the proposed algorithm is a new method for calculating the first eigenvector.

Equations

$$A\circ B=(a_{ij}{b}_{ij})$$

$$\alpha A\circ \beta B=\alpha \beta (A\circ B)$$

$$A\circ B=(\mathbf{u}{\mathbf{v}}^T)\circ (\mathbf{x}{\mathbf{y}}^T)=(\mathbf{u}\circ \mathbf{x})(\mathbf{v}\circ \mathbf{y}{)}^T.$$

$$\underset{i=1,\dots ,n}{\min }\left\{r_i(A)\right\}\le \rho (A)\le \underset{i=1,\dots ,n}{\max }\left\{r_i(A)\right\}$$

$$\underset{j=1,\dots ,n}{\min }\left\{c_j(A)\right\}\le \rho (A)\le \underset{j=1,\dots ,n}{\max }\left\{c_j(A)\right\}$$

$$A=\left(\begin{array}[]{ccc}0&1&0\\ 3&0&3\\ 0&2&0\end{array}\right)$$

$$B=D_{\mathbf{x}}^{-1}A{D}_{\mathbf{x}}$$

$$I_n+B$$

$$(I_n+B{)}^{n-1}$$

$$\underset{i=1,\dots ,n}{\min }\left\{r_i(D_{\mathbf{r}}^{-1}A{D}_{\mathbf{r}})\right\}\le \rho (A)\le \underset{i=1,\dots ,n}{\max }\left\{r_i(D_{\mathbf{r}}^{-1}A{D}_{\mathbf{r}})\right\}$$

$$\underset{j=1,\dots ,n}{\min }\left\{c_j(D_{\mathbf{c}}^{-1}A{D}_{\mathbf{c}})\right\}\le \rho (A)\le \underset{j=1,\dots ,n}{\max }\left\{c_j(D_{\mathbf{c}}^{-1}A{D}_{\mathbf{c}})\right\}$$

$$D_{\mathbf{r}}^{-1}A{D}_{\mathbf{r}}$$

$$=$$

$$x_{ij}=\frac{r_j(A)}{r_i(A)}\text{mbox};i,j=1,2\dots ,n$$

$$A\circ X=\left(\begin{array}[]{cc}a_{11}&a_{12}x\\ a_{21}/x&a_{22}\end{array}\right)$$

$$x=\frac{S-a_{11}}{a_{12}}\text{mbox};\frac{1}{x}=\frac{S-a_{22}}{a_{21}}$$

$$S^2-(a_{11}-a_{22})S+a_{11}{a}_{22}-a_{12}{a}_{21}=0$$

$$S=\left(a_{11}+a_{22}+\sqrt{(a_{11}-a_{22}{)}^2+4a_{12}{a}_{21}}\right)/2$$

$$\sum _{j=1}^n{a}_{ij}{x}_j=\lambda x_i$$

$$r_i=\frac{1}{x_i}\left(\sum _{j=1}^n{a}_{ij}{x}_j\right)$$

$$X=\left(\begin{array}[]{ccccc}1&x_{2}/x_{1}&x_{3}/x_{1}&\ldots&x_{n}/x_{1}\\ x_{1}/x_{2}&1&x_{3}/x_{2}&\ldots&x_{n}/x_{2}\\ x_{1}/x_{3}&x_{2}/x_{3}&1&\ldots&x_{n}/x_{3}\\ \vdots&\vdots&\vdots&\ddots&\vdots\\ x_{1}/x_{n}&x_{2}/x_{n}&x_{3}/x_{n}&\ldots&1\end{array}\right)$$

$$A^{(t)}=D_{{\mathbf{r}}^{(t-1)}}^{-1}{A}^{(t-1)}{D}_{{\mathbf{r}}^{(t-1)}}\text{mbox};t=1,2,\dots$$

$$a_{ij}^{(t)}$$

$$r_i^{(t)}(A^{(t)})$$

$$A^{(1)}$$

$$A^{(2)}$$

$$A^{(t)}$$

$$B=A^{(t)}=A\circ X$$

$$X$$

$$x_{ij}$$

$$X=\mathbf{x}{\mathbf{y}}^T$$

$$\underset{i}{\min }\left\{r_i^{(t)}(A^{(t)})\right\}$$

$$\underset{i}{\max }\left\{r_i^{(t)}(A^{(t)})\right\}$$

$${\xi }^{(t)}$$

Full text

Keywords: primitive matrix; Perron root; Markov chain; stochastic matrix.

MSC(2010): 15A18 15A48 15A03 65F10 65F15 65C40

1 Introduction

Given a nonnegative matrix, the problem of computing the first eigenvalue and eigenvector is considered in this paper. For two matrices $A=(a_{ij})$ and $B=(b_{ij})$ with the same number of rows and columns, their Hadamard product is a matrix of elementwise products:

[TABLE]

For scalars $\alpha$ and $\beta$:

[TABLE]

If $A$ and $B$ are rank one matrices, i.e. $A=\mathbf{u}\mathbf{v}^{T}$ and $B=\mathbf{x}\mathbf{y}^{T}$ then

[TABLE]

Many properties for Hadamard product are given in [1],[2, chapter 5].

If $A=(a_{ij})\in\mathbb{R}^{n\times n}$, then $A$ is called positive if $a_{ij}>0$, and nonnegative if $a_{ij}\geq 0$. Perron [3] showed that the spectral radius of a positive matrix $A$ is a simple eigenvalue [4, page 667] that dominates all other eigenvalues in modulus. This eigenvalue, denoted $\rho(A)$, is called the Perron root and the associated normalized positive vector is called the Perron vector. Nonnegative matrices are frequently encountered in real life applications [5, 6]. Frobenius [7] extented Perron’s work on positive matrices to nonnegative matrices. The spectral radius of a nonnegative matrix $A$ is positive if it is irreducible, i.e. $(I_{n}+A)^{n-1}$ is a positive matrix [8, page 534], [4, page 672]. The dominant eigenvalue of an irreducible matrix is unique if it is primitive [8, page 540], [4, page 674]. A nonnegative matrix is primitive if $A^{m}$ is a positive matrix for some non nul $m$ [8, page 540], [4, page 678]. Wielandt, [9], showed that a nonnegative matrix $A$ of order $n$ is primitive if $A^{n^{2}-2n+2}$ is a positive matrix [8, page 543]. To verify primitivity of a nonnegative matrix using Frobenius or Wielandt formula leads to huge calculations especially when $n$ is high. It is shown in [8, page 544] that only some power calculations of the matrix are necessary.

This paper is on the calculation of the Perron root. The power method is generally used to obtain the eigenvalue with the maximum modulus and associated eigenvector [8, page 545], [10, page 330], [4, page 533]. The convergence rate of the power method depends on the ratio of the second eigenvalue to the first [10, page 330], [4, page 533]. More iterations will be required when the modulus of the second highest eigenvalue is close to that of the first. Here, an iterative algorithm is proposed for calculating the Perron root for primitive matrices. This algorithm is based on successive improvement of bounds for the Perron root. There are many research works on localization of the Perron root for nonnegative matrices [11, 12, 13, 14, 15]. Frobenius carried out the following bounds [8, page 521]:

[TABLE]

where $r_{i}(A)=\sum_{j=1}^{n}a_{ij}$ and $c_{j}(A)=\sum_{i=1}^{n}a_{ij}$ are the row and column sums of $A$, respectively. In (4) and (5), equalities occur when $\rho(A)$ is equal to the row or column sums. The column sums of the matrix in (6) are both equal to $3$.

[TABLE]

The matrix $A$ in (6) is imprimitive and its eigenvalues are: $3$, $-3$ and [math]. This example shows that equality in (4) or (5) can occur for an imprimitive matrix.

Let $\textbf{x}=(x_{1},x_{2},\ldots,x_{n})$ a vector with only positive values, $x_{i}>0$, and $D_{\mathbf{x}}$ a diagonal matrix formed with x. The matrix $B$ defined by:

[TABLE]

is diagonally similar to $A$ [16], and we have:

Lemma 1.1**.**

The matrices $A$ and $B$ in (7) have the same eigenvalues, and

• a)

if $A$ is irreducible, then $B$ is irreducible,

• b)

if $A$ is primitive, then $B$ is primitive,

Proof.

a) if $A$ is irreducible then $(I_{n}+A)^{n-1}$ is a positive matrix. From (7), we have:

[TABLE]

$D_{\mathbf{x}}$ has only positive values and $(I_{n}+A)^{n-1}$ is a positive matrix. The matrix $(I_{n}+B)^{n-1}$ is then positive that implies irreducibility of the matrix $B$.

b) if $A$ is a primitive matrix then there exists an integer $m$ such that $A^{m}$ is positive. Using (7) we have $B^{m}=D_{\mathbf{x}}^{-1}A^{m}D_{\mathbf{x}}$, that implies $B^{m}$ is positive and the result follows. ∎

Using an improvement of bounds in (4) and (5) by Minc [17], relation (7) and the uniqueness of eigenvalue with a maximum modulus for a primitite matrix, an iterative algorithm is proposed to obtaining the Perron root.

2 Methods

Lemma 2.1**.**

Let $A=(a_{ij})\in\mathcal{R}^{n\times n}$ a nonnegative matrix. If matrix $A$ has a row (column) with only zero entries, then $A$ cannot be a primitive matrix.

Proof.

see an exercise on irreducible matrices in [8, page 522]. ∎

For a primitive matrix, from the Lemma 2.1 and relations (4)-(5), we have the following two observations. The minimum value of the row (column) sums for a primitive matrix is greater than zero. The maximum value of the row (column) sums for a primitive matrix is greater than or equal to the Perron root.

Let us note $D_{\mathbf{r}}$ and $D_{\mathbf{c}}$ diagonal matrices formed with the row sums $\mathbf{r}=(r_{1}(A),r_{2}(A),\ldots,r_{n}(A))$ and the column sums $\mathbf{c}=(c_{1}(A),c_{2}(A),\ldots,c_{n}(A))$ of A, respectively. The Frobenius bounds (4) and (5) have been improved by Minc [17, page 27]:

[TABLE]

In (10) and (11), equalities hold when the row or column sums are the same and correspond to the Perron root. Using the row sums relation, (10) allows to write:

[TABLE]

where $X$ is a positive matrix formed with:

[TABLE]

The unicity of the Perron root for a primitive matrix and (28) suggest that the components of the matrix $X$ can be chosen to have the same row sums for $A\circ X$. For a second order nonnegative matrix ($n=2$), we have:

[TABLE]

where $x=r_{2}(A)/r_{1}(A)$, $r_{1}(A)=a_{11}+a_{12}$ and $r_{2}(A)=a_{21}+a_{22}.$

If the row sums in (30) are the same and equal to $S$, an expression can be obtained for the parameter $x$:

[TABLE]

To have a value for $x$, $a_{12}$ and $a_{21}$ should be nonzero. Relation (31) allows to have a second order equation which resolution leads to a value for $S$:

[TABLE]

The solution of (32) with the maximum modulus is:

[TABLE]

A second order nonnegative matrix $A$ is primitive if $a_{12}$, and $a_{21}$ are both nonzero, on the one hand. On the other hand, at least $a_{11}$ or $a_{22}$ should be nonzero. Hence, for a second order primitive matrix, explicit expressions can be obtained for matrix $X$ (parameter $x$) and the Perron root, (33). However, a direct search for components of the matrix $X$ in (28) becomes difficult when $n>2$.

Lemma 2.2**.**

Let $A=(a_{ij})$ a square matrix of order $n$, $\mathbf{y}=\alpha\mathbf{x}$ a vector where $\alpha$ is a nonzero scalar and $\mathbf{x}$ is the eigenvector associated with eigenvalue $\lambda$ of $A$. If all components of $\mathbf{x}$ have nonzero value, the row sums of the matrix $D_{\mathbf{y}}^{-1}AD_{\mathbf{y}}$ are the same and equal to the eigenvalue $\lambda$ of $A$.

A similar result is obtained using $A^{T}$ or the column sums.

Proof.

Using the definition of the eigenvalue, the component $i$ of $A\mathbf{x}=\lambda\mathbf{x}$ is:

[TABLE]

The component $i$ of the row sums of the matrix $D_{\mathbf{y}}^{-1}AD_{\mathbf{y}}$ is:

[TABLE]

By replacing the right hand expression of (34) in (35) the result follows. ∎

Compared to the Frobenius bounds (4) and (5), the bounds in (10) and 11 are based on a modification of the initial matrix. This process can be repeated to further sharpen the bounds. From the unicity of the Perron root for a primitive matrix, a repeative improvement of the Minc bounds will lead to equalities of the row (column) sums. The main result of this paper is the following.

Theorem 2.1**.**

Let $A$ be a primitive matrix of order $n$. There exists a positive rank one matrix $X$ of the form

[TABLE]

such that the matrix $B=A\circ X$ is similar to $A$. In addition, the row (column) sums of $B$ are the same and equal to the Perron root of $A$.

Proof.

The row sums are used for the proof, the column sums can also be used in a similar way.

Let us note $A^{(0)}$ the initial matrix and its row sums vector as $\mathbf{r}^{(0)}$. Relation (10) allows to write:

[TABLE]

From (28) and (29), the components of matrix $A^{(t)}$ and its row sums are:

[TABLE]

At iteration $t$, $A^{(t)}$ is similar to $A^{(t-1)}$. From the Minc relation (10), the bounds with $A^{(t)}$ are improved compared to those with $A^{(t-1)}$. Hence, when $t\to\infty$, equalies hold in (10) for primitive matrix and the row sums $r_{i}^{(t)}(A^{(t)}),i=1,2,\ldots,n$, have the same value which is equal to $\rho(A)$, thank to Lemma 2.2.

From (37) and the notation in (28), we can write:

[TABLE]

At the convergence iteration $t$, $X^{(t)}$ is formed with only $1$. Then, from (42) we have

[TABLE]

where:

[TABLE]

Relation (43) is another form of (17), then matrix $B$ is similar to matrix $A$. Since the row sums of matrix $B$ are the same, equalities occur in (4) and lead to the Perron root.

The matrix $X$ in (36) can be write as a product of two vectors:

[TABLE]

where $\mathbf{x}=(1,x_{1}/x_{2},\ldots,x_{1}/x_{n})^{T}$ and $\mathbf{y}=(1,x_{2}/x_{1},\ldots,x_{n}/x_{1})^{T}$, i.e. the first column and first row of the matrix $X$, respectively. $X$ is then a rank one matrix, and is also positive because formed with the row sums of a primitive matrix, see Lemma 2.1. ∎

Corollary 2.1**.**

Vector $\mathbf{y}$ allowing to obtain the matrix $X$ in (46) is in the space spanned by the Perron vector of a primitive matrix $A$.

Proof.

Use (17), (28) and Lemma 2.2. ∎

2.1 Convergence of the proposed algorithm

Theorem 2.2**.**

The iterative algorithm based on a successive improvement of the Minc bounds is convergent for a primitive matrix.

Proof.

At iteration $t$, the row sum vectors associated with matrices $A^{(t)}$ and $A^{(t-1)}$ are $\mathbf{r}^{(t)}(A^{(t)})$ and $\mathbf{r}^{(t-1)}(A^{(t-1)})$, respectively. From the Minc theorem [17, page 27], we have:

[TABLE]

Let us define two decreasing sequences as follows:

[TABLE]

From (48) and (50), we have

[TABLE]

where $0<c^{(t)}\leq 1$. Let $\zeta^{(0)}$ denotes the initial value obtained using (50). From relation (52) we have:

[TABLE]

Since $c^{(i)}$, $i=1,2,\ldots,t$, are positive numbers not all equal to $1$, we have

[TABLE]

A similar reasoning using (47) and (49) leads to $\lim_{t\to\infty}\xi^{(t)}=0$. Hence, when the number of iteration goes to infinity, relations (49) and (50) show that the row sums obtained with the algorithm converges to the Perron root. Referring to Lemma 2.2, it is like a vector with only ones is used to obtaining to Perron root. ∎

Relation (53) can be used to estimate the minimum number of iterations required by the algorithm before convergence when an error level $\alpha$ is set. Assuming that $E(c^{(i)})=c$ is the mean of the $c^{(i)}$ coefficients, relation (53) becomes $\zeta^{(t)}=c^{t}\zeta^{(0)}$ and we have

[TABLE]

Corollary 2.2**.**

For primitive matrices, the convergence rate of the proposed algorithm is similar to that of the power method and depends on the magnitude of the second highest eigenvalue.

Proof.

The matrix $A$ can be write as the sum of rank one matrices using its eigenvalues and eigenvectors:

[TABLE]

From (57), (46) and (3) we have:

[TABLE]

Using (56) and (57), an expression for power $k$ of matrix $A$ is

[TABLE]

This expression allows to have another one similar to (58). Then, the row sums at the first step of the algorithm using $A^{k}$ are:

[TABLE]

where $\delta_{i}=(\mathbf{v}_{i}\circ\mathbf{y}^{(0)})^{T}\mathbf{1}$, $\mathbf{z}_{i}=(\mathbf{u}_{i}\circ\mathbf{x}^{(0)})$, $i=1,2,\ldots,n$,

$\mathbf{x}^{(0)}=\left(1,\frac{r_{1}^{(0)}}{r_{2}^{(0)}},\frac{r_{1}^{(0)}}{r_{3}^{(0)}},\ldots,\frac{r_{1}^{(0)}}{r_{n}^{(0)}}\right)^{T}$ and $\mathbf{y}^{(0)}=\left(1,\frac{r_{2}^{(0)}}{r_{1}^{(0)}},\frac{r_{3}^{(0)}}{r_{1}^{(0)}},\ldots,\frac{r_{n}^{(0)}}{r_{1}^{(0)}}\right)^{T}$

From corollary 2.1, vector $(A\circ X)\mathbf{1}\in span\left\{(A^{k}\circ X^{(0)})\mathbf{1}\right\}$ then,

[TABLE]

∎

2.2 Algorithm and implementation

Relations (38) and (39) allow to obtain an algorithm for computing the matrix $B$ in Theorem 2.1. A convergence test is based on the difference between the maximum and the minimum values of the row or column sums, i.e. the range value.

[TABLE]

Another convergence test can be based on the examination of the minimum and maximum row sum values, i.e. by using the decreasing sequences $\xi^{(t)}$ and $\zeta^{(t)}$ defined in Theorem 2.2.

2.2.1 Algorithm A (using row sums): Perron root only

Initialization

• •

set: $t\leftarrow 0$, $a_{ij}^{(t)}\leftarrow a_{ij}$, calculate the row sums using (39)

• •

set stopping rules: $eps$ (the acceptable error), $maxIter$ (the maximum number of iterations), compute the initial error value using (63) 2. 2.

while ($error>eps$ and $t<maxIter$)

• •

update matrix: (38)

• •

calculate row sums: (39)

• •

compute error: (63)

• •

increase iteration number: $t\leftarrow t+1$

Remark 2.1**.**

As mentioned, Theorem 2.1 is also valid using column sums. In the implementation, one can compute the row and column sums and perform the next step using the sum where the initial error is the lowest.

Remark 2.2**.**

From an iteration to the next, the error should decrease by an amount that depends on the convergence rate. Otherwise, we must stop the algorithm because the matrix does not seem to be primitive. This observation can be used as an indirect test for primitivity of a matrix.

Indeed, if the spectral radius is not simple (case of an irreducible imprimitive matrix) there may be at least two eigenvectors associated with eigenvalues having the same modulus.

Remark 2.3**.**

With the proposed algorithm, the diagonal entries of $A$ remain unchanged, only the off-diagonal components are modified. The Gerschgorin discs [8, page 388] associated with a matrix allow to illustrate this.

Let us consider an example:

[TABLE]

The off-diagonal components of $A$ are modified in such a way all discs cross the same highest point on the x-axis ($4$ for this example). To show this, let us write:

[TABLE]

where $D_{A}$ is a diagonal matrix formed with the diagonal elements of $A$ and $P$ is matrix $A$ where the diagonal elements are set to zero. From (43), we have:

[TABLE]

Hence, the row sums of matrix $B$ are given by:

[TABLE]

where $\mathbf{p}_{i.}$ and $\mathbf{x}_{i.}$ are vectors formed with row $i$ of matrices $P$ and $X$, respectively.

Remark 2.4**.**

Compared to the power method, there is no initial vector to set. The results obtained using the power method are the Perron root and the associated eigenvalue. Only the Perron root is obtained using this algorithm. However, since the row (column) sums of the matrix $B$ in (43) are the same, a vector $\mathbf{1}$ formed with only ones is an eigenvector of $B$ ($B^{T}$).

[TABLE]

Remark 2.5**.**

The proposed algorithm uses another matrix in comparison with the power method. At each iteration, the total numbers of multiplications and additions of the matrix-vector multiplication by the power method are equal to the total number of operations for the proposed algorithm. Hence, using the power method, the additional arithmetic operations used for calculating the eigenvector and the eigenvalue are extra computational load compared to the proposed algorithm.

2.2.2 Algorithm B (using row sums): Perron root and vector

Instead of the algorithm A (2.2.1), another one can consist in searching for a vector $\mathbf{y}$ similar to the Perron vector. For this purpose, the matrix $B$ is initially equal to $A$ and the vector $\mathbf{y}$ is set to $\mathbf{1}$, then, the row sums of $B$ are calculated. At iteration $t$, a vector $\mathbf{y}^{(t)}$ is formed using row sums and the matrix $B$ is updated. At the convergence, we should have $\mathbf{y}^{(t)}\approx\mathbf{1}$.

Initialization

• •

set: $t\leftarrow 0$, $B\leftarrow A$, $\mathbf{y}\leftarrow\mathbf{1}$ and compute row sums $r_{i}^{(0)}$ of $B$,

• •

calculate initial error: $\max(r_{i}^{(0)})-\min(r_{i}^{(0)})$.

• •

set stopping rules: $eps$ and $maxIter$ (see Algorithm A) 2. 2.

while ($error>eps$ and $t<maxIter$)

• •

form $\mathbf{y}^{(t)}$ using row sums $r_{i}^{(t)}$

• •

update $\mathbf{y}$: $\mathbf{y}\circ\mathbf{y}^{(t)}$

• •

form $\mathbf{x}$ ($1/\mathbf{y}$) and update matrix $B$: $A\circ(\mathbf{x}*\mathbf{y}^{T})$

• •

compute error: $\max|\mathbf{y}^{(t)}-\mathbf{1}|$

• •

increase iteration number: $t\leftarrow t+1$

• •

calculate row sums $r_{i}^{(t)}$ of $B$

The convergence test of this algorithm consists to have only ones for vector $\mathbf{y}$. Instead, the convergence test can be based on the examination of the minimum and maximum row sum values of the decreasing sequences $\xi^{(t)}$ and $\zeta^{(t)}$ defined in Theorem 2.2. The code in appendix A.1 is the R implementations of algorithm B using this test. An input matrix should be square, nonnegative and each row sum should be greater than zero.

2.2.3 Application to row-stochastic matrices

Let us consider Markov chains modeling a dynamic system with finite discrete $n$ states. At each time $t$, this system is in one state. When the system is in state $i$, the move to state $j$ or to stay in state $i$ is controlled by the probability $p_{ij}\geq 0$, $\sum_{j=1}^{n}p_{ij}=1$. The probabilities are organized in a transition matrix $P$ which is nonnegative. Interestingly, the power $k$ of matrix $P$ is also a transition matrix which element $p_{ij}$ corresponds to the probability to move from state $i$ to state $j$ at time $t+k$. The Markov chains are used in: (a) biological sequences analysis [18, 19], (b) internet traffic or search engines [20, 21, 22], (c) …The state occupied by the Markov chain process at time $t+k$ depends on the nature of transition matrix $P$: reducible/irreducible, imprimitive/primitive. In some applications, the observed transition matrix is modified to be primitive [20]. For this particular case, the power of the modified matrix $\hat{P}$ converges to a matrix formed with the same vector $\mathbf{u}$ verifying: $\mathbf{u}^{T}\hat{P}=\mathbf{u}^{T}$. The vector $\mathbf{u}$ is the left-hand vector of the matrix $\hat{P}$. In the search engines based on the Markov chains, the components of the vector $\mathbf{u}$ allow to hierarchize the information. Applying algorithm B (2.2.2) to the transposed of the row-stochastic (modified transition) matrix lead to a vector $\mathbf{y}$ which components are then normalized in a way that their sum is $1$. This normalized vector corresponds to $\mathbf{u}$.

3 Results and conclusions

All calculations were performed on the same computer (a laptop equipped with i7-66600U processor, 16 GB of RAM, under Microsoft Windows 10) and R version 3.6.2. The default error level was arbitrarily set to 1.0E-8.

The first example is:

[TABLE]

The range values for the row and the column sums are $4$ and $2.5$ respectively. The algorithm is performed using the column sums.

Figure 1 presents the Gerschgorin discs for all iterations. The Perron root for this example is $5.739952$, the proposed algorithm and the power method require $17$ and $19$ iterations, respectively. Figure 1 shows that the major modifications of the off-diagonal elements of the matrix A are done during the first five iterations.

For all of the tests performed, the algorithm proposed and the power method have close number of iterations. Worse results, in term of the number of iterations, were obtained using a tridiagonal matrix. Let $T(n;c,a,b)$ a tridiagonal matrix of order $n$, where $a$ is the value for the diagonal components, $b$ is the value for the upper diagonal components and $c$ is the value for the under diagonal components. An explicit expression relating eigenvalues of matrix $T$ is available [23]:

[TABLE]

The ratio $\lambda_{2}/\lambda_{1}$ for $T$ is near $1$ when $n$ is high. For $n=50$, $a=3$, $b=2$ and $c=1$ the first two eigenvalues of matrix $T$ are: $5.823063$ and $5.806989$. The proposed algorithm took $5,890$ (algorithm A) or $5,174$ (algorithm B) iterations to calculate the Perron root. The power method need $5,159$ iterations.

Except the cases where the modulus of the second eigenvalue is near to that of the first, the algorithm proposed converges after few iterations, especially when the first eigenvalue is largely dominant. The proposed method has been succesfully used for a matrix of order $15,515$ that results from high-throughput biological data.

With a convergence rate similar to that of the classic power method, the proposed algorithm for computing the Perron root is computationally less demanding. But, it applies to only primitive matrices. However, it can be used as a low cost primitivity test compared to the matrix power calculations involved in the Frobenius and Wielandt tests.

Acknowledgements

This work was supported by funds from CNRS, INSERM and University of Strasbourg.

Author is grateful to a referee for the valuable comments and suggestions.

Appendix A R code using row sums

A.1 Algorithm B

This function computes iteratively the Perron root

and the eigenvector of matrix A using row sums

A: nonnegative square matrix (all row sums are > 0)

tol: error level used (stopping criterion)

maxIter: maximum number of iterations (stooping criterion)

Returned

B: matrix which has the same row sums (B = A o X)

pfr: minimum and maximum row sum, that defines to the Perron root

y: vector (leading to have the eigenvector and matrix X)

iter: number of iterations performed

rmin: sequences with minimum row sum values for iterations

rmax: sequences with maximum row sum values for iterations

calcPRc <- function(A, tol=1.0e-8, maxIter=50) { n <- nrow(A); m <- ncol(A) ri <- apply(A, 1, sum) ko <- (sum(A<0) || (min(ri)==0)) if ((n != m) || (ko)) { stop("calcPRc(): for nonnegative primitive matrices") } y <- rep(1,n) iter <- 1; B <- A rmin <- c(); erMin <- rmin[iter] <- min(ri) rmax <- c(); erMax <- rmax[iter] <- max(ri) erIter <- ((erMin > tol) || (erMax > tol)) while (erIter && (iter < maxIter)) { yt <- ri/ri[1]; y <- yyt B <- A * ((1/y) %% t(y)) ri <- apply(B, 1, sum) iter <- iter + 1 rmin[iter] <- ri.min <- min(ri) rmax[iter] <- ri.max <- max(ri) erMin <- rmin[iter] - rmin[iter-1] erMax <- rmax[iter-1] - rmax[iter] erIter <- ((erMin > tol) || (erMax > tol)) } pfr <- c(ri.min, ri.max) list(B=B, pfr=pfr, y=y, iter=iter-1, rmin=rmin, rmax=rmax) }