
==== Front
Heliyon
Heliyon
Heliyon
2405-8440
Elsevier

S2405-8440(24)11879-8
10.1016/j.heliyon.2024.e35848
e35848
Research Article
Penalized estimation in parametric frailty model
Ahelali Marwan H. malhilaly@ut.edu.sa
a
Alamri Osama Abdulaziz oalmughamisi@ut.edu.sa
b
Sirohi Anu bsc.ashori@gmail.com
c⁎
a Department of Statistic, University of Tabuk, Tabuk-71491, Kingdom of Saudi Arabia
b Department of Statistics, University of Tabuk, Tabuk-71491, Kingdom of Saudi Arabia
c Department of Statistics, AIAS, Amity University, Noida, India
⁎ Corresponding author. bsc.ashori@gmail.com
08 8 2024
30 8 2024
08 8 2024
10 16 e3584810 11 2023
1 8 2024
5 8 2024
© 2024 Published by Elsevier Ltd.
2024

https://creativecommons.org/licenses/by-nc-nd/4.0/ This is an open access article under the CC BY-NC-ND license (http://creativecommons.org/licenses/by-nc-nd/4.0/).
Frailty model examines the effect of observable and non-observable factors on time to event data. Presence of collinearity produces unstable estimates of parameters. Therefore, this research focus on the penalized estimation of frailty model and proposed the new estimator which is the extension of ridge and principal component estimators. Simulation is run to reveal the performance of proposed estimator. Moreover, the technique is applied on NFHS (National Family Health Survey) data to examine the infant mortality in India.

Keywords

Frailty model
Collinearity
Ridge estimator
Principal component estimator
Infant mortality
==== Body
pmc1 Introduction

In survival analysis frailty model is a random effect model which analyzes the measured and unmeasured heterogeneity simultaneously [21]. It is the extension of proportional hazard regression model which was introduced by Cox [4]. Maximum likelihood estimation is used for the estimation of frailty model but, estimates become unstable, if explanatory variables are highly correlated to each other. In several real life situations, explanatory variables are not independent to each other therefore, alternative inference procedures have been proposed.

Hoerl and Kennard [5], Hoerl et al. [6], introduced the ridge estimator to overcome the problem of collinearity in linear regression model then Schafer et al. [18] extent this approach in logistic regression. Xue et al. [23] implemented the same procedure in Cox regression model furthermore, Albalawi et al. [3] extended it into the frailty model. Accuracy of continuous biomarker is checked through different estimators in different sampling situations, Mahdizadeh and Zamanzade [8], [9], Zamanzade et al. [10], Mahdizadeh and Zamanzade [12], [11]. Some other estimators are also developed in survival analysis, Ahmad et al. [2] and Michimae and Emura [13]. Ozkale and Arican [16] developed an estimator in logistic regression which is the combined form of ridge and principal component estimator, Sirohi and Rai [19], extent this estimator in proportional hazard regression model. But, this approach is not developed in frailty models.

This research generalizes the Ozkale and Arican's approach in frailty model. Proposed estimator is compared with maximum likelihood (ML), ridge and principal component estimators in respect of different situations. Further, the approach is utilized to examine the infant mortality in India.

2 Estimation in frailty model

Frailty model is a popular mixed model in survival analysis as it considers measurable and non-measurable factors simultaneously. Mathematical formulation of this model is represented as,(1) λ(τi|xi,νi)=νiλ0(τi)exp(γ′xi)

In equation (1), ith(i=1,2,...,n) individual's hazard is represented as λ(τi|xi,νi) at time τi. Explanatory variables (observable) are represented as the vector xi′=(xi1,xi2,...,xiq) and unobserved variable is represented through νi which is named as frailty variable. Base-line hazard is represented through the function λ0(τi) and regression coefficients are represented through the vector γ′=(γ1,γ2,...,γq).

Research is considering the parametric form of frailty model therefore it will have an assumption of the distribution regarding survival time which is taken as Weibull distribution with one parameter η, hence λ0(τi)=ητiη−1,η>0. Likelihood functional form of model-(1) is,(2) L(γ)=Πi=1n(νiλ0(τi)eγ′xi)Δiexp(−νiΨ0(τi)eγ′xi)

In equation (2), Ψ0(τi) represents the base-line cumulative hazard and Δi is the notation of censoring indicator having the value one (event has been occurred) and zero (observation is right censored). Oakes [14] and Sastry [17] incorporated an assumption that the frailty variable is measured through the variance ρ2 of gamma distribution (f(νi)=(1/ρ)e−νi/ρ), this paper also has the same assumption. Following this the frailty term is integrated out and resulted likelihood function is given below,(3) L(γ)=∏i=1n(λ0(τi)eγ′xi1+ρ2Ψ0(τi)eγ′xi)Δi(1+ρ2Ψ0(τi)eγ′xi)−1/ρ2

On applying log, equation (3) is transformed in log-likelihood function i.e. ℓ(γ) function,(4) ℓ(γ)=∑i=1n[Δi[log(λ0(τi))+γ′xi−log(1+ρ2Ψ0(τi)eγ′xi)]−(1/ρ2)log(1+ρ2Ψ0(τi)eγ′xi)]

Newtom-Raphson method is utilized to maximize the equation (4),(5) γˆ(s)=γˆ(s−1)+P(γˆ(s−1))−1Q(γˆ(s−1))

where s representing the iteration, Q(γˆ(s−1)) stands for first-order derivative and P(γˆ(s−1)) stands for the negative of second-order derivative of equation (5) in respect of γ.Q(γˆ(s−1))=∑i=1n[Δixi−{(1+Δiρ2)Ψ0(τi)xieγ′(s−1)xi}/{1+ρ2Ψ0(τi)eγ′(s−1)xi}]

(6) P(γˆ(s−1))=X′VX

In equation (6), V is the notation of diagonal matrix (n x n) and (i,i)th element of it is defined as, vii=τiηeγ′ˆ(s−1)xi(1+Δiρ2)/(1+ρ2τiηeγ′ˆ(s−1)xi).

If collinearity presents in the data set, variance matrix raises and leads the unstable estimates. To overcome collinearity this research generalizes the approach of Ozkale and Arican [16] in parametric frailty model. This approach firstly performs principal component estimation then utilizes ridge estimation on the reduced model. Smith and Marx [20] and Aguilera et al. [1] spectrally decomposed the information matrix ζ=X′VˆX and developed a reduced model. On the basis of this approach, ζ is represented as ζ=X′VˆX=UΩU′ here U=[U1...Uq] is an orthogonal matrix of order q x q and U′ζU=A′VˆA=Ω=diag(ωj) is the diagonal matrix (q x q) of eigenvalues of ζ(ω1=ωmax≥ω2≥...≥ωq=ωmin). Orthogonal linear spans with maximum variance of the columns of matrix X is considered as principal component and represented as Aj=XUj, where Uj stand for the eigenvectors of ζ which is associated with their eigenvalues ωj,(j=1,2,...,q). These eigenvalues are the variances of principal components. To overcome the problem of collinearity, principal component frailty model is developed on a reduced set of principal components.

Aguilera et al. [1] and Sirohi and Rai [19] expressed logistic regression and Cox proportional hazard regression in matrix form by using logit transformation respectively. In continuation this paper writes frailty model in terms of matrix by using logit transformation.L=Xγ=XUU′γ=Aβ

Partitioned form of the matrix A and vector β is given as A=[ArAq−r] and β=[βr′βq′−r]′ here Aq−r(r≤q) includes the principal component which are going to be deleted. Therefor U and Ω is partitioned as, U=[UrUq−r] and[Ωr00Ωq−r]

respectively, where Ar′VˆAr=Ur′X′VˆXUr=Ωr and Aq−r′VˆAq−r=Uq−r′X′VˆXUq−r=Ωq−r. On incorporating decomposition, logit transformation can be represented as, L=Aβ=Arβr+Aq−rβq−r. Collinearity is addressed by principal component regression, then the model L=Xγ=XUrUr′γ+XUq−rUq−r′γ=Arβr+Aq−rβq−r is shrinked into the model L=Arβr where βr=Ur′γ and Ar=XUr. Shrinked frailty model is written as,λ(τi|a(r)i,νi)=νiλ0(τi)exp(βr′a(r)i)

where, a(r)i′=[ai1,ai2,...,air] stands for the ith observations of Ar and βr′=[β1,β2,...,βr] is the vectors of parameters. Log likelihood function of shrinked frailty model is defined as,(7) ℓ(βr)=∑i=1n[Δi[log(λ0(τi))+βr′a(r)i−log(1+ρ2Ψ0(τi)eβr′a(r)i)]−(1/ρ2)log(1+ρ2Ψ0(τi)eβr′a(r)i)]

Newton-Raphson method is utilized to maximize equation (7),βrˆ(s)=βrˆ(s−1)+P(βrˆ(s−1))−1Q(βrˆ(s−1))

where Q(βˆ(s−1)) again stands for the first order derivative of equation (10), P(βˆ(s−1))−1 stands for the negative of second order derivative. Estimated parameters are converted into the original form as, γˆr=Urβˆr=UrUr′γˆ and named as principal component frailty (PCF) estimator.

On following the approach of Ozkale and Arican [16], penalty term is introduced into the log-likelihood function of shrinked frailty model.ℓk(βr)=ℓ(βr)−k12‖β‖r2

On maximizing penalized likelihood functionβrˆk(s)=βrˆk(s−1)+Pk(βrˆk(s−1))−1Qk(βrˆk(s−1))

where,Qk(βrˆk(s−1))=Q(βrˆk(s−1))−kβrˆk(s−1)

Pk(βrˆk(s−1))=P(βrˆk(s−1))+kIr

Ir is representing the identity matrix (r x r). On using the Cessie and Houwelingen [24] approach,βˆrk(s)=(Ar′VˆAr+kIr)−1(Ar′VˆAr)βˆr(s)

which is expressed in the original form as,γˆrk(s)=Ur(Ur′X′VˆXUr+kIr)−1Ur′X′VˆXUrUr′γˆ(s)

and named as r−k frailty estimator. If we take r→q, then r−k converts into ridge frailty estimator and if k→0, then r−k converts into principal component (PC) frailty estimator. If r→q and k→0 then r−k frailty estimator tends to Maximum likelihood (ML) frailty estimator.

3 Comparison of estimators

Asymptotic bias of r−k frailty estimator is,E(γˆrk−γ)=[Ur(Ur′X′VˆXUr+kIr)−1Ur′X′VˆXUrUr′−Iq]γ

On taking Iq=UrUr′+Uq−rUq−r′(8) E(γˆrk−γ)=[−kUr(Ur′X′VˆXUr+kIr)−1Ur′−Uq−rUq−r′]γ

Asymptotic covariance is defined as,(9) Cov(γˆrk)=Ur(Ur′X′VˆXUr+kIr)−1Ur′X′VˆXUrUr′Cov(γˆ)UrUr′X′VˆXUrUr′(Ur′X′VˆXUr+kIr)−1

On the basis of equation (8) and (9), MSE is defined as,MSE(γˆrk)=Ur(Ur′X′VˆXUr+kIr)−1Ur′X′VˆXUrUr′Cov(γˆ)UrUr′X′VˆXUrUr′(Ur′X′VˆXUr+kIr)−1+[−kUr(Ur′X′VˆXUr+kIr)−1Ur′−Uq−rUq−r′]γγ′[−kUr(Ur′X′VˆXUr+kIr)−1Ur′−Uq−rUq−r′]

we know that, (X′VˆX)−1 can be rewritten as, U(Ω)−1U′=Ur(Ωr)−1Ur′+Uq−r(Ωq−r)−1Uq−r′MSE(γˆrk)=Ur(Ωr+kIr)−1Ωr(Ωr+kIr)−1Ur′+[−kUr(Ωr+kIr)−1Ur′−Uq−rUq−r′]γγ′[−kUr(Ωr+kIr)−1Ur′−Uq−rUq−r′]

If r=q, then MSE of ridge frailty estimator is defined as,MSE(γˆk)=(X′VˆX+kIq)−1X′VˆX(X′VˆX+kIq)−1+k2(X′VˆX+kIq)−1γγ′(X′VˆX+kIq)−1

MSE(γˆk)=U(Ω+kIq)−1Ω(Ω+kIq)−1U′+k2U(Ω+kIq)−1U′γγ′U(Ω+kIq)−1U′

If k=0, then MSE of PC frailty estimator is defined as,MSE(γˆr)=UrΩr−1Ur′+Uq−rUq−r′γγ′Uq−rUq−r′

If k=0 and r=q, then MSE of ML frailty estimator is defined as,MSE(γˆ)=(X′VˆX)−1=UΩ−1U′

Three choices of k i.e. 1/γˆ′γˆ, q/γˆ′γˆ, (q+1)/γˆ′γˆ are taken on the basis of the approach of Xue et al. [23], Sirohi and Rai [19] and Albalawi et al. [3].

4 Simulation study

Estimators are compared through simulation study. Survival function is defined as,S(τi)=exp(−Ψ0(τi)eγ′xi)

Distributions function (G) assumed to follow uniform distribution i.e. U(0,1) then (1−G)∼U(0,1). Therefore,G=exp(−Ψ0(τi)eγ′xi)∼U(0,1)

τ=Ψ0−1(−log(G)/eγ′xi)

On the basis of Weibull distribution assumption with one parameter η=0.6, Ψ0−1(.)=t1/η and survival time is defined as,(10) τ=(−log(G)/eγ′xi)1/η

Equation (10) is utilized to generate the survival time. Size of samples are taken as, 50,150,250,350,450,5000. Censored time is generated through uniform distribution with in the range zero to six and gave approximate 39% censoring rate. Number of explanatory variables are taken as 3,5,7,9,11 among these, two variables are generated on the basis of the values of correlation coefficients (p)=0.60,0.80,0.85,0.90,0.95,0.99. Rest of the variables are generated through standard normal distribution. Initial values of the parameters are taken as, 0.150,0.263,0.270 for (q=3), 0.150,0.263,0.270,0.331,0.121 for (q=5), 0.150,0.263,0.270,0.331, 0.121,0.167,0.216 for (q=7), 0.150,0.263,0.270,0.331,0.121,0.167, 0.216,0.450,0.163 for q=9 and 0.150,0.263,0.270,0.331,0.121,0.167, 0.216,0.450,0.163,0.034,0.045 for q=11. Replication number of the experiment is taken as 1000 for each category. Convergence criteria of the estimated parameter are taken as, 10−4 and iteration is terminated at ‖γ(s)−γ(s−1)≤10−5‖. Averaged MSE are computed to compare the r−k, PC, ridge and ML frailty estimators. Findings of simulation are given in Table 1, Table 2, Table 3, Table 4, Table 5, Table 6.Table 1 Averaged Mean Square error (n = 50).

Table 1p		γˆ	γˆr	γˆk	γˆrk	γˆk	γˆrk	γˆk	γˆrk	
k=1/γˆ′γˆ	k=q/γˆ′γˆ	k=(q+1)/γˆ′γˆ	
0.6	q=3	0.14633	0.07671	0.02600	0.01817	0.00906	0.00647	0.00650	0.00463	
0.8		0.22965	0.07704	0.02985	0.02284	0.01119	0.00903	0.00831	0.00665	
0.85		0.28459	0.07627	0.03275	0.02384	0.01270	0.01025	0.00964	0.00781	
0.9		0.40190	0.07376	0.03705	0.02706	0.01493	0.01287	0.01155	0.01008	
0.95		0.74523	0.07901	0.04296	0.03207	0.01791	0.01632	0.01411	0.01303	
0.99		3.50290	0.07731	0.11773	0.05038	0.04229	0.03537	0.03369	0.03118	
0.6	q=5	0.23208	0.17267	0.06643	0.05289	0.01524	0.01187	0.01230	0.00955	
0.8		0.32511	0.21391	0.07460	0.06139	0.01829	0.01524	0.01494	0.01243	
0.85		0.38361	0.23469	0.06910	0.05673	0.01760	0.01480	0.01456	0.01222	
0.9		0.50007	0.26194	0.07257	0.06133	0.01824	0.01535	0.01499	0.01256	
0.95		0.86207	0.30457	0.08149	0.07072	0.02304	0.02021	0.01933	0.01687	
0.99		3.83314	0.64316	0.16453	0.13585	0.05427	0.05606	0.04820	0.04968	
0.6	q=7	0.32498	0.26469	0.12499	0.10629	0.02543	0.02159	0.02190	0.01856	
0.8		0.32511	0.21391	0.07460	0.06139	0.01829	0.01524	0.01494	0.01243	
0.85		0.48510	0.35109	0.13811	0.11950	0.02827	0.02419	0.02440	0.02082	
0.9		0.62312	0.41296	0.14015	0.12315	0.02987	0.02610	0.02595	0.02261	
0.95		1.01744	0.52072	0.17053	0.15094	0.04442	0.03910	0.03940	0.03452	
0.99		4.14629	0.91224	0.32952	0.24949	0.09358	0.08720	0.08557	0.07994	
0.6	q=9	0.44262	0.37200	0.19725	0.17430	0.02586	0.02314	0.02244	0.02007	
0.8		0.53915	0.43754	0.19416	0.17381	0.02522	0.02272	0.02190	0.01971	
0.85		0.64351	0.48788	0.19252	0.17210	0.02493	0.02240	0.02166	0.01945	
0.9		0.74140	0.53799	0.18991	0.17298	0.02510	0.02240	0.02181	0.01945	
0.95		1.16508	0.69922	0.18161	0.16925	0.02502	0.02244	0.02178	0.01951	
0.99		1.25678	1.40652	0.16999	0.16204	0.02448	0.02189	0.02129	0.01903	
0.6	q=11	0.56273	0.48590	0.29244	0.26055	0.04794	0.04294	0.04343	0.03888	
0.8		0.67444	0.57390	0.29711	0.26812	0.04638	0.04258	0.04207	0.03861	
0.85		2.97974	0.59370	0.28707	0.25988	0.04941	0.04405	0.04540	0.04031	
0.9		0.90243	0.70097	0.32626	0.29458	0.05716	0.05270	0.05217	0.04811	
0.95		1.32056	0.91030	0.30733	0.28484	0.05338	0.04916	0.04871	0.04485	
0.99		4.88150	1.53459	0.22865	0.22249	0.02692	0.02462	0.02394	0.02188	
p: Correlation coefficient, q: No. of covariates, γ: regression coefficients.

Table 2 Averaged Mean Square error (n = 150).

Table 2p		γˆ	γˆr	γˆk	γˆrk	γˆk	γˆrk	γˆk	γˆrk	
k=1/γˆ′γˆ	k=q/γˆ′γˆ	k=(q+1)/γˆ′γˆ	
0.6	q=3	0.04658	0.02347	0.01540	0.01098	0.00561	0.00533	0.00539	0.00411	
0.8		0.07363	0.02312	0.01478	0.01111	0.00551	0.00594	0.00537	0.00476	
0.85		0.09266	0.02353	0.01851	0.01325	0.00727	0.00788	0.00714	0.00651	
0.9		0.12885	0.02374	0.01747	0.01298	0.00674	0.00752	0.00638	0.00616	
0.95		0.24191	0.02375	0.02042	0.01468	0.00786	0.00945	0.00765	0.00802	
0.99		1.13137	0.02320	0.02879	0.01602	0.01089	0.01120	0.01027	0.00984	
0.6	q=5	0.06969	0.05266	0.02367	0.01901	0.00540	0.00508	0.00518	0.00416	
0.8		0.09770	0.06668	0.02668	0.02166	0.00616	0.00595	0.00603	0.00491	
0.85		0.11495	0.06927	0.02447	0.02062	0.00529	0.00517	0.00507	0.00422	
0.9		0.15275	0.08189	0.02524	0.02168	0.00600	0.00585	0.00573	0.00483	
0.95		0.26590	0.09880	0.03113	0.02714	0.00878	0.00883	0.00852	0.00752	
0.99		1.16635	0.15840	0.05196	0.04614	0.01713	0.01830	0.01640	0.01614	
0.6	q=7	0.09340	0.07793	0.04007	0.03447	0.00646	0.00690	0.00672	0.00589	
0.8		0.12208	0.09536	0.04053	0.03522	0.00692	0.00716	0.00704	0.00614	
0.85		0.14211	0.10466	0.04623	0.04034	0.00891	0.00920	0.00926	0.00800	
0.9		0.18032	0.11751	0.04661	0.04162	0.00900	0.00924	0.00904	0.00801	
0.95		0.29266	0.15590	0.04488	0.04182	0.00927	0.00944	0.00938	0.00823	
0.99		1.21029	0.33298	0.06146	0.05988	0.01950	0.01946	0.01871	0.01763	
0.6	q=9	0.11925	0.10372	0.09174	0.08064	0.02898	0.02694	0.02755	0.02456	
0.8		0.14839	0.12322	0.09868	0.08642	0.02822	0.02655	0.02695	0.02420	
0.85		0.16791	0.13546	0.10071	0.08837	0.02814	0.02643	0.02678	0.02411	
0.9		0.20558	0.15491	0.10210	0.08999	0.02766	0.02621	0.02653	0.02393	
0.95		0.32227	0.19464	0.09844	0.09024	0.02671	0.02576	0.02582	0.02353	
0.99		1.25403	0.44267	0.08287	0.08668	0.02650	0.02515	0.02513	0.02299	
0.6	q=11	0.14611	0.13075	0.11523	0.10384	0.03276	0.03141	0.03187	0.02903	
0.8		0.17556	0.15281	0.12315	0.11072	0.03175	0.03063	0.03089	0.02828	
0.85		0.19471	0.16495	0.12558	0.11277	0.03149	0.03063	0.03084	0.02831	
0.9		0.23580	0.18629	0.12816	0.11578	0.03136	0.03044	0.03056	0.02811	
0.95		0.35143	0.23044	0.12456	0.11570	0.03059	0.02981	0.02984	0.02756	
0.99		1.29604	0.52297	0.10660	0.10362	0.03023	0.02942	0.02956	0.02722	
p: Correlation coefficient, q: No. of covariates, γ: regression coefficients.

Table 3 Averaged Mean Square error (n = 250).

Table 3p		γˆ	γˆr	γˆk	γˆrk	γˆk	γˆrk	γˆk	γˆrk	
k=1/γˆ′γˆ	k=q/γˆ′γˆ	k=(q+1)/γˆ′γˆ	
0.6	q=3	0.02767	0.01343	0.01251	0.00828	0.00518	0.00517	0.00563	0.00433	
0.8		0.04367	0.01360	0.01279	0.00886	0.00489	0.00576	0.00529	0.00484	
0.85		0.05466	0.01359	0.01224	0.00876	0.00469	0.00574	0.00497	0.00485	
0.9		0.07670	0.01386	0.01262	0.00934	0.00498	0.00650	0.00527	0.00561	
0.95		0.14325	0.01363	0.01283	0.00955	0.00505	0.00690	0.00538	0.00605	
0.99		0.67649	0.01369	0.02535	0.01087	0.00838	0.00870	0.00794	0.00799	
0.6	q=5	0.04102	0.03124	0.01564	0.01263	0.00332	0.00353	0.00360	0.00290	
0.8		0.05716	0.03784	0.01669	0.01397	0.00376	0.00414	0.00410	0.00345	
0.85		0.06845	0.04087	0.01677	0.01405	0.00388	0.00415	0.00417	0.00346	
0.9		0.09039	0.04457	0.01702	0.01507	0.00406	0.00452	0.00437	0.00378	
0.95		0.15758	0.05588	0.01921	0.01817	0.00542	0.00625	0.00575	0.00532	
0.99		0.68968	0.09352	0.02970	0.02498	0.00931	0.01056	0.00921	0.00937	
0.6	q=7	0.05514	0.04618	0.02813	0.02417	0.00531	0.00605	0.00607	0.00527	
0.8		0.07165	0.05554	0.02752	0.02387	0.00505	0.00571	0.00568	0.00497	
0.85		0.07165	0.05554	0.02752	0.02387	0.00505	0.00571	0.00568	0.00497	
0.9		0.10553	0.07059	0.03064	0.02736	0.00622	0.00694	0.00699	0.00609	
0.95		0.17167	0.09437	0.03142	0.02896	0.00731	0.00780	0.00779	0.00691	
0.99		0.71156	0.17449	0.04322	0.04153	0.01410	0.01478	0.01398	0.01350	
0.6	q=9	0.06908	0.06045	0.05903	0.05197	0.02627	0.02410	0.02518	0.02248	
0.8		0.08619	0.07143	0.06633	0.05717	0.02561	0.02377	0.02464	0.02213	
0.85		0.09713	0.07878	0.06907	0.05951	0.02521	0.02349	0.02419	0.02186	
0.9		0.11942	0.08766	0.07261	0.06160	0.02461	0.02330	0.02370	0.02167	
0.95		0.18848	0.11539	0.07549	0.06576	0.02399	0.02310	0.02331	0.02151	
0.99		0.73688	0.20692	0.06098	0.06788	0.02310	0.02277	0.02259	0.02121	
0.6	q=11	0.08378	0.07525	0.07285	0.06569	0.02999	0.02833	0.02927	0.02668	
0.8		0.10084	0.08769	0.08063	0.07189	0.02956	0.02807	0.02885	0.02641	
0.85		0.10084	0.08769	0.08063	0.07189	0.02956	0.02807	0.02885	0.02641	
0.9		0.13570	0.10851	0.08841	0.07825	0.02850	0.02768	0.02808	0.02604	
0.95		0.20394	0.13637	0.09231	0.08169	0.02796	0.02725	0.02763	0.02566	
0.99		0.75606	0.25660	0.07763	0.07236	0.02715	0.02682	0.02677	0.02527	
p: Correlation coefficient, q: No. of covariates, γ: regression coefficients.

Table 4 Averaged Mean Square error (n = 350).

Table 4p		γˆ	γˆr	γˆk	γˆrk	γˆk	γˆrk	γˆk	γˆrk	
k=1/γˆ′γˆ	k=q/γˆ′γˆ	k=(q+1)/γˆ′γˆ	
0.6	q=3	0.01969	0.00973	0.01075	0.00690	0.00494	0.00498	0.00563	0.00438	
0.8		0.03121	0.00969	0.01046	0.00687	0.00424	0.00506	0.00473	0.00447	
0.85		0.03892	0.00968	0.01042	0.00699	0.00409	0.00528	0.00463	0.00468	
0.9		0.05474	0.00984	0.01112	0.00741	0.00433	0.00587	0.00485	0.00529	
0.95		0.10205	0.00968	0.01060	0.00741	0.00423	0.00599	0.00469	0.00544	
0.99		0.48025	0.00960	0.01736	0.00793	0.00594	0.00675	0.00601	0.00632	
0.6	q=5	0.02933	0.02239	0.01367	0.01102	0.00339	0.00371	0.00387	0.00313	
0.8		0.04056	0.02661	0.01301	0.01094	0.00302	0.00347	0.00348	0.00291	
0.85		0.04842	0.02916	0.01401	0.01198	0.00355	0.00407	0.00407	0.00345	
0.9		0.06392	0.03471	0.01365	0.01202	0.00346	0.00395	0.00395	0.00335	
0.95		0.11198	0.03955	0.01549	0.01453	0.00447	0.00519	0.00495	0.00446	
0.99		0.48918	0.07351	0.02276	0.01930	0.00744	0.00873	0.00748	0.00780	
0.6	q=7	0.03926	0.03283	0.02180	0.01868	0.00468	0.00549	0.00561	0.00485	
0.8		0.05066	0.03915	0.02174	0.01903	0.00418	0.00498	0.00497	0.00437	
0.85		0.05868	0.04333	0.02246	0.01964	0.00460	0.00543	0.00552	0.00480	
0.9		0.07441	0.04960	0.02258	0.02014	0.00476	0.00565	0.00562	0.00501	
0.95		0.12183	0.06420	0.02252	0.02148	0.00523	0.00607	0.00600	0.00539	
0.99		0.50435	0.13682	0.02670	0.02821	0.00876	0.00932	0.00901	0.00847	
0.6	q=9	0.04914	0.04305	0.04389	0.03862	0.02334	0.02120	0.02245	0.02003	
0.8		0.06045	0.05049	0.04985	0.04287	0.02279	0.02098	0.02193	0.01976	
0.85		0.06823	0.05557	0.05279	0.04513	0.02243	0.02082	0.02165	0.01960	
0.9		0.08443	0.06186	0.05724	0.04747	0.02197	0.02084	0.02130	0.01961	
0.95		0.13225	0.08329	0.06197	0.05214	0.02118	0.02025	0.02040	0.01904	
0.99		0.51562	0.14699	0.05256	0.05356	0.01996	0.01998	0.01948	0.01880	
0.6	q=11	0.05899	0.05299	0.05343	0.04814	0.02690	0.02505	0.02621	0.02385	
0.8		0.07060	0.06136	0.05994	0.05316	0.02626	0.02480	0.02565	0.02357	
0.85		0.07853	0.06716	0.06328	0.05591	0.02600	0.02470	0.02550	0.02347	
0.9		0.09489	0.07790	0.06833	0.05999	0.02549	0.02421	0.02487	0.02297	
0.95		0.14257	0.09859	0.07413	0.06425	0.02476	0.02388	0.02427	0.02269	
0.99		0.53127	0.20607	0.06491	0.06370	0.02400	0.02354	0.02352	0.02238	
p: Correlation coefficient, q: No. of covariates, γ: regression coefficients.

Table 5 Averaged Mean Square error (n = 450).

Table 5p		γˆ	γˆr	γˆk	γˆrk	γˆk	γˆrk	γˆk	γˆrk	
k=1/γˆ′γˆ	k=q/γˆ′γˆ	k=(q+1)/γˆ′γˆ	
0.6	q=3	0.01517	0.00752	0.00888	0.00564	0.00435	0.00440	0.00503	0.00399	
0.8		0.02423	0.00747	0.00961	0.00590	0.00418	0.00505	0.00472	0.00471	
0.85		0.03026	0.00757	0.00960	0.00599	0.00395	0.00513	0.00449	0.00477	
0.9		0.04258	0.00764	0.00991	0.00607	0.00394	0.00522	0.00437	0.00487	
0.95		0.07937	0.00761	0.00906	0.00616	0.00365	0.00544	0.00409	0.00512	
0.99		0.37265	0.00749	0.01468	0.00646	0.00510	0.00585	0.00516	0.00562	
0.6	q=5	0.02250	0.01739	0.01127	0.00908	0.00292	0.00332	0.00349	0.00283	
0.8		0.03138	0.02081	0.01094	0.00919	0.00266	0.00315	0.00315	0.00267	
0.85		0.03773	0.02346	0.01241	0.01058	0.00329	0.00394	0.00396	0.00339	
0.9		0.04968	0.02556	0.01184	0.01059	0.00312	0.00382	0.00368	0.00327	
0.95		0.08649	0.03264	0.01330	0.01217	0.00393	0.00477	0.00441	0.00413	
0.99		0.38060	0.06167	0.01673	0.01599	0.00579	0.00715	0.00617	0.00639	
0.6	q=7	0.03012	0.02517	0.01691	0.01451	0.00353	0.00427	0.00432	0.00377	
0.8		0.03909	0.03080	0.01790	0.01551	0.00377	0.00459	0.00464	0.00407	
0.85		0.04519	0.03366	0.01825	0.01598	0.00387	0.00463	0.00471	0.00411	
0.9		0.05754	0.03893	0.01846	0.01637	0.00406	0.00484	0.00492	0.00431	
0.95		0.09440	0.05021	0.01968	0.01844	0.00496	0.00584	0.00585	0.00526	
0.99		0.39060	0.09945	0.02743	0.02679	0.00841	0.00949	0.00894	0.00872	
0.6	q=9	0.03767	0.03303	0.03453	0.03037	0.02063	0.01866	0.01998	0.01780	
0.8		0.04669	0.03906	0.04003	0.03429	0.02047	0.01871	0.01971	0.01776	
0.85		0.05272	0.04303	0.04289	0.03642	0.02016	0.01852	0.01940	0.01754	
0.9		0.06505	0.04990	0.04732	0.03944	0.01965	0.01831	0.01900	0.01732	
0.95		0.10210	0.05866	0.05368	0.04184	0.01882	0.01835	0.01821	0.01736	
0.99		0.39855	0.12554	0.04793	0.04636	0.01769	0.01776	0.01715	0.01681	
0.6	q=11	0.04535	0.04087	0.04201	0.03793	0.02392	0.02213	0.02331	0.02123	
0.8		0.05451	0.04741	0.04790	0.04234	0.02375	0.02227	0.02319	0.02131	
0.85		0.06050	0.05192	0.05090	0.04485	0.02347	0.02199	0.02286	0.02103	
0.9		0.07283	0.05805	0.05571	0.04759	0.02299	0.02195	0.02251	0.02098	
0.95		0.11005	0.07355	0.06326	0.05230	0.02223	0.02157	0.02178	0.02061	
0.99		0.40885	0.14253	0.05949	0.05771	0.02119	0.02107	0.02082	0.02014	
p: Correlation coefficient, q: No. of covariates, γ: regression coefficients.

Table 6 Averaged Mean Square error (n = 5000).

Table 6p		γˆ	γˆr	γˆk	γˆrk	γˆk	γˆrk	γˆk	γˆrk	
k=1/γˆ′γˆ	k=q/γˆ′γˆ	k=(q+1)/γˆ′γˆ	
0.6	q=3	0.001365	0.000664	0.001284	0.000648	0.001631	0.000672	0.001213	0.000704	
0.8		0.002167	0.000666	0.001896	0.000653	0.001914	0.000700	0.001550	0.000749	
0.85		0.002713	0.000667	0.002260	0.000655	0.002039	0.000708	0.001696	0.000761	
0.9		0.003806	0.000670	0.002863	0.000660	0.002186	0.000723	0.001813	0.000783	
0.95		0.007097	0.000667	0.003989	0.000658	0.002217	0.000729	0.001744	0.000794	
0.99		0.033488	0.000668	0.005104	0.000660	0.001648	0.000728	0.001165	0.000792	
0.6	q=5	0.002004	0.001505	0.001798	0.001384	0.001405	0.001270	0.001471	0.001270	
0.8		0.002797	0.001828	0.002255	0.001633	0.001466	0.001536	0.001497	0.001538	
0.85		0.003339	0.002078	0.002506	0.001788	0.001482	0.001613	0.001503	0.001607	
0.9		0.004424	0.002370	0.002831	0.001973	0.001421	0.001772	0.001442	0.001750	
0.95		0.007693	0.002688	0.003268	0.002222	0.001327	0.002107	0.001334	0.002066	
0.99		0.033815	0.004412	0.003763	0.002919	0.001269	0.002523	0.001224	0.002466	
0.6	q=7	0.001366	0.000663	0.001285	0.000648	0.001636	0.000674	0.001207	0.000708	
0.8		0.003466	0.002719	0.002961	0.002418	0.001816	0.001901	0.001997	0.001876	
0.85		0.004002	0.002985	0.003212	0.002580	0.001785	0.001976	0.001960	0.001942	
0.9		0.005088	0.003449	0.003664	0.002878	0.001795	0.002084	0.001954	0.002035	
0.95		0.008358	0.004516	0.004320	0.003372	0.001704	0.002277	0.001847	0.002208	
0.99		0.034542	0.006974	0.004916	0.004268	0.001677	0.003032	0.001731	0.002854	
0.6	q=9	0.003311	0.002911	0.003285	0.002889	0.003244	0.002731	0.003077	0.002714	
0.8		0.004101	0.003411	0.004040	0.003370	0.003787	0.003090	0.003586	0.003060	
0.85		0.004640	0.003764	0.004544	0.003703	0.004092	0.003310	0.003866	0.003270	
0.9		0.005735	0.004354	0.005537	0.004249	0.004571	0.003633	0.004311	0.003575	
0.95		0.008994	0.005485	0.008278	0.005243	0.005312	0.004131	0.005004	0.004046	
0.99		0.035112	0.011353	0.022172	0.009318	0.004923	0.005363	0.004694	0.005174	
0.6	q=11	0.003966	0.003578	0.003939	0.003554	0.003883	0.003344	0.003677	0.003325	
0.8		0.004758	0.004146	0.004699	0.004102	0.004406	0.003738	0.004172	0.003707	
0.85		0.005307	0.004526	0.005216	0.004463	0.004700	0.003964	0.004447	0.003924	
0.9		0.006385	0.005111	0.006203	0.005007	0.005128	0.004267	0.004848	0.004213	
0.95		0.009655	0.006703	0.009002	0.006420	0.005776	0.004878	0.005458	0.004788	
0.99		0.035771	0.013124	0.023745	0.010864	0.005249	0.006013	0.005067	0.005847	
p: Correlation coefficient, q: No. of covariates, γ: regression coefficients.

Result of simulation shows that the r−k frailty estimators have smaller values of MSE as compared to ML, PC and ridge frailty estimator for all categories of p,n,q and k except n=5000 and these MSE are minimum for k=(q+1)/γˆ′γˆ. But ridge frailty estimator for k=(q+1)/γˆ′γˆ,q/γˆ′γˆ, big sample sizes (n=5000), high correlation (p=0.99) and (q≥5) has smallest MSE as compared to ML, PC and r−k frailty estimators.

5 Application

To show the importance of proposed estimator in real life, this research utilized data from the third phase of NFHS (National Family Health survey) India. But this data is confidential and not allowed to print. Data includes 11,581 infants out of which 10,448 are right censored. Here, dependent variable is taken as the survival time of infant in months and eight explanatory variables are considered namely, Family Size (FS), Iron Folic Acid (IFA), Mother Education (ME), Mother Age (MA), Tetanus Toxoid (TT), Antenatal Care (ANC), Marital Status (MS), Breastfeeding (BF). Out of these TT and ANC explain each other upto 60% i.e. p=0.6.

Sometime correlation coefficient misleads collinearity therefore, this research utilizes the Ozkale [15]'s approach to compute condition number. The calculated condition number is 36.45 and it is greater than 10. This is the indication of collinearity [7], [22].

Frailty model with different estimation techniques is applied on the data set. Table 7, compares all the estimators in terms of MSE and revealed that the r−k frailty estimator has smaller MSE as compared to the other estimators and the value is minimum for k=(q+1)/γˆ′γˆ. Therefore, Table 8, shows the estimated values of the parameters corresponding to different covariates at k=(q+1)/γˆ′γˆ. If hazard ratio (exp(γˆ)) is less than one, there is a decreased risk of death otherwise there is an increased risk of death.Table 7 Mean square Error of estimators in NFHS data.

Table 7Method	MSE(γ)	
MLE	4.214	
γˆk	3.124	
γˆk (k=1/γˆ′γˆ)	3.025	
γˆrk (k=1/γˆ′γˆ)	2.548	
γˆk (k=q/γˆ′γˆ)	3.012	
γˆrk (k=q/γˆ′γˆ)	2.037	
γˆk (k=(q+1)/γˆ′γˆ)	2.127	
γˆrk (k=(q+1)/γˆ′γˆ)	1.982	

Table 8 Estimated parameters in NFHS Data.

Table 8Category	exp(γˆ)	exp(γˆr)	exp(γˆk)	exp(γrˆk)	
FS					
One to Three	2.342⁎	2.331⁎	2.519⁎	2.723⁎	
Four to Six#	-	-	-	-	
Greater than Six	0.612⁎	0.659⁎	0.834⁎	0.921⁎	
IFA					
No#	-	-	-	-	
Yes	1.153⁎	1.112⁎	1.261⁎	1.312⁎	
ME					
Not Educated	10.631⁎	10.695⁎	10.732⁎	10.921⁎	
Primary	10.230⁎	10.123⁎	10.342⁎	10.639⁎	
Secondary	7.461⁎	7.731⁎	7.923⁎	8.036⁎	
Higher#	-	-	-	-	
MA (in Years)					
< 18	1.356⁎	1.389⁎	1.493⁎	1.562⁎	
18 - 24#	-	-	-	-	
25 - 29	1.412⁎	1.459⁎	1.515⁎	1.623⁎	
30 - 34	1.631⁎	1.672⁎	1.722⁎	1.917⁎	
Greater than 35	1.245⁎	1.213⁎	1.316⁎	1.441⁎	
TT					
No	1.121	1.113	1.261	1.273	
One	1.012⁎	1.013⁎	1.142⁎	1.241⁎	
Greater than One#	-	-	-	-	
ANC					
No	1.216	1.312	1.161	1.215	
One to Two	1.362⁎	1.316⁎	1.318⁎	1.531⁎	
Greater than Two#	-	-	-	-	
Marital status					
unmarried	0.413	0.512	0.612	0.712	
0 - 4 Years#	-	-	-	-	
5 - 9 Years	1.312⁎	1.415⁎	1.561⁎	1.712⁎	
>=10 Years	1.726⁎	1.816⁎	1.912⁎	2.112⁎	
Breastfeeding					
No#	-	-	-	-	
Yes	7.881E-05⁎	7.761E-05⁎	7.912E-05⁎	7.325E-05⁎	
Statistical significance: ⁎ at 5% level of significance.

Reference Category: #.

All the variables have significant impact on infant survival and interpretation is given as: If size of the family is small (one to three members in the family), the risk of infant death is 2.3 times higher as compared to the infant who is having the family of 4−6 members. If infant is breastfeeding, risk of infant death will be 99% lower as compared to the infant who is not breastfeeding. Uneducated mother has 9.9 times higher risk of infant death as compared to highly educated mother. In the same manner we can see the impact of other variables. Estimated value of ρ2 is 2.85 which is non zero that means unobserved variables also have impact on infant death.

6 Discussion and conclusion

Frailty model is a mixture model in survival analysis which considers time to event data with observable and un-observable variables simultaneously. In case of collinearity estimated parameters of frailty becomes unstable. Therefore this research proposed the biased (r−k) frailty estimator and compared it with the ML, PC and ridge frailty estimator through simulation study. Findings suggest that the proposed estimator performs better as compared to other estimators for all categories of p,q and k at the sample sizes 50, 150, 250, 350 and 450. Precision of the estimator is more accurate at k=(q+1)/γˆ′γˆ.

Further, the developed approach is applied to examine the infant mortality in NFHS data, India where, collinearity was present in the data set. Here again it is found that the r−k frailty estimator performed better than the other estimators and results are more benefited at k=(q+1)/γˆ′γˆ. Frailty model suggests that breastfeeding has a very positive significant impact on infant mortality. To make less infant mortality, people should believe on joint family system, avoid teenage pregnancy and promote women education. Meanwhile mother should take proper antenatal checkup and IFA tablets.

This research considers the parametric frailty model, in future semi-parametric model can also be considered. Researchers can also work on the estimation of ridge parameter.

CRediT authorship contribution statement

Marwan H. Ahelali: Formal analysis, Conceptualization. Osama Abdulaziz Alamri: Conceptualization, Formal analysis. Anu Sirohi: Methodology, Software.

Declaration of Competing Interest

The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.
==== Refs
References

1 Aguilera A.M. Escabias M. Valderrama M.J. Using principal components for estimating logistic regression with high-dimensional multicollinear data Comput. Stat. Data Anal. 50 8 2006 1905 1924
2 Ahmad S. Aslam M. Ahmad S. Extending the Liu estimator for the Cox proportional hazards regression model with multicollinearity Commun. Stat., Simul. Comput. 2023 1 14
3 Albalawi O. Sirohi A. Rai P.K. Alanzi A.R.A. Parametric frailty analysis in presence of collinearity: an application to assessment of infant mortality Mathematics 10 13 2022 1 10
4 Cox D.R. Regression models and life tables J. R. Stat. Soc. B 34 1972 187 202
5 Hoerl A.E. Kennard R.W. Ridge regression: biased estimation for nonorthogonal problems Technometrics 12 1 1970 55 67
6 Hoerl A.E. Kannard R.W. Baldwin K.F. Ridge regression: some simulations Commun. Stat., Theory Methods 4 2 1975 105 123
7 Mackinnon M.J. Puterman A.L. Collinearity in generalized linear models Commun. Stat., Theory Methods 18 9 1989 3463 3472
8 Mahdizadeh M. Zamanzade E. On estimating the area under the ROC curve in ranked set sampling Stat. Methods Med. Res. 31 8 2022 1500 1514 35549545
9 Mahdizadeh M. Zamanzade E. Smooth estimation of the area under the ROC curve in multistage ranked set sampling Stat. Pap. 62 4 2021 1753 1776
10 Zamanzade E. Mahdizadeh M. Samawi H.M. Efficient estimation of cumulative distribution function using moving extreme ranked set sampling with application to reliability AStA Adv. Stat. Anal. 104 3 2020 485 502
11 Mahdizadeh M. Zamanzade E. Dynamic reliability estimation in a rank-based design Probab. Math. Stat. 39 2019 1 18
12 Mahdizadeh M. Zamanzade E. Reliability estimation in multistage ranked set sampling REVSTAT Stat. J. 15 4 2017 565 581
13 Michimae H. Emura T. Bayesian ridge regression for survival data based on a vine copula-based prior AStA Adv. Stat. Anal. 107 4 2023 755 784
14 Oakes D. A model for association in bivariate survival data J. R. Stat. Soc., Ser. B, Stat. Methodol. 44 3 1982 414 422
15 Ozkale M.R. The red indicator and corrected VIFs in generalized linear models Commun. Stat., Simul. Comput. 50 12 2021 4144 4170
16 Ozkale M.R. Arican E. A new biased estimator in logistic regression model Statistics 50 2 2016 233 253
17 Sastry N. A nested frailty model for survival data, with an application to the study of child survival in Northeast Brazil J. Am. Stat. Assoc. 92 438 1997 426 435 12155403
18 Schafer R.L. Roi L.D. Wolfe R.A. A ridge logistic estimator Commun. Stat., Theory Methods 13 1 1984 99 113
19 Sirohi A. Rai P.K. Some r−k class proportional hazard regression models in the presence of collinearity: an evidence from Indian infant mortality Commun. Stat., Simul. Comput. 52 10 2023 4880 4890
20 Smith E.P. Marx B.D. III-conditioned information matrices, generalized linear models and estimation of the effects of acid rain Environmetrics 1 1 1990 57 71
21 Vaupel J. Manton K. Stallard E. The impact of heterogeneity in individual frailty on the dynamics of mortality Demography 16 3 1979 439 454 510638
22 Weissfeld L.A. Sereika S.M. A multicollinearity diagnostic for generalized linear models Commun. Stat., Theory Methods 20 4 1991 1183 1198
23 Xue X. Kim M.Y. Roy E.S. Cox regression analysis in presence of collinearity: an application to assessment of health risks associated with occupational radiation exposure Lifetime Data Anal. 13 2007 333 350 17661175
24 Cessie S.L. Houwelingen J.C.V. Ridge estimators in logistic regression Applied Statistics 41 1992 191 201
