Home A matrix variate inverse Lomax distribution
Article Open Access

A matrix variate inverse Lomax distribution

  • Daya K. Nagar , Alejandro Roldán-Correa and Saralees Nadarajah EMAIL logo
Published/Copyright: September 11, 2025

Abstract

This study presents the matrix variate inverse Lomax distribution as a generalization of the univariate inverse Lomax distribution and shows that this distribution can be derived by using matrix variate gamma distributions. We study several properties such as cumulative distribution function, marginal distribution of sub-matrix, triangular factorization, moment generating function, Kullback-Liebler divergence, entropies and expected values of the inverse Lomax matrix. We also derive maximum likelihood estimators of the parameters and check their finite sample performance by simulations. Some of these results are expressed in terms of special functions of matrix arguments and zonal polynomials.

1 Introduction

The inverse Lomax distribution, also called the inverse Pareto Distribution of the second kind, with shape parameter α > 0 and location parameter σ > 0 is given by the probability density function

(1) α σ α w α 1 1 + w σ ( α + 1 ) , w > 0 , α > 0 .

For σ = 1 , the inverse Lomax distribution slides to a standard inverse Lomax distribution or inverse Pareto distribution of the second kind. We will write w IL ( σ , α ) if w follows the inverse Lomax distribution given by the density (1) and if σ = 1 , then we simply write w IL ( α ) .

By transforming v = 1 w in the probability density function of the inverse Lomax distribution, the Lomax distribution is defined by

(2) α λ 1 + v λ ( α + 1 ) , v > 0 ,

where λ = 1 σ . The Lomax distribution, named after K. S. Lomax, is a heavy-tail probability distribution often used in business, economics, and actuarial modeling. A short notation to designate that v has the two parameter Lomax distribution is w IL ( λ , α ) and the standard Lomax distribution or Pareto Distribution of the second kind with λ = 1 is denoted by w IL ( α ) .

The inverse Lomax distribution is one of the significant lifetime models, and it is used in economic sciences, geography, econometrics, and clinical fields; the inverse Lomax distribution was used by Kleiber [19] to obtain Lorenz ordering of order statistics. Singh et al. [28,29] estimated the reliability for the inverse Lomax distribution based on censored type II observations. Rahman et al. [23] discussed estimation and prediction of an inverse Lomax model using Bayesian approach under various loss functions. Bantan et al. [7] estimated the entropy for the inverse Lomax distribution under multiple censored data. Al-Omar et al. [5] estimated the reliability of the inverse Lomax distribution using extreme ranked set sampling. Jan and Ahmad [17] performed Bayesian analysis of the inverse Lomax distribution using approximation techniques. Yadav et al. [31] performed Bayesian estimation of P ( Y < X ) for the inverse Lomax distribution under progressive type-II censoring scheme. Alam et al. [2] studied step stress partially accelerated life test under adaptive type-II progressive hybrid censoring for the inverse Lomax distribution. Ahmed et al. [1] investigated a competing risks model with type-II generalized hybrid censored data from the inverse Lomax distribution. Hora et al. [15] provided classical and Bayesian inference for the inverse Lomax distribution under adaptive progressive type-II censored data with COVID-19 application. Yadav et al. [33] provided characterizations and estimation of the inverse Lomax distribution under progressively type-II censoring. Ramadan et al. [24] estimated multi stress-strength reliability based on progressive first failure with lifetime following the inverse Lomax distribution. Amer and Shalabi [6] constructed an optimal design of step stress partially accelerated life test under progressive type-II censored data with random removal for the inverse Lomax distribution. Tolba et al. [30] estimated stress-strength reliability using the inverse Lomax lifetime distribution with mechanical engineering applications. Al-Muhayfith and Al-Bossly [4] provided inference of type-I hybrid censored samples from the inverse Lomax distribution. Yadav et al. [32] considered hybrid censoring of the inverse Lomax distribution with application to survival data.

Although the literature offers a wealth of results on inverse Lomax distribution [3,9,10,14,18,21] nothing seems to have been done to define and investigate matrix variate inverse Lomax distributions.

Therefore, the aim of this study is to define matrix variate inverse Lomax distribution and examine several of its properties. The contents of this study are organized as follows. Some known definitions and results needed to define the matrix variate inverse Lomax distribution and study its properties are stated in Section 2. The definition of the matrix variate inverse Lomax distribution and its derivation are given in Section 3. Several properties such as invariance, marginal and conditional distributions of sub-matrices are derived in Section 4. Section 5 deals with the Cholesky decomposition of the inverse Lomax matrix and several independent properties. A collection of expected values of functions of inverse Lomax matrix are given in Section 6. Sections 7 and 8 deal with the Kullback-Liebler divergence and entropies. Section 9 gives the joint probability density of eigenvalues of the matrix distributed as inverse Lomax while Section 10 explores maximum likelihood estimation of parameters. Finally, Section 11 checks finite sample performance of the maximum likelihood estimators by simulations.

2 Preliminaries

We will need the following standard notations (cf. Gupta and Nagar [12]). Let A = ( a i , j ) be an m × m matrix. Then, A T denotes the transpose of A ; tr ( A ) = trace of A = a 1,1 + + a m , m ; etr ( A ) = exp ( tr ( A ) ) ; det ( A ) = determinant of A ; diag ( A ) = diagonal matrix of order m whose entries are the m leading diagonal elements of A = ( δ i , j a i , j ) , δ i , j is the Kronecker delta; A = norm of A ; A > 0 means that A is symmetric positive definite and A 1 2 denotes the unique symmetric positive definite square root of A > 0 . The sub matrices A ( α ) and A ( α ) , 1 α m , of the matrix A are defined as A ( α ) = ( a i , j ) , 1 i , j α , and A ( α ) = ( a i , j ) , α i , j m , respectively.

The well-known multivariate gamma function, denoted by Γ m ( a ) , is defined by

(3) Γ m ( a ) = X > 0 etr ( X ) det ( X ) a ( m + 1 ) 2 d X = π m ( m 1 ) 4 i = 1 m Γ a i 1 2 , Re ( a ) > m 1 2 ,

where the differential element d X stands for the exterior product of m ( m + 1 ) 2 differentials d x 1,1 , , d x 1 , m , , d x m 1 , m 1 , d x m 1 , m , d x m , m .

A number of Jacobian of transformations are available in the studies by Muirhead [20] and Gupta and Nagar [12]. A few that we will need in this study are

X = U U T d X = 2 m i = 1 m u i , i m + 1 i d U , X = V V T d X = 2 m i = 1 m v i , i i d V , X = A Y A T d X = det ( A A T ) ( m + 1 ) 2 d A , X = ( I m + Y ) 1 Y d X = det ( I m + Y ) ( m + 1 ) d Y ,

where X and Y are m × m real symmetric positive definite matrices of functionally independent variables, U = ( u i , j ) and V = ( v i , j ) are m × m lower and upper triangular matrices with positive diagonal elements, respectively, and A is a constant nonsingular matrix of order m .

The multivariate generalization of the beta function is defined by

(4) B m ( a , b ) = 0 I m det ( X ) a ( m + 1 ) 2 det ( I m X ) b ( m + 1 ) 2 d X = Γ m ( a ) Γ m ( b ) Γ m ( a + b ) = B m ( b , a ) ,

where Re ( a ) > ( m 1 ) 2 and Re ( b ) > ( m 1 ) 2 . Further, by using the matrix transformation X = ( I m + Y ) 1 Y in (4) with the Jacobian J ( X Y ) = det ( I m + Y ) ( m + 1 ) , one can easily establish the identity

(5) B m ( a , b ) = Y > 0 det ( Y ) a ( m + 1 ) 2 det ( I m + Y ) ( a + b ) d Y .

The beta type 1 and beta type 2 families of distributions are defined by the probability density functions [18]

(6) { B ( α , β ) } 1 u α 1 ( 1 u ) β 1 , 0 < u < 1

and

(7) { B ( α , β ) } 1 v α 1 ( 1 + v ) ( α + β ) , v > 0 ,

respectively, where α > 0 , β > 0 , and B ( α , β ) is the usual beta function.

If a random variable u has the probability density function (6), then we will write u B1 ( α , β ) , and if the probability density function of the random variable v is given by (7), then v B2 ( α , β ) . Further, it is well known that v ( 1 + v ) B1 ( α , β ) , 1 ( 1 + v ) B1 ( β , α ) , u ( 1 u ) B2 ( α , β ) , and ( 1 u ) u B1 ( β , α ) .

Next we define matrix variate generalizations of (6) and (7) (Gupta and Nagar [1113]).

Definition 2.1

An m × m random symmetric positive definite matrix U is said to have a matrix variate beta type 1 distribution with parameters ( α , β ) , denoted as U B1 ( m , α , β ) , if its probability density function is given by

(8) det ( U ) α ( m + 1 ) 2 det ( I m U ) β ( m + 1 ) 2 B m ( α , β ) , 0 < U < I m ,

where α > ( m 1 ) 2 and β > ( m 1 ) 2 .

Definition 2.2

An m × m random symmetric positive definite matrix V is said to have a matrix variate beta type 2 distribution with parameters ( α , β ) , denoted as V B2 ( m , α , β ) , if its probability density function is given by

(9) det ( V ) α ( m + 1 ) 2 det ( I m + V ) ( α + β ) B m ( α , β ) , V > 0 ,

where α > ( m 1 ) 2 and β > ( m 1 ) 2 .

Using definitions and suitable transformations, we can establish that V ( I m + V ) 1 B1 ( m , α , β ) , ( I m + V ) 1 B1 ( m , β , α ) , U ( I m U ) 1 B2 ( m , α , β ) , and ( I m U ) U 1 B1 ( m , β , α ) .

3 Matrix variate inverse Lomax distribution

First, we define different forms (regular and inverse) of matrix variate Lomax distribution or matrix variate Pareto distribution of the second kind.

Definition 3.1

An m × m random symmetric positive definite matrix X is said to have a matrix variate inverse Lomax distribution or matrix variate Pareto distribution of the second kind with parameters Σ and α , denoted as X IL ( m , Σ , α ) , if its probability density function is given by

(10) f MIL ( X ; Σ , α ) = Γ m [ α + ( m + 1 ) 2 ] Γ m ( α ) Γ m [ ( m + 1 ) 2 ] det ( Σ ) α det ( X ) α ( m + 1 ) 2 det ( I m + Σ 1 X ) α ( m + 1 ) 2 , X > 0 ,

where Σ is an m × m symmetric positive definite matrix and α > ( m 1 ) 2 .

Definition 3.2

An m × m random symmetric positive definite matrix U is said to have a matrix variate Lomax distribution with parameters Λ and β , denoted as U L ( m , Λ , β ) , if its probability density function is given by

(11) Γ m [ β + ( m + 1 ) 2 ] Γ m ( β ) Γ m [ ( m + 1 ) 2 ] det ( Λ ) ( m + 1 ) 2 det ( I m + Λ 1 U ) β ( m + 1 ) 2 , U > 0 ,

where Λ is an m × m symmetric positive definite matrix and β > ( m 1 ) 2 .

The standard forms of matrix variate inverse Lomax distribution and matrix variate Lomax distribution can be obtained by substituting Σ = I m and Λ = I m in Definition 3.1 and Definition 3.2, respectively, and in these cases, we will write X IL ( m , α ) and U L ( m , β ) .

From Definition 3.1, we see that if W IL ( m , α ) , then for an m × m symmetric positive definite constant matrix Σ , Σ 1 2 W Σ 1 2 IL ( m , Σ , α ) and if X IL ( m , Σ , α ) , then Σ 1 2 X Σ 1 2 IL ( m , Σ , α ) . Further, If W IL ( m , α ) , then W 1 L ( m , Σ 1 , α ) .

For m = 1 , the matrix variate inverse Lomax distribution (matrix variate inverse Pareto distribution of second kind) and the matrix variate Lomax distribution (matrix variate Pareto distribution of second kind) reduce to their respective univariate forms.

The matrix variate inverse Lomax distribution can be derived by using independent gamma matrices. A random matrix Y is said to have a matrix variate gamma distribution with parameters Ψ ( > 0 ) , and κ ( > ( m 1 ) 2 ) , denoted by Y Ga ( m , κ , Ψ ) , if its probability density function is given by

etr ( Ψ 1 Y ) det ( Y ) κ ( m + 1 ) 2 Γ m ( κ ) det ( Ψ ) κ , Y > 0 .

Theorem 3.1

Let Y 1 and Y 2 be independent, Y 1 Ga ( m , ( m + 1 ) 2 , Ψ ) and Y 2 Ga ( m , α , I m ) . Then, Y 1 1 2 Y 2 Y 1 1 2 IL ( m , Ψ 1 , α ) .

Proof

The joint probability density function of Y 1 and Y 2 is given by

etr [ ( Ψ 1 Y 1 + Y 2 ) ] det ( Y 2 ) α ( m + 1 ) 2 det ( Ψ ) ( m + 1 ) 2 Γ m [ ( m + 1 ) 2 ] Γ m ( α ) , Y 1 > 0 , Y 2 > 0 .

Transforming W = Y 1 1 2 Y 2 Y 1 1 2 and V = Y 1 with the Jacobian J ( Y 1 , Y 2 V , W ) = det ( V ) ( m + 1 ) 2 in the joint probability density of Y 1 and Y 2 , we obtain the joint probability density of W and V as

(12) etr [ V 1 2 ( Ψ 1 + W ) V 1 2 ] det ( W ) α ( m + 1 ) 2 det ( V ) α + ( m + 1 ) 2 ( m + 1 ) 2 det ( Ψ ) ( m + 1 ) 2 Γ m [ ( m + 1 ) 2 ] Γ m ( α ) ,

where W > 0 and V > 0 . Now, the desired result is obtained by integrating V in (12) by using (3).□

The cumulative distribution function of W , W IL ( m , Σ , α ) , can be obtained as

P ( W < Ω ) = Γ m [ α + ( m + 1 ) 2 ] det ( Σ ) α Γ m ( α ) Γ m [ ( m + 1 ) 2 ] 0 Ω det ( W ) α ( m + 1 ) 2 det ( I m + Σ 1 W ) α ( m + 1 ) 2 d W = Γ m [ α + ( m + 1 ) 2 ] Γ m ( α ) Γ m [ ( m + 1 ) 2 ] 0 Ω 1 det ( Z ) α ( m + 1 ) 2 d Z ,

where Ω 1 = ( I m + Σ 1 2 Ω Σ 1 2 ) 1 Σ 1 2 Ω Σ 1 2 and the last line has been obtained by substituting Y = Σ 1 2 W Σ 1 2 and Z = ( I m + Y ) 1 Y with the Jacobian J ( W Z ) = J ( W Y ) J ( Y Z ) = det ( Σ ) ( m + 1 ) 2 det ( I m Z ) ( m + 1 ) . Next, substituting G = Ω 1 1 2 Z Ω 1 1 2 with the Jacobian J ( Z G ) = det ( Ω 1 ) ( m + 1 ) 2 above, we obtain

P ( W < Ω ) = Γ m [ α + ( m + 1 ) 2 ] Γ m ( α ) Γ m [ ( m + 1 ) 2 ] det ( Ω 1 ) α 0 I m det ( G ) α ( m + 1 ) 2 d G = det ( ( Σ + Ω ) 1 Ω ) α ,

where we have used (4).

The moment generating function of X IL ( m , Σ , α ) can be derived as

M X ( Z ) = Γ m [ α + ( m + 1 ) 2 ] Γ m ( α ) Γ m [ ( m + 1 ) 2 ] det ( Σ ) α X > 0 etr ( Z X ) det ( X ) α ( m + 1 ) 2 det ( I m + Σ 1 X ) α ( m + 1 ) 2 d X ,

where Z ( m × m ) = ( ( 1 + δ i , j ) z i , j 2 ) . Now, evaluating the above integral, we obtain

M X ( Z ) = Γ m [ α + ( m + 1 ) 2 ] Γ m [ ( m + 1 ) 2 ] Ψ ( α , 0 ; Σ Z ) ,

where the confluent hypergeometric function Ψ of m × m symmetric matrix argument Y is defined by the integral [20, p. 472]

Ψ ( a , c ; Y ) = 1 Γ m ( a ) R > 0 etr ( R Y ) det ( R ) a ( m + 1 ) 2 det ( I m + R ) c a ( m + 1 ) 2 d R

valid for Re ( Y ) > 0 and Re ( a ) > ( m 1 ) 2 .

4 Properties

In this section, we give several properties of the matrix variate inverse Lomax distribution defined in Section 3.

Theorem 4.1

Let X IL ( m , Σ , α ) and A be an m × m constant nonsingular matrix. Then, the probability density function of Y = A X A T is

Γ m [ α + ( m + 1 ) 2 ] det ( A Σ A T ) α Γ m ( α ) Γ m [ ( m + 1 ) 2 ] det ( Y ) α ( m + 1 ) 2 det ( I m + ( A Σ A T ) 1 Y ) α ( m + 1 ) 2 , Y > 0 .

That is, A X A T IL ( m , A Σ A T , α ) and ( A Σ A T ) 1 2 A X A T ( A Σ A T ) 1 2 IL ( m , α ) .

Corollary 4.1.1

Let X IL ( m , Σ , α ) and A T A = Σ 1 , then A X A T IL ( m , α ) .

In the next theorem, we show that the matrix variate inverse Lomax distribution is orthogonally invariant.

Theorem 4.2

Let W IL ( m , α ) and H be an orthogonal matrix order m whose elements are either constants or random variables distributed independent of W. Then, the distribution of W is invariant under the transformation W H W H T . Further, if H is a random matrix, then W and H are distributed independently.

Proof

First, let H be a constant orthogonal matrix. Then, from Theorem 4.1, H W H T IL ( m , α ) . If, however, H is a random orthogonal matrix, then H W H T H IL ( m , α ) . Since this distribution does not depend on H we have independence between H and H W H T , and H W H T IL ( m , α ) .□

Theorem 4.3

Suppose that W IL ( m , Σ , α ) , then the mode of the distribution is given by

W = α ( m + 1 ) 2 m + 1 Σ .

Proof

Taking logarithm on both sides of the density of the matrix variate inverted Lomax distribution given by (10), we obtain

ln f MIL ( W ; Σ , α ) = constant + α m + 1 2 ln det ( W ) α + m + 1 2 ln det ( Σ + W ) .

Now, differentiating (refer equation 17.46 of [27, p. 368]) the above expression with respect to W gives

f MIL ( W ; Σ , α ) W = 2 M diag ( M ) ,

where

M = α m + 1 2 W 1 α + m + 1 2 ( Σ + W ) 1 .

Equating the derivative to 0 yields M = 0 and

α m + 1 2 W 1 = α + m + 1 2 ( Σ + W ) 1 .

The final result is obtained from the above equation after re-arrangement of terms.□

Theorem 4.4

Let W = W 1,1 W 1,2 W 2,1 W 2,2 , where W 1,1 is a q × q matrix. Define W 1,1 2 = W 1,1 W 1,2 W 2,2 1 W 2,1 and W 2,2 1 = W 2,2 W 2,1 W 1,1 1 W 1,2 . If W IL ( m , α ) , then (i) W 1,1 and W 2,2 1 are independent, W 1,1 IL ( q , α ) and W 2,2 1 B2 ( m q , α q 2 , ( m + 1 ) 2 ) , (ii) W 2,2 and W 1,1 2 are independent, W 2,2 IL ( m q , α ) and W 1,1 2 B2 ( q , α ( m q ) 2 , ( m + 1 ) 2 ) .

Proof

From the partition of W , we have

det ( W ) = det ( W 1,1 ) det ( W 2,2 1 ) ,

(13) det ( I m + W ) = det ( I q + W 1,1 ) det ( I m q + W 2,2 1 + W 2,1 W 1,1 1 ( I q + W 1,1 ) 1 W 1,2 ) .

Now, making the transformation W 1,1 = W 1,1 , R = W 2,1 W 1,1 1 2 , and W 2,2 1 = W 2,2 W 2,1 W 1,1 1 W 1,2 = W 2,2 R R T with the Jacobian J ( W 1,1 , W 2,2 , W 2,1 W 1,1 , W 2,2 1 , R ) = det ( W 1,1 ) ( m q ) 2 in the density of W , we obtain the joint density of W 1,1 , W 2,2 1 , and R as

Γ m [ α + ( m + 1 ) 2 ] Γ m ( α ) Γ m [ ( m + 1 ) 2 ] det ( W 1,1 ) α ( q + 1 ) 2 det ( I q + W 1,1 ) α ( m + 1 ) 2 det ( W 2,2 1 ) α ( m + 1 ) 2 det ( I m q + W 2,2 1 + R ( I q + W 1,1 ) 1 R T ) α ( m + 1 ) 2 .

Further, transforming Y = ( I m q + W 2,2 1 ) 1 2 R ( I q + W 1,1 ) 1 2 with the Jacobian J ( R Y ) = det ( I m q + W 2,2 1 ) q 2 det ( I q + W 1,1 ) ( m q ) 2 , the joint density of W 1,1 , W 2,2 1 , and Y is derived as

Γ m [ α + ( m + 1 ) 2 ] Γ m ( α ) Γ m [ ( m + 1 ) 2 ] det ( W 1,1 ) α ( q + 1 ) 2 det ( I m + W 1,1 ) α ( q + 1 ) 2 det ( W 2,2 1 ) α ( m + 1 ) 2 det ( I m q + W 2,2 1 ) α ( m q + 1 ) 2 det ( I m q + Y Y T ) α ( m + 1 ) 2 ,

where W 1,1 > 0 , W 2,2 1 > 0 , and Y R ( m q ) × q . From the above factorization, we see that W 1,1 and W 2,2 1 are independent, W 1,1 IL ( q , α ) and W 2,2 1 B2 ( m q , α q 2 , ( m + 1 ) 2 ) . The second part is similar.□

Theorem 4.5

Let A be a q × m constant matrix of rank q ( m ) . If X IL ( m , Σ , α ) , then

[ ( A Σ 1 A T ) 1 2 A X 1 A T ( A Σ 1 A T ) 1 2 ] 1 B2 q , α m q 2 , m + 1 2

and

( A Σ A T ) 1 2 A X A T ( A Σ A T ) 1 2 IL ( q , α ) .

Proof

Let B = A Σ 1 2 and W = Σ 1 2 X Σ 1 2 , where Σ 1 2 is a symmetric positive definite square root of Σ . Then, W IL ( m , α ) and

( A X 1 A T ) 1 = ( A Σ 1 2 Σ 1 2 X 1 Σ 1 2 Σ 1 2 A T ) 1 = ( B W 1 B T ) 1 .

Write B = M ( I q , 0 ) G , where M is a q × q nonsingular matrix and G is an m × m orthogonal matrix. Now,

( B W 1 B T ) 1 = M I q 0 G W 1 G T I q 0 T M T 1 = ( M T ) 1 I q 0 Y 1 I q 0 1 M 1 = ( M T ) 1 ( Y 1,1 ) 1 M 1 ,

where Y = Y 1,1 Y 1,2 Y 2,1 Y 2,2 = G W G T IL ( m , α ) , Y 1,1 is a q × q matrix, and Y 1,1 = ( Y 1,1 Y 1,2 Y 2,2 1 Y 2,1 ) 1 = Y 1,1 2 1 . From Theorem 4.4, Y 1,1 2 B2 ( q , α ( m q ) 2 , ( m + 1 ) 2 ) and therefore, Z = ( M T ) 1 Y 1,1 2 M 1 has the probability density function proportional to

Γ q [ α + ( q + 1 ) 2 ] Γ q [ α ( m q ) 2 ] Γ q [ ( m + 1 ) 2 ] det ( M M T ) ( m + 1 ) 2 det ( Z ) α ( m + 1 ) 2 det ( ( M M T ) 1 + Z ) α ( q + 1 ) 2 , Z > 0 .

Now, noting that M M T = B B T = A Σ 1 A T and transforming S = ( A Σ 1 A T ) 1 2 Z ( A Σ 1 A T ) 1 2 with the Jacobian J ( Z S ) = det ( A Σ 1 A T ) ( q + 1 ) 2 in the above density, we obtain the desired result. The proof of the second part is similar.□

Corollary 4.5.1

If X IL ( m , Σ , α ) and a R m , a 0 , then

a T Σ 1 a a T X 1 a B2 α m 1 2 , m + 1 2 ,

a T X 1 a a T Σ 1 a B2 m + 1 2 , α m 1 2 ,

a T X 1 a a T ( Σ 1 + X 1 ) a B1 m + 1 2 , α m 1 2 ,

a T Σ 1 a a T ( Σ 1 + X 1 ) a B1 α m 1 2 , m + 1 2 ,

and

a T X a a T Σ a IL ( α ) .

Let X = ( x i , j ) , X 1 = ( x i , j ) , Σ = ( σ i , j ) , and Σ 1 = ( σ i , j ) and a = ( a 1 , , a m ) T . Then, by setting a i = 1 and a k = 0 , for k i , 1 j , k m , in Corollary 4.5.1, we obtain

σ i , i x i , i B2 α m 1 2 , m + 1 2

and

x i , i σ i , i IL ( α ) .

Since, for i = 1 , , m , x i , i σ i , i IL ( α ) , the joint distribution of x 1,1 σ 1,1 , , x m , m σ m , m is a multivariate inverse Lomax or multivariate inverse Pareto distribution of the second kind. Also, note that the joint distribution of σ 1,1 x 1,1 , , σ m , m x m , m is a multivariate beta type 2 distribution.

In Corollary 4.5.1, the distributions of a T Σ 1 a a T X 1 a , a T X 1 a a T Σ 1 a , a T X 1 a a T ( Σ 1 + X 1 ) a , a T Σ 1 a a T ( Σ 1 + X 1 ) a , and a T X a a T Σ a do not depend on a . Thus, if y is an m × 1 random vector, independent of X , and P ( y 0 ) = 1 , then we see that

y T Σ 1 y y T X 1 y B2 α m 1 2 , m + 1 2 ,

y T X 1 y y T Σ 1 y B2 m + 1 2 , α m 1 2 ,

y T X 1 y y T ( Σ 1 + X 1 ) y B1 m + 1 2 , α m 1 2 ,

y T Σ 1 y y T ( Σ 1 + X 1 ) y B1 α m 1 2 , m + 1 2 ,

and

y T X y y T Σ y IL ( α ) .

5 Cholesky decomposition of the inverse Lomax matrix

The Cholesky decomposition, also known as Cholesky factorization, is a process of writing a Hermitian positive-definite matrix into the product of a lower (upper) triangular matrix and its conjugate transpose. In this section, we will consider Cholesky decompositions of inverted Lomax matrix and develop distributional result for associated triangular matrices.

Theorem 5.1

If W IL ( m , α ) and W = U U T , where U = ( u i , j ) is a lower triangular matrix with positive diagonal elements, then u 1,1 , , u m , m are all independent, u i , i 2 B2 ( α ( i 1 ) 2 , ( i + 1 ) 2 ) , i = 1 , , m .

Proof

Making the transformation W = U U T with the Jacobian J ( W U ) = 2 m i = 1 m u i , i m + 1 i in (10), the density of U can be derived as

(14) 2 m Γ m [ α + ( m + 1 ) 2 ] Γ m ( α ) Γ m [ ( m + 1 ) 2 ] det ( U U T ) α ( m + 1 ) 2 det ( I m + U U T ) α ( m + 1 ) 2 i = 1 m u i , i m + 1 i ,

where u i , i > 0 , i = 1 , , m and < u i , j < for 1 j < i m . Now, partition U as

U = u 1,1 0 u 1 U ( 2 ) ,

where u 1 is an ( m 1 ) × 1 vector and U ( 2 ) is an ( m 1 ) × ( m 1 ) lower triangular matrix. Then,

det ( I m + U U T ) = det 1 + u 1,1 2 u 1,1 u 1 T u 1,1 u 1 I m 1 + u 1 u 1 T + U ( 2 ) U ( 2 ) T = ( 1 + u 1,1 2 ) det ( I m 1 + U ( 2 ) U ( 2 ) T ) 1 + 1 1 + u 1,1 2 u 1 T ( I m 1 + U ( 2 ) U ( 2 ) T ) 1 u 1 .

Now, make the transformation

y 1 = 1 ( 1 + u 1,1 2 ) 1 2 ( I m 1 + U ( 2 ) U ( 2 ) T ) 1 2 u 1

with the Jacobian J ( u 1 y 1 ) = ( 1 + u 1,1 2 ) ( m 1 ) 2 det ( I m 1 + U ( 2 ) U ( 2 ) T ) 1 2 in (14) to obtain the joint density of u 1,1 , U ( 2 ) and y 1 as

2 m Γ m [ α + ( m + 1 ) 2 ] Γ m ( α ) Γ m [ ( m + 1 ) 2 ] u 1,1 2 α 1 ( 1 + u 1,1 2 ) α 1 det ( U ( 2 ) U ( 2 ) T ) α ( m + 1 ) 2 det ( I m 1 + U ( 2 ) U ( 2 ) T ) α m 2 i = 2 m u i , i m + 1 i ( 1 + y 1 T y 1 ) α ( m + 1 ) 2 .

From the above factorization, we see that u 1,1 , U ( 2 ) , and y 1 are all independent, u 1,1 2 B2 ( α , 1 ) and the density of U ( 2 ) is proportional to

(15) det ( U ( 2 ) U ( 2 ) T ) α ( m + 1 ) 2 det ( I m 1 + U ( 2 ) U ( 2 ) T ) α m 2 i = 2 m u i , i m + 1 i .

Next, partition U ( 2 ) as

U ( 2 ) = u 2,2 0 u 2 U ( 3 ) ,

where u 2 is an ( m 2 ) × 1 vector and U ( 3 ) is an ( m 2 ) × ( m 2 ) lower triangular matrix. Then,

det ( I m + U ( 2 ) U ( 2 ) T ) = det 1 + u 2,2 2 u 2,2 u 2 T u 2,2 u 2 I m 2 + u 2 u 2 T + U ( 3 ) U ( 3 ) T = ( 1 + u 2,2 2 ) det ( I m 2 + U ( 3 ) U ( 3 ) T ) 1 + 1 1 + u 2,2 2 u 2 T ( I m 2 + U ( 3 ) U ( 3 ) T ) 1 u 2 .

Now, make the transformation

y 2 = 1 ( 1 + u 2,2 2 ) 1 2 ( I m 2 + U ( 3 ) U ( 3 ) T ) 1 2 u 2

with the Jacobian J ( u 2 y 2 ) = ( 1 + u 2,2 2 ) ( m 2 ) 2 det ( I m 2 + U ( 3 ) U ( 3 ) T ) 1 2 in (15) to obtain the joint density of u 2,2 , U ( 3 ) , and y 2 proportional to

u 2,2 2 α 2 ( 1 + u 2,2 2 ) α 1 det ( U ( 3 ) U ( 3 ) T ) α ( m + 1 ) 2 det ( I m 1 + U ( 3 ) U ( 3 ) T ) α ( m 1 ) 2 i = 3 m u i , i m + 1 i ( 1 + y 2 T y 2 ) α m 2 .

From the above factorization, we see that u 2,2 , U ( 3 ) , and y 2 are all independent, u 2,2 2 B2 ( α 1 2,3 2 ) . Continuing further with similar steps, we can show that u 1,1 2 , , u m , m 2 are independent, u i , i 2 B2 ( α ( i 1 ) 2 , ( i + 1 ) 2 ) , i = 1 , , m .□

Corollary 5.1.1

If W IL ( m , α ) , then the distribution of det ( W ) is same as the distribution of the product of m independent beta type 2 variables. That is, det ( W ) i = 1 m z i , where z i B2 ( α ( i 1 ) 2 , ( i + 1 ) 2 ) , i = 1 , , m .

Corollary 5.1.2

If W IL ( m , α ) , then

det ( W ( 1 ) ) det ( W ( 0 ) ) , det ( W ( 2 ) ) det ( W ( 1 ) ) , , det ( W ( m ) ) det ( W ( m 1 ) ) ( det ( W ( 0 ) ) 1 )

are independently distributed. Further, for i = 1 , , m , det ( W ( i ) ) det ( W ( i 1 ) ) B2 ( α ( i 1 ) 2 , ( i + 1 ) 2 ) .

Theorem 5.2

If W IL ( m , α ) and W = V V T , where V = ( v i , j ) is an upper triangular matrix with positive diagonal elements, then v 1,1 , , v m , m are all independent, v i , i 2 B2 ( α ( m i ) 2 , ( m i + 2 ) 2 ) , i = 1 , , m .

Proof

Making the transformation W = V V T with the Jacobian J ( W V ) = 2 m i = 1 m v i , i i in (10), the density of V can be derived as

(16) 2 m Γ m [ α + ( m + 1 ) 2 ] Γ m ( α ) Γ m [ ( m + 1 ) 2 ] det ( V V T ) α ( m + 1 ) 2 det ( I m + V V T ) α ( m + 1 ) 2 i = 1 m v i , i i ,

where v i , i > 0 , i = 1 , , m and < v i , j < for 1 i < j m . Now, partition V as

V = V ( m 1 ) v 1 0 v m , m ,

where v 1 is an ( m 1 ) × 1 vector and V ( m 1 ) is an ( m 1 ) × ( m 1 ) upper triangular matrix. Then,

det ( I m + V V T ) = det I m 1 + v 1 v 1 T + V ( m 1 ) ( V ( m 1 ) ) T v m , m v 1 v m , m v T 1 + v m , m 2 = ( 1 + v m , m 2 ) det ( I m 1 + V ( m 1 ) ( V ( m 1 ) ) T ) 1 + 1 1 + v m , m 2 v 1 T ( I m 1 + V ( m 1 ) ( V ( m 1 ) ) T ) 1 v 1 .

Now, make the transformation

y 1 = 1 ( 1 + v m , m 2 ) 1 2 ( I m 1 + V ( m 1 ) ( V ( m 1 ) ) T ) 1 2 v 1

with the Jacobian J ( v 1 y 1 ) = ( 1 + v m , m 2 ) ( m 1 ) 2 det ( I m 1 + V ( m 1 ) ( V ( m 1 ) ) T ) 1 2 in (16) to obtain the joint density of V ( m 1 ) , y 1 , and v m m as

2 m Γ m [ α + ( m + 1 ) 2 ] Γ m ( α ) Γ m [ ( m + 1 ) 2 ] v m , m 2 α 1 ( 1 + v m , m 2 ) α 1 det ( V ( m 1 ) ( V ( m 1 ) ) T ) α ( m + 1 ) 2 det ( I m 1 + V ( m 1 ) ( V ( m 1 ) ) T ) α m 2 i = 1 m 1 v i , i i ( 1 + y 1 T y 1 ) α ( m + 1 ) 2 .

From the above factorization, we see that V ( m 1 ) , v m , m , and y 1 are all independent, v m , m 2 B2 ( α , 1 ) and the density of V ( m 1 ) is proportional to

(17) det ( V ( m 1 ) ( V ( m 1 ) ) T ) α ( m + 1 ) 2 det ( I m 1 + V ( m 1 ) ( V ( m 1 ) ) T ) α m 2 i = 1 m 1 v i , i i .

Now, partition V ( m 1 ) as

V ( m 1 ) = V ( m 2 ) v 2 0 v m 1 , m 1 ,

where v 2 is an ( m 2 ) × 1 vector and V ( m 2 ) is an ( m 2 ) × ( m 2 ) upper triangular matrix. Then,

det ( I m + V ( m 1 ) ( V ( m 1 ) ) T ) = det I m 2 + v 2 v 2 T + V ( m 2 ) ( V ( m 2 ) ) T v m 1 , m 1 v 2 v m 1 , m 1 v 2 T 1 + v m 1 , m 1 2 = ( 1 + v m 1 , m 1 2 ) det ( I m 2 + V ( m 2 ) ( V ( m 2 ) ) T ) 1 + 1 1 + v m 1 , m 1 2 v 2 T ( I m 2 + V ( m 2 ) ( V ( m 2 ) ) T ) 1 v 2 .

Now, make the transformation

y 2 = 1 ( 1 + v m 1 , m 1 2 ) 1 2 ( I m 2 + V ( m 2 ) ( V ( m 2 ) ) T ) 1 2 v 2

with the Jacobian J ( v 2 y 2 ) = ( 1 + v m 1 , m 1 2 ) ( m 2 ) 2 det ( I m 2 + V ( m 2 ) ( V ( m 2 ) ) T ) 1 2 in (17) to obtain the joint density of V ( m 2 ) , y 2 , and v m 1 , m 1 proportional to

v m 1 , m 1 2 α 2 ( 1 + v m 1 , m 1 2 ) α 1 det ( V ( m 2 ) ( V ( m 2 ) ) T ) α ( m + 1 ) 2 det ( I m 2 + V ( m 2 ) ( V ( m 2 ) ) T ) α ( m 1 ) 2 i = 1 m 2 v i , i i ( 1 + y 2 T y 2 ) α m 2 .

From the above factorization, we see that V ( m 2 ) , v m 1 , m 1 , and y 2 are all independent, v m 1 , m 1 2 B2 ( α 1 2 , 3 2 ) . Repeating the argument given above on the density function of V ( m 2 ) = V ( m 3 ) v 3 0 v m 2 , m 2 , we observe that v m 2 , m 2 2 B2 ( α 1 , 2 ) and is independent of v m m and v m 1 , m 1 . Continuing further with the same argument, we obtain the desired result.□

Corollary 5.2.1

If W IL ( m , α ) , then the distribution of det ( W ) is the same as the distribution of the product of m independent beta type 2 variables. That is, det ( W ) i = 1 m z i , where z i B2 ( α ( m i ) 2 , ( m i + 2 ) 2 ) , i = 1 , , m .

Corollary 5.2.2

If W IL ( m , α ) , then

det ( W ( 1 ) ) det ( W ( 2 ) ) , det ( W ( 2 ) ) det ( W ( 3 ) ) , , det ( W ( m ) ) det ( W ( m + 1 ) ) ( det ( W ( m + 1 ) ) 1 )

are independently distributed. Further, for i = 1 , , m , det ( W ( i ) ) det ( W ( i + 1 ) ) B2 ( α ( m i ) 2 , ( m i + 2 ) 2 ) .

6 Expected values

The expected value of X 1 , for X IL ( m , Σ , α ) , can easily be obtained from Corollary 4.5.1. For any fixed a R m , a 0 ,

E a T X 1 a a T Σ 1 a = E ( v ) ,

where v B2 ( ( m + 1 ) 2 , α ( m 1 ) 2 ) . Hence, for all c R m ,

a T E ( X 1 ) a = a T Σ 1 a E ( v ) = m + 1 2 α m 1 a T Σ 1 a , α > m + 1 2 .

Since the above result is valid for all a R m , we obtain

E ( X 1 ) = m + 1 2 α m 1 Σ 1 , α > m + 1 2 .

Theorem 6.1

Let W IL ( m , α ) , then

E det ( W ) r det ( I m + W ) s = Γ m [ α + ( m + 1 ) 2 ] Γ m ( α + r ) Γ m [ ( m + 1 ) 2 + s r ] Γ m ( α ) Γ m [ ( m + 1 ) 2 ] Γ m [ α + s + ( m + 1 ) 2 ] ,

where Re ( s r ) 0 and Re ( r ) 0 .

Proof

By definition

E det ( W ) r det ( I m + W ) s = Γ m [ α + ( m + 1 ) 2 ] Γ m ( α ) Γ m [ ( m + 1 ) 2 ] W > 0 det ( W ) α + r ( m + 1 ) 2 d W det ( I m + W ) α + ( m + 1 ) 2 + s .

Now, evaluating the above integral using (5), we obtain the result.□

Corollary 6.1.1

If W IL ( m , α ) , then

E det ( W ) r det ( I m + W ) r = Γ m [ α + ( m + 1 ) 2 ] Γ m ( α + r ) Γ m ( α ) Γ m [ α + r + ( m + 1 ) 2 ] , Re ( r ) 0

and

E [ det ( I m + W ) s ] = Γ m [ α + ( m + 1 ) 2 ] Γ m [ ( m + 1 ) 2 + s ] Γ m [ ( m + 1 ) 2 ] Γ m [ α + ( m + 1 ) 2 + s ] , Re ( s ) 0 .

By writing multivariate gamma functions in terms of ordinary gamma functions, expressions E det ( W ) r det ( I m + W ) r and E [ det ( I m + W ) s ] can be simplified as

(18) E det ( W ) r det ( I m + W ) r = j = 1 m Γ [ α + ( m + 1 ) 2 ( j 1 ) 2 ] Γ ( α + r ( j 1 ) 2 ) Γ [ α + ( m + 1 ) 2 ( j 1 ) 2 + r ] Γ [ α ( j 1 ) 2 ]

and

(19) E [ det ( I m + W ) s ] = j = 1 m Γ [ α + ( m + 1 ) 2 ( j 1 ) 2 ] Γ [ ( m + 1 ) 2 + s ( j 1 ) 2 ] Γ [ α + ( m + 1 ) 2 + s ( j 1 ) 2 ] Γ [ ( m + 1 ) 2 ( j 1 ) 2 ] .

Substituting r , s = 1 , 2 in (18) and (19), the first and second order moments of det ( W ) det ( I m + W ) and det ( I m + W ) 1 can be calculated as

E det ( W ) det ( I m + W ) = ( 2 α m + 1 ) m ( 2 α + 2 ) m ,

E det ( W ) 2 det ( I m + W ) 2 = ( 2 α m + 1 ) m ( 2 α m + 3 ) m ( 2 α + 2 ) m ( 2 α + 4 ) m ,

E [ det ( I m + W ) 1 ] = ( 2 ) m ( 2 α + 2 ) m ,

and

E [ det ( I m + W ) 2 ] = ( 2 ) m ( 4 ) m ( 2 α + 2 ) m ( 2 α + 4 ) m ,

where the Pochhammer notation ( a ) k is defined by ( a ) k = a ( a + 1 ) ( a + k 1 ) , k = 1 , 2 , with ( a ) 0 = 1 .

Finally, we give the following result used in the derivation of Kullback-Liebler divergence.

Theorem 6.2

If W IL ( m , Σ , α ) , then

(20) E ln det ( W ) det ( I m + Σ 1 W ) = ln det ( Σ ) + i = 1 m ψ α i 1 2 ψ α + i + 1 2 ,

where ψ is the digamma function ψ ( z ) = d d z ln Γ ( z ) .

Proof

By definition,

E ln det ( W ) det ( I m + Σ 1 W ) = Γ m [ α + ( m + 1 ) 2 ] Γ m ( α ) Γ m [ ( m + 1 ) 2 ] det ( Σ ) α W > 0 ln det ( W ) det ( I m + Σ 1 W ) det ( W ) α ( m + 1 ) 2 det ( I m + Σ 1 W ) α + ( m + 1 ) 2 d W .

Now, writing

ln det ( W ) det ( I m + Σ 1 W ) det ( W ) det ( I m + Σ 1 W ) α = d d α det ( W ) det ( I m + Σ 1 W ) α

in the above expression, we have

E ln det ( W ) det ( I m + Σ 1 W ) = Γ m [ α + ( m + 1 ) 2 ] Γ m ( α ) Γ m [ ( m + 1 ) 2 ] det ( Σ ) α d d α W > 0 det ( W ) α ( m + 1 ) 2 det ( I m + Σ 1 W ) α + ( m + 1 ) 2 d W = Γ m [ α + ( m + 1 ) 2 ] det ( Σ ) α Γ m ( α ) d d α det ( Σ ) α Γ m ( α ) Γ m [ α + ( m + 1 ) 2 ] ,

where the last line has been obtained by using (5). Now, we obtain

d d α det ( Σ ) α Γ m ( α ) Γ m [ α + ( m + 1 ) 2 ] = det ( Σ ) α Γ m ( α ) Γ m [ α + ( m + 1 ) 2 ] ln det ( Σ ) + i = 1 m ψ α i 1 2 ψ α + i + 1 2 .

Now, substituting appropriately we obtain the desired result.□

7 Kullback-Liebler divergence

Assume two matrix variate inverse Lomax distributions L 1 and L 2 specified by the probability laws L 1 : W IL ( m , Σ , α 1 ) and L 2 : W IL ( m , Σ , α 2 ) . The Kullback-Liebler divergence of L 1 from L 2 is given by

KL ( L 1 ; L 2 ) = W > 0 ln f MIL ( W ; Σ , α 1 ) f MIL ( W ; Σ , α 2 ) f MIL ( W ; Σ , α 1 ) d W ,

where f MIL ( W ; Σ , α 1 ) and f MIL ( W ; Σ , α 2 ) are matrix variate inverse Lomax densities. Substituting appropriately from (10), the above expression can be re-written as

KL ( L 1 ; L 2 ) = W > 0 ln Γ m [ α 1 + ( m + 1 ) 2 ] Γ m ( α 2 ) Γ m [ α 2 + ( m + 1 ) 2 ] Γ m ( α 1 ) det ( Σ 1 W ) det ( I m + Σ 1 W ) α 1 α 2 f MIL ( W ; Σ , α 1 ) d W = ln Γ m [ α 1 + ( m + 1 ) 2 ] Γ m ( α 2 ) Γ m [ α 2 + ( m + 1 ) 2 ] Γ m ( α 1 ) det ( Σ ) α 2 α 1 + ( α 1 α 2 ) W > 0 ln det ( W ) det ( I m + Σ 1 W ) f MIL ( W ; Σ , α 1 ) d W = ln Γ m [ α 1 + ( m + 1 ) 2 ] Γ m ( α 2 ) Γ m [ α 2 + ( m + 1 ) 2 ] Γ m ( α 1 ) det ( Σ ) α 2 α 1 + ( α 1 α 2 ) E W ln det ( W ) det ( I m + Σ 1 W ) ,

where W IL ( m , Σ , α 1 ) . Finally, application of (20) in the above expression yields

KL ( L 1 ; L 2 ) = ln Γ m [ α 1 + ( m + 1 ) 2 ] Γ m ( α 2 ) Γ m [ α 2 + ( m + 1 ) 2 ] Γ m ( α 1 ) det ( Σ ) α 2 α 1 + ( α 1 α 2 ) ln det ( Σ ) + i = 1 m ψ α 1 i 1 2 ψ α 1 + i + 1 2 .

8 Entropy

The variation in uncertainty of a random variable is measured by its entropy. In this section, we derive expressions for Shannon and Rényi entropies for the matrix variate inverse Lomax distribution.

The well-known Shannon entropy, denoted by H S H ( f ) , introduced in Shannon [26] is defined by

(21) H S H ( f ) = f ( x ) ln f ( x ) d x .

One of the main extensions of the Shannon entropy was defined by Rényi [25]. This generalized entropy measure is given by

(22) H R ( η , f ) = ln G ( η ) 1 η

for η > 0 and η 1 , where

G ( η ) = f η ( x ) d x .

The additional parameter η is used to describe complex behavior in probability models and the associated process under study. Rényi entropy is monotonically decreasing in η , while Shannon entropy (21) is obtained from (22) for η 1 . For details, refer [22,34,35].

Theorem 8.1

For the matrix variate inverse Lomax distribution defined by the probability density function (10), the Rényi and the Shannon entropies are given by

H R ( η , f ) = 1 1 η η ln Γ m [ α + ( m + 1 ) 2 ] Γ m ( α ) Γ m [ ( m + 1 ) 2 ] ( η 1 ) ( m + 1 ) 2 ln det ( Σ ) + ln Γ m η α ( η 1 ) ( m + 1 ) 2 + ln Γ m ( 2 η 1 ) ( m + 1 ) 2 ln Γ m η α + m + 1 2

and

H S H ( f ) = ln Γ m [ α + ( m + 1 ) 2 ] det ( Σ ) ( m + 1 ) 2 Γ m ( α ) Γ m [ ( m + 1 ) 2 ] i = 1 m α m + 1 2 ψ α i 1 2 + ( m + 1 ) ψ m i + 2 2 α + m + 1 2 ψ α + m i + 2 2 ,

respectively.

Proof

For η > 0 and η 1 , using the density given in (10), we have

(23) G ( η ) = Γ m [ α + ( m + 1 ) 2 ] det ( Σ ) α Γ m ( α ) Γ m [ ( m + 1 ) 2 ] η W > 0 det ( W ) η ( α ( m + 1 ) 2 ) det ( I m + Σ 1 W ) η ( α + ( m + 1 ) 2 ) d W = Γ m [ α + ( m + 1 ) 2 ] det ( Σ ) α Γ m ( α ) Γ m [ ( m + 1 ) 2 ] η Γ m [ η α ( η 1 ) ( m + 1 ) 2 ] Γ m [ ( 2 η 1 ) ( m + 1 ) 2 ] det ( Σ ) η α + ( η 1 ) ( m + 1 ) 2 Γ m [ η ( α + ( m + 1 ) 2 ) ] ,

where the last line has been obtained by using beta integral (5). Now, taking logarithm of G ( η ) given in (23) and using (22), we obtain H R ( η , f ) . The Shannon entropy can be obtained from H R ( η , f ) by taking η 1 and using L’Hopital’s rule.□

9 Eigenvalues of the inverse Lomax matrix

In this section, we derive densities of eigenvalues of the random matrix distributed as inverse Lomax. We begin by stating the following result, which is essential in determining the density of eigenvalues.

Theorem 9.1

Let A be a positive definite random matrix of order m with the probability density function f(A). Then, the joint probability density function of eigenvalues l 1 , l 2 , , l m of A is given by

(24) π m 2 2 Γ m ( m 2 ) i < j m ( l i l j ) O ( m ) f ( H L H ) [ d H ] , ( l 1 > l 2 > > l m > 0 ) ,

where L = diag ( l 1 , l 2 , , l m ) and [ d H ] is the unit invariant Haar measure on the group of orthogonal matrices.

If A and B are symmetric positive definite and symmetric matrices, respectively, of order m with A 1 2 B A 1 2 < 1 , then we can expand F 0 1 ( α ; A H B H ) = det ( I m A H B H ) α in terms of zonal polynomials as

(25) det ( I m A H B H ) α = k = 0 κ k ( α ) κ C κ ( A H B H ) k ! ,

where κ k denotes summation over all ordered partitions κ , κ = ( k 1 , , k m ) , k 1 k m 0 and k 1 + + k m = k , ( α ) κ is the generalized hypergeometric coefficient and C κ ( X ) denotes the zonal polynomial of the symmetric matrix X corresponding to the ordered partition κ . Now, integrating (25) over the group of orthogonal matrices

(26) O ( m ) det ( I m A H B H ) α [ d H ] = k = 0 κ k ( α ) κ k ! O ( m ) C κ ( A H B H ) [ d H ] = k = 0 κ k ( α ) κ k ! C κ ( A ) C κ ( B ) C κ ( I m ) = F 0 ( m ) 1 ( α ; A , B ) ,

where F 0 ( m ) 1 ( α ; A , B ) is the hypergeometric function of two matrix arguments and the integral has been evaluated by using the fundamental property of zonal polynomials given by

(27) O ( m ) C κ ( A H B H ) [ d H ] = C κ ( A ) C κ ( B ) C κ ( I m ) ,

where A is a symmetric positive definite and B is a symmetric matrix of order m . Also, if one of the argument matrices is proportional to the identity matrix F 0 ( m ) 1 ( α ; A , B ) reduces to a one argument function. That is, if A = λ I m , then

F 0 ( m ) 1 ( α ; λ I m , B ) = F 0 ( m ) 1 ( α ; λ B ) = det ( I m λ B ) α .

Proof of Theorem 9.1 and several other results such as (27) can be found in [8,16,20].

Theorem 9.2

If W IL ( m , Σ , α ) , then the joint probability density function of the eigenvalues w 1 , , w m of W is given by

(28) π m 2 2 Γ m [ α + ( m + 1 ) 2 ] det ( Σ ) ( m + 1 ) 2 Γ m ( α ) Γ m [ ( m + 1 ) 2 ] Γ m ( m 2 ) i < j m ( w i w j ) i = 1 m [ w i α ( m + 1 ) 2 ( 1 + w i ) α ( m + 1 ) 2 ] F 0 ( m ) 1 α + m + 1 2 ; I m Σ , ( I m + L ) 1

for I m Σ < 1 and for I m Σ 1 < 1 ,

(29) π m 2 2 Γ m [ α + ( m + 1 ) 2 ] det ( Σ ) α Γ m ( α ) Γ m [ ( m + 1 ) 2 ] Γ m ( m 2 ) i < j m ( w i w j ) i = 1 m [ w i α ( m + 1 ) 2 ( 1 + w i ) α ( m + 1 ) 2 ] F 0 ( m ) 1 α + m + 1 2 ; I m Σ 1 , L ( I m + L ) 1 ,

where 0 < w m < < w 1 < and L = diag ( w 1 , , w m ) .

Proof

The probability density function of W is given by (10). Applying Theorem 9.1, we obtain the joint probability density function of the eigenvalues w 1 , , w m of W as

(30) Γ m [ α + ( m + 1 ) 2 ] det ( Σ ) α Γ m ( α ) Γ m [ ( m + 1 ) 2 ] π m 2 2 Γ m ( m 2 ) i < j m ( w i w j ) i = 1 m w i α ( m + 1 ) 2 O ( m ) det ( I m + Σ 1 H L H ) α ( m + 1 ) 2 [ d H ] .

Writing

det ( I m + Σ 1 H L H ) = det ( Σ ) 1 det ( Σ + H L H T ) = det ( Σ ) 1 det ( I m + H L H T ( I m Σ ) ) = det ( Σ ) 1 det ( I m + L ) det ( I m ( I m Σ ) H ( I m + L ) 1 H T )

if I m Σ < 1 and

det ( I m + Σ 1 H L H ) = det ( I m + H L H T H L H T + Σ 1 H L H T ) = det ( I m + H L H T ( I m Σ 1 ) H L H T ) = det ( I m + L ) det ( I m ( I m Σ 1 ) H L ( I m + L ) 1 H T )

if I m Σ 1 < 1 in the above integral and applying (26), we evaluate the above integral as

(31) O ( m ) det ( I m + Σ 1 H L H ) α ( m + 1 ) 2 [ d H ] = det ( Σ ) α + ( m + 1 ) 2 det ( I m + L ) α ( m + 1 ) 2 F 0 ( m ) 1 α + m + 1 2 ; I m Σ , ( I m + L ) 1

and

(32) O ( m ) det ( I m + Σ 1 H L H ) α ( m + 1 ) 2 [ d H ] = det ( I m + L ) α ( m + 1 ) 2 F 0 ( m ) 1 α + m + 1 2 ; I m Σ 1 , L ( I m + L ) 1

for I m Σ < 1 and I m Σ 1 < 1 , respectively. Finally, substituting (31) and (32) in (30), we obtain the desired result.□

Corollary 9.2.1

If W IL ( m , α ) , then the joint probability density function of the eigenvalues w 1 , , w m of W is given by

π m 2 2 Γ m [ α + ( m + 1 ) 2 ] Γ m ( α ) Γ m [ ( m + 1 ) 2 ] Γ m ( m 2 ) i < j m ( w i w j ) i = 1 m [ w i α ( m + 1 ) 2 ( 1 + w i ) α ( m + 1 ) 2 ] ,

where 0 < w m < < w 1 < .

Since the density over its support set integrates to one, from (28) and (29), we obtain

0 < w m < < w 1 < i < j m ( w i w j ) i = 1 m [ w i α ( m + 1 ) 2 ( 1 + w i ) α ( m + 1 ) 2 ] F 0 ( m ) 1 α + m + 1 2 ; I m Σ , ( I m + L ) 1 i = 1 m d w i = Γ m ( α ) Γ m [ ( m + 1 ) 2 ] Γ m ( m 2 ) π m 2 2 Γ m [ α + ( m + 1 ) 2 ] det ( Σ ) ( m + 1 ) 2

for I m Σ < 1 and for I m Σ 1 < 1 ,

0 < w m < < w 1 < i < j m ( w i w j ) i = 1 m [ w i α ( m + 1 ) 2 ( 1 + w i ) α ( m + 1 ) 2 ] F 0 ( m ) 1 α + m + 1 2 ; I m Σ 1 , L ( I m + L ) 1 i = 1 m d w i = Γ m ( α ) Γ m [ ( m + 1 ) 2 ] Γ m ( m 2 ) π m 2 2 Γ m [ α + ( m + 1 ) 2 ] det ( Σ ) α .

10 Estimation

Let W 1 , , W n be a random sample from the matrix variate inverse distribution defined by the density (10). The log likelihood function, denoted by l ( α , Σ ) , is given by

l ( α , Σ ) = n ln Γ m α + m + 1 2 ln Γ m ( α ) ln Γ m m + 1 2 + m + 1 2 ln det ( Σ ) + α m + 1 2 i = 1 n ln det ( W i ) α + m + 1 2 i = 1 n ln det ( Σ + W i ) .

Now, differentiating l ( α , Σ ) appropriately, we obtain

(33) l ( α , Σ ) α = n i = 1 m ψ α + m i + 2 2 ψ α i 1 2 + i = 1 n ln det ( W i ) i = 1 n ln det ( Σ + W i )

and (using equation 17.46 of [27, p. 368])

(34) l ( α , Σ ) Σ = n ( m + 1 ) 2 [ 2 Σ 1 diag ( Σ 1 ) ] α + m + 1 2 i = 1 n [ 2 ( Σ + W i ) 1 diag ( Σ + W i ) 1 ] = 2 B diag ( B ) ,

where

B = n ( m + 1 ) 2 Σ 1 α + m + 1 2 i = 1 n ( Σ + W i ) 1 .

Setting (34) to zero results in B = 0 and

(35) n ( m + 1 ) 2 Σ ^ 1 = α ^ + m + 1 2 i = 1 n ( Σ ^ + W i ) 1 .

Equating (33) to zero yields

(36) n i = 1 m ψ α ^ + m i + 2 2 ψ α ^ i 1 2 = i = 1 n ln det ( W i ) + i = 1 n ln det ( Σ ^ + W i ) .

Equations (35) and (36) can be easily solved since in-built routines for digamma function are widely available even in some pocket calculators. Simulations in Section 11 show that the solutions of (35) and (36) are unique for a wide range of initial values and for a wide range of sample sizes.

11 Simulation study

In this section, we perform simulations to assess the finite sample performance of the maximum likelihood estimators in Section 10 in terms of biases and mean squared errors. The assessment is based on the following scheme.

  1. set initial values for α and Σ .

  2. simulate a random sample of size n from (10).

  3. estimate the maximum likelihood estimates of α and Σ .

  4. repeat steps (ii) and (iii) 10,000 times, yielding the estimates ( α ^ i , Σ ^ i ) for i = 1, 2 , , 10,000 .

  5. compute the biases as

    bias ( α ^ ) = 1 10,000 i = 1 n ( α ^ i α ) , bias ( Σ ^ ) = 1 10,000 i = 1 n tr ( Σ ^ i Σ ) .

  6. compute the mean squared errors as

    mse ( α ^ ) = 1 10,000 i = 1 n ( α ^ i α ) 2 , mse ( Σ ^ ) = 1 10,000 i = 1 n [ tr ( Σ ^ i Σ ) ] 2 .

  7. repeat steps (ii)–(vi) for n = 10 , 11 , , 100 .

Plots of the biases vs n are shown in Figure 1. Plots of the mean squared errors vs n are shown in Figure 2.

Figure 1 
               Biases for 
                     
                        
                        
                           
                              
                                 α
                              
                              
                                 ^
                              
                           
                        
                        \widehat{\alpha }
                     
                   and 
                     
                        
                        
                           
                              
                                 Σ
                              
                              
                                 ^
                              
                           
                        
                        \widehat{\Sigma }
                     
                   vs 
                     
                        
                        
                           n
                        
                        n
                     
                  .
Figure 1

Biases for α ^ and Σ ^ vs n .

Figure 2 
               Mean squared errors for 
                     
                        
                        
                           
                              
                                 α
                              
                              
                                 ^
                              
                           
                        
                        \widehat{\alpha }
                     
                   and 
                     
                        
                        
                           
                              
                                 Σ
                              
                              
                                 ^
                              
                           
                        
                        \widehat{\Sigma }
                     
                   vs 
                     
                        
                        
                           n
                        
                        n
                     
                  .
Figure 2

Mean squared errors for α ^ and Σ ^ vs n .

The biases for α ^ and Σ ^ approach zero with increasing n . The biases appear negative for Σ ^ . The biases appear smaller in magnitude for α ^ . The mean squared errors for α ^ and Σ ^ decrease to zero with increasing n . The mean squared errors appear smaller for Σ ^ . It is pleasing that both the biases and mean squared errors appear reasonably small for all n 100 .

In our simulations, we took α = 2 and Σ = 1 0.1 0.1 1 . The results of simulations were similar for a wide range of other values of α and Σ . In particular, the biases always approached zero with increasing n , the mean squared errors always decreased to zero with increasing n , and both the biases and mean squared errors always appeared reasonably small for all n 100 .

Acknowledgments

The authors would like to thank the Editor and the two referees for careful reading and comments which greatly improved the article.

  1. Funding information: This research was supported by the Sistema Universitario de Investigación, Universidad de Antioquia [Project No. 2023-66310].

  2. Author contributions: All authors contributed equally to the manuscript.

  3. Conflict of interest: The authors have no conflicts of interest.

  4. Code availability: Not applicable.

  5. Consent to participate: Not applicable.

  6. Consent for publication: Not applicable.

  7. Data availability statement: Not applicable.

References

[1] S. M. Ahmed, M. Abdalla and A. A. Farghal, Statistical inference for competing risks model with type-II generalized hybrid censored of inverse Lomax distribution with applications, Alexandr. Eng. J. 118 (2025), 132–146, https://doi.org/10.1016/j.aej.2025.01.004. Search in Google Scholar

[2] I. Alam, M. Kamal, A. Rahman, Nayabuddin, A study on step stress partially accelerated life test under adaptive type-II progressive hybrid censoring for inverse Lomax distribution, Int. J. Reliabil. Safety 18 (2024), 1–18, https://doi.org/10.1504/IJRS.2024.139194. Search in Google Scholar

[3] R. Alharbi, R. Garg, I. Kumar, A. Kumari, and R. Aldallal, On estimation of P(Y<X) for inverse Pareto distribution based on progressively first failure censored data, PLoS One 18 (2023), e0287473, https://doi.org/10.1371/journal.pone.0287473. Search in Google Scholar PubMed PubMed Central

[4] F. E. Al-Muhayfith and A. Al-Bossly, Statistical inference of type-I hybrid censored inverse Lomax samples, Therm. Sci. 26 (2022), S339–S351, https://doi.org/10.2298/TSCI22S1339A. Search in Google Scholar

[5] A. I. Al-Omar, A. S. Hassan, N. Alotaibi, M. Shrahili, and H. F. Nagy, Reliability estimation of inverse Lomax distribution using extreme ranked set sampling, Adv. Math. Phys. 2021 (2021), no. 1, 4599872, https://doi.org/10.1155/2021/4599872. Search in Google Scholar

[6] Y. M. Amer and R. M. Shalabi, Optimal design of step stress partially accelerated life test under progressive type-II censored data with random removal for inverse Lomax distribution, Am. J. Appl. Math. Stat. 7 (2019), 171–177, https://doi.org/10.12691/ajams-7-5-3. Search in Google Scholar

[7] R. A. R. Bantan, M. Elgarhy, C. Chesneau and F. Jamal, Estimation of entropy for inverse Lomax distribution under multiple censored data, Entropy 22 (2020), 601, https://doi.org/10.3390/e22060601. Search in Google Scholar PubMed PubMed Central

[8] A. G. Constantine, Some noncentral distribution problems in multivariate analysis, Ann. Math. Stat. 34 (1963), 1270–1285, https://doi.org/10.1214/aoms/1177703863. Search in Google Scholar

[9] S. Dankunprasert, U. Jaroengeratikun, and T. Talangtam, The properties of inverse Pareto distribution and its application to extreme events, Thailand Stat. 19 (2021), 1–13. Search in Google Scholar

[10] R. Garg, M. Kumari, R. K. Sahoo, and A. Kumari, Stress-strength reliability estimation of multicomponent system with non-identical strength components from inverse Pareto distribution, Life Cycle Reliabil. Safety Eng. 13 (2024), 351–363, https://doi.org/10.1007/s41872-024-00258-6. Search in Google Scholar

[11] A. K. Gupta and D. K. Nagar, Properties of matrix variate beta type 3 distribution, Int. J. Math. Math. Sci. 2009 (2009), 308518, https://doi.org/10.1155/2009/308518. Search in Google Scholar

[12] A. K. Gupta and D. K. Nagar, Matrix Variate Distributions, Chapman and Hall/CRC, Boca Raton, 2000. https://djvu.online/file/Jovkf6iNhv3V5. Search in Google Scholar

[13] A. K. Gupta and D. K. Nagar, Matrix variate beta distribution, Int. J. Math. Math. Sci. 23 (2000), no. 7, 449–459. 10.1155/S0161171200002398Search in Google Scholar

[14] L. Guo and W. Gui, Bayesian and classical estimation of the inverse Pareto distribution and its application to strength-stress models, Am. J. Math. Manag. Sci. 37 (2018), 80–92, https://doi.org/10.1080/01966324.2017.1383217. Search in Google Scholar

[15] R. Hora, N. C. Kabdwal, and P. Srivastava, Classical and Bayesian inference for the inverse Lomax distribution under adaptive progressive type-II censored data with COVID-19 application, J. Reliabil. Stat. Stud. 15 (2022), 505–534, https://doi.org/10.13052/jrss0974-8024.1525. Search in Google Scholar

[16] A. T. James, Distributions of matrix variate and latent roots derived from normal samples, Ann. Math. Stat. 35 (1964), 475–501, 10.1214/aoms/1177703550. Search in Google Scholar

[17] U. Jan and S. P. Ahmad, Bayesian analysis of inverse Lomax distribution using approximation techniques, Math. Theory Model. 7 (2017). Search in Google Scholar

[18] N. L. Johnson, S. Kotz, and N. Balakrishnan, Continuous Univariate Distributions-2, 2nd edition, John Wiley and Sons, New York, 1994. Search in Google Scholar

[19] C. Kleiber, Lorenz ordering of order statistics from log-logistic and related distributions, J. Stat. Plan. Infer. 120 (2004), 13–19, https://doi.org/10.1016/S0378-3758(02)00495-0. Search in Google Scholar

[20] R. J. Muirhead, Aspects of Multivariate Statistical Theory, John Wiley and Sons, Inc., New York, 1982, https://doi.org/10.1002/9780470316559. Search in Google Scholar

[21] G. Mustafa, M. Ijaz, and F. Jamal, Order statistics of inverse Pareto distribution, Comput. J. Math. Stat. Sci. 1 (2022), 51–62, https://doi.org/10.21608/cjmss.2022.272724. Search in Google Scholar

[22] S. Nadarajah and K. Zografos, Expressions for Rényi and Shannon entropies for bivariate distributions, Inform. Sci. 170 (2005), 173–189, https://doi.org/10.1016/j.ins.2004.02.020. Search in Google Scholar

[23] J. Rahman, M. Aslam, and S. Ali, Estimation and prediction of inverse Lomax model via Bayesian approach, Caspian J. Appl. Sci. Res. 2 (2013), 43–56. Search in Google Scholar

[24] D. A. Ramadan, E. M. Almetwally, and A. H. Tolba, Statistical inference for multi stress-strength reliability based on progressive first failure with lifetime inverse Lomax distribution and analysis of transformer insulation data, Quality Reliabil Eng. Int. 39 (2023), 2558–2581, https://doi.org/10.1002/qre.3362. Search in Google Scholar

[25] A. Rényi, On measures of entropy and information, in: Proc. 4th Berkeley Sympos. Math. Statist. and Prob., Vol. I, Univ. California Press, Berkeley, Calif., 1961, pp. 547–561. Search in Google Scholar

[26] C. E. Shannon, A mathematical theory of communication, Bell Syst. Tech. J. 27 (1948), 379–423, 623–656, https://doi.org/10.1002/j.1538-7305.1948.tb01338.x.Search in Google Scholar

[27] G. A. F. Seber, A Matrix Handbook for Statisticians, John Wiley and Sons, Hoboken, NJ, 2008. 10.1002/9780470226797Search in Google Scholar

[28] S. K. Singh, U. Singh, and A. Singh Yadav, Reliability estimation for inverse Lomax distribution under type-II censored data using Markov chain Monte Carlo method, Int. J. Math. Stat. 17 (2016), 128–146. Search in Google Scholar

[29] S. K. Singh, U. Singh, and A. S. Yadav, Bayesian estimation of Lomax distribution under type-II hybrid censored data using Lindley’s approximation method, Int. J. Data Sci. 2 (2017), 352–268, 10.1504/IJDS.2017.10009048. Search in Google Scholar

[30] A. H. Tolba, D. A. Ramadan, E. M. Almetwally, T. M. Jawa, and N. Sayed-Ahmed, Statistical inference for stress-strength reliability using inverse Lomax lifetime distribution with mechanical engineering applications, Therm Sci. 26 (2022), S303–S326, 10.2298/TSCI22S1303T. Search in Google Scholar

[31] A. S. Yadav, S. K. Singh, and U. Singh, Bayesian estimation of P(Y<X) for inverse Lomax distribution under progressive type-II censoring scheme, Int J Syst Assurance Eng Manag. 10 (2019), 905–917, https://doi.org/10.1007/s13198-019-00820-x. Search in Google Scholar

[32] A. S. Yadav, S. K. Singh, and U. Singh, On hybrid censored inverse Lomax distribution: Application to the survival data, Statistica 76 (2016), 185–203, https://doi.org/10.6092/issn.1973-2201/5993. Search in Google Scholar

[33] A. S. Yadav, S. Shukla, M. Saha, D. Koley, and N. Jaiswal, On progressively type-II censored inverse Lomax distribution: Characterizations, estimation and application to cancer data, REVSTAT-Stat. J. (2024), https://doi.org/10.57805/revstat.vi.696. Search in Google Scholar

[34] K. Zografos, On maximum entropy characterization of Pearson’s type II and VII multivariate distributions, J. Multivariate Anal. 71 (1999), 67–75, https://doi.org/10.1006/jmva.1999.1824. Search in Google Scholar

[35] K. Zografos and S. Nadarajah, Expressions for Rényi and Shannon entropies for multivariate distributions, Stat. Probabil. Lett. 71 (2005), 71–84, https://doi.org/10.1016/j.spl.2004.10.023. Search in Google Scholar

Received: 2025-02-22
Revised: 2025-06-08
Accepted: 2025-07-02
Published Online: 2025-09-11

© 2025 the author(s), published by De Gruyter

This work is licensed under the Creative Commons Attribution 4.0 International License.

Downloaded on 1.10.2025 from https://www.degruyterbrill.com/document/doi/10.1515/spma-2025-0039/html?lang=en
Scroll to top button