This paper is devoted to the technique of construction of confidence ellipsoids for coefficients of linear regression in the case, when statistical data are derived from a mixture with finite number of components and the mixing probabilities (the concentrations of components) are different for different observations. These mixing probabilities are assumed to be known, but the distributions of the components are unknown. (Such mixture models are also known as mixtures with varying concentrations, see [1] and [7]).

The problem of estimation of regression coefficients by such mixed observations was considered in the parametric setting in [4] and [5]. The authors of these papers assume that the distributions of regression errors and regressors are known up to some unknown parameters. The models considered in these papers are called finite mixtures of regression models. Some versions of maximum likelihood estimates are used for the estimation of unknown parameters of distributions and regression coefficients under these models. The EM-algorithm is used for computation of the estimates. (This algorithm is also implemented in R package mixtools, [2]. See [13] for the general theory of EM-algorithm and its application to mixture models).

In [6] a nonparametric approach was proposed under which no parametric models on the distributions of components are assumed. A weighted least squares technique is used to derive estimates for regression coefficients. Consistency and asymptotic normality of the estimates are demonstrated.

Note that in [6, 10, 11] a nonparametric approach to the analysis of mixtures with varying concentrations was developed in the case when the concentrations of components (mixing probabilities) are known. Some examples of real life data analysis under this assumption were considered in [10, 11].

Namely, in [11] an application to the analysis of the Ukrainian parliamentary elections-2006 was considered. Here the observed subjects were respondents of the Four Wave World Values survey (conduced in Ukraine in 2006) and the mixture components were the populations of different electoral behavior adherents. The mixing probabilities were obtained from the official results of voting in 27 regions of Ukraine.

In [10] an application to DNA microarray data analysis was presented. Here the subjects were nearly 3000 of genes of the human genome. They were divided into two components by the difference of their expression in two types of malignantly transformed tissues. The concentrations were defined as a posteriori probabilities for the genes to belong to a given component, calculated by observations on the genes expression in sample tissues.

In this paper we will show how to construct confidence sets (ellipsoids) for regression coefficients under both parametric and nonparametric approaches. Quality of obtained ellipsoids is compared via simulations.

The rest of the paper in organized as follows. In Section 2 we present a formal description of the model. Nonparametric and parametric estimates of regression coefficients and their asymptotic properties are discussed in Sections 3 and 4. Estimation of asymptotic covariances of these estimates is considered in Section 5. The confidence ellipsoids are constructed in Section 6. Results of simulations are presented in Section 7. In Section 8 we present a toy example of application to a real life sociological data. Section 9 contains concluding remarks.

We consider the model of mixture with varying concentrations. It means that each observed subject O belongs to one of M different subpopulations (mixture components). The number of component which the subject belongs to is denoted by κ(O)∈{1,2,…,M} . This characteristic of the subject is not observed. The vector of observed variables of O will be denoted by ξ(O) . It is considered as a random vector with the distribution depending on the subpopulation which O belongs to. A structural linear regression model will be used to describe these distributions. (See [14] for general theory of linear regression).

That is, we consider one variable Y=Y(O) in ξ(O)=(Y(O),X1(O),…,Xd(O))T as a response and all other ones X(O)=(X1(O),…,Xd(O))T as regressors in the model (1) Y(O)= ∑i=1dbiκ(O)Xi(O)+ε(O), where bik , i=1,…,d , k=1,…,M are unknown regression coefficients for the k-th component of the mixture, ε(O) is the error term. Denote by bk=(b1k,…,bdk)T the vector of the k-th component’s coefficients. We consider ε(O) as a random variable and assume that E[ε(O)|κ(O)=m]=0,m=1,…,M, and σm2=Var[ε(O)|κ(O)=m]<∞. ( σm2 are unknown).

It is also assumed that the regression error term ε(O) and regressors X(O) are conditionally independent for fixed κ(O)=m , m=1,…,M .

The observed sample Ξn consists of values ξj=(Yj,XjT)T=ξ(Oj) , j=1,…,n , where O1 ,…, On are independent subjects which can belong to different components with probabilities pjm=P{κ(Oj)=m},m=1,…,M;j=1,…,n. (all mixing probabilities pjm are known).

To describe completely the probabilistic behavior of the observed data we need to introduce the distributions of ε(O) and X(O) for different components. Let us denote FX,m(A)=P{X(O)∈A|κ(O)=m}for any measurableA⊆Rd, and Fε,m(A)=P{ε(O)∈A|κ(O)=m}for any measurableA⊆R. The corresponding probability densities fX,m and fε,m(x) are defined by FX,m(A)=∫AfX,m(x)dx,Fε,m(A)=∫Afε,m(x)dx (for all measurable A).

The distribution of observed ξj is a mixture of distributions of components with the mixing probabilities pjm , e.g. P{Xj∈A}= ∑m=1MpjmFXj,m(A) and the probability density fj(y,x) of ξj=(Yj,XjT)T at a point (y,xT)T∈Rd+1 is fj(y,x)= ∑m=1MpjmfX,m(x)fε,m(y−xTbm). In what follows we will discuss two approaches to the estimation of the parameters of interest bk , for a fixed k∈{1,…,M} .

The first one is the nonparametric approach. Under this approach we do not need to know the densities fX,m and fε,m . Moreover we even do not assume the existence of these densities. The estimates are based on some modification of the least squares tehnique proposed in [6].

In the second, parametric approach we assume that the densities of components are known up to some unknown nuisance parameters ϑm∈Θ⊆RL : (2) fX,m(x)=f(x;ϑm),fε,m(x)=fε(x;ϑm). In the most popular parametric normal mixture model these densities are normal, i.e. (3) fε,m∼N(0,σm2),fX,m(x)∼N(μm,Σm), where μm∈Rd is the mean of X for the m-th component and Σm∈Rd×d is its covariance matrix. All the parameters are usually unknown. So, in this case the unknown nusance parameters are ϑm=(μm,Σm,σm2),m=1,…,M.

Let us consider the nonparametric approach to the estimation of the regression coefficients developed in [6]. It is based on the minimization of weighted least squares Jk;n(b)=def1n ∑j=1naj;nk(Yj− ∑i=1dbiXji)2, over all possible b=(b1,…,bd)T∈Rd .

Here ak=(a1;nk,…,an;nk) are the minimax weights for estimation of the k-th component’s distribution. They are defined by (4) aj;nk=1detΓn ∑m=1M(−1)k+mγmk;npjm, where γmk;n is the (mk) -th minor of the matrix Γn=(1n ∑j=1npjlpji)l,i=1M, see [6, 10] for details.

Define X=def(Xji)j=1,…,n;i=1,…,d to be the n×d -matrix of observed regressors, Y=def(Y1,…,YN)T be the vector of observed responses, A=defdiag(a1;nk,…,an;nk) be the diagonal weights matrix for estimation of k-th component. Then the stationarity condition ∂Jk;n(b)∂b=0 has the unique solution in b, (5) bˆLS(k,n)=def(XTAX)−1XTAY, if the matrix XTAX is nonsingular.

Note that the weight vector ak defined by (4) contains negative weights, so bˆLS(k,n) is not allways the point of minimum of Jk;n(b) . But in what follows we will consider bˆLS(k,n) as a generalized least squares estimate for bk .

The asymptotic behavior of bˆLS(k,n) as n→∞ was investigated in [6]. To describe it we will need some additional notation.

Let us denote by D(m)=defE[X(O)XT(O)|κ(O)=m] the matrix of second moments of regressors for the m-th component.

The consistency conditions for the estimator bˆLS(k,n) are given by the following theorem.

Let us introduce the following notation. (6) Dis(m)=defE[Xi(O)Xs(O)|κ(O)=m],Lis(m)=def(E[Xi(O)Xs(O)Xq(O)Xl(O)|κ(O)=m])l,q=1d,Mis(m,p)=def(Dil(m)Dsq(p))l,q=1d,αs,q(k)=limn→∞1n ∑j=1n(aj;nk)2pjspjq (if this limit exists), αs(k)=limn→∞1n ∑j=1n(aj;nk)2pjs= ∑q=1Mαs,q(k). The following theorem provides conditions for the asymptotic normality and describes the dispersion matrix of the estimator bˆLS(k,n) . Theorem 2

Assume that 1.

E[(Xi(O))4|κ(O)=m]<∞ and E[(ε(O))4|κ(O)=m]<∞ for all m=1,…,M , i=1,…,d .

Then n(bˆLS(k,n)−b(k))⟶WN(0,V) , where (7) V=defD−1ΣD−1, (8) Σ=(Σil)il=1d,Σil= ∑s=1Mαsk(Dil(s)σs2+(bs−bk)TLil(s)(bs−bk))− ∑s=1M ∑m=1Mαs,mk(bs−bk)TMil(s,m)(bm−bk).

In this section we discuss the parametric approach to the estimation of bk based on papers [4, 5]. We will assume that the representation (2) holds with some unknown ϑm∈Θ⊆RL , m=1,…,M .

Then the set of all unknown parameters τ=(bk,β,ϑ) consists of β=(bm,m=1,…,M,m≠k) and ϑ=(ϑ1,…,ϑM). Here bk is our parameter of interest, β and ϑ are the nuisance parameters.

In this model the log-likelihood for the unknown τ by the sample Ξn can be defined as L(τ)= ∑j=1nL(ξj,pj,τ), where pj=(pj1,…,pjM)T , L(ξj,pj,τ)=ln( ∑m=1MpjmfX,m(Xj;ϑm)fε,m(Yj−XjTbm)). The general maximum likelihood estimator (MLE) τˆnMLE=(bˆk,MLE,βˆMLE,ϑˆMLE) for τ is defined as τˆnMLE=argmaxτL(τ), where the maximum is taken over all possible values of τ. Unfortunately, this estimator is not applicable to most common parametric mixture models, since the log-likelihood L(τ) usually is not bounded on the set of all possible τ.

For example, it is so in the normal mixture model (3). Really, in this model L(τ)→∞ as σ12→0 and Y1−X1Tb1=0 with all other parameters being arbitrary fixed.

The usual way to cope with this problem is to use the one-step MLE, which can be considered as one iteration of the Newton–Raphson algorithm of approximate calculation of MLE, starting from some pilot estimate (see [15], section 4.5.3). Namely, let τˆn(0) be some pilot estimate for τ. Let us consider τ as a vector of dimension P=d×M+M×L and denote its entries by τi : τ=(τ1,…,τP). Denote the gradient of L(τ) by sn(τ)=∂L(τ)∂τ=(∂L(τ)∂τ1,…,∂L(τ)∂τP)T and the Hessian of L(τ) by γn(τ)=(∂L(τ)∂τiτl)i,l=1P. Then the one-step estimator for τ starting from τ(0)ˆ is defined as τˆnOS=τˆ(0)−(γn(τˆ(0)))−1sn(τˆ(0)). Theorem 4.19 in [15] provides general conditions under which τˆnOS constructed by an i.i.d. sample is consistent, asymptotically normal and asymptotically efficient.1¹

Note that in our setting the sample Ξn is not an i.i.d. sample. But one can consider it as i.i.d. if the vectors of concentrations (pj1,…,pjM) are generated by some stochastic mechanism as i.i.d. vectors. See [12] for an example of such stochastic concentrations models.

So, if the assumptions of theorem 4.19 (or other analogous statement) hold, there is no need to use an iterative procedure to derive an estimate with asymptotically optimal performance. But on samples of moderate size τˆnOS can be not good enough.

Another popular way to obtain a stable estimate for τ is to use some version of EM-algorithm. A general EM-algorithm is an iterative procedure for approximate calculation of maximum likelihood estimates when information on some variables is missed. We describe here only the algorithm which calculates EM estimates τˆnEM under the normal mixture model assumptions (3), cf. [4, 5].

The algorithm starts from some pilot estimate τˆ(0)=(bˆm(0),σˆm2(0),μˆm(0),Σˆm(0),m=1,…,M) for the full set of the model parameters.

Then for i=1,2,… the estimates are iteratively recalculated in the following way.

Assume that on the i-iteration estimates bˆ(i) , σˆm2(i) , μˆm(i) , Σˆm(i) , m=1,…,M are obtained. Then the i-th stage weights are defined as (9) wjm(i)=wjm(ξj,τˆ(i))=pjmfX,m(Xj;μˆm(i),Σˆm(i))fε,m(Yj−XjTbˆm(i);σˆm2(i))∑l=1MpjlfX,l(Xj;μˆl(i),Σˆl(i))fε,l(Yj−XjTbˆl(i);σˆl2(i)) (note that wjm(i) is the posterior probability P{κj=m|ξj} calculated for τ=τˆ(i) ).

Let w¯m=∑j=1nwjm(i) . Then the estimators of the i+1 iteration are defined as μˆm(i+1)=1w¯m ∑j=1nwjm(i)Xj,Σˆm(i+1)=1w¯m ∑j=1nwjm(i)(Xj−μˆm(i))(Xj−μˆm(i))T,bˆm(i+1)=( ∑j=1nwjm(i)XjXjT)−1 ∑j=1nwjm(i)YjXj,σˆm2(i+1)=1w¯m ∑j=1nwjm(i)(Yj−XjTbˆm(i))2. The iterations are stopped when some stopping condition is fulfilled. For example, it can be ‖τˆ(i+1)−τˆ(i)‖<δ, where δ is a prescribed target accuracy.

It is known that this procedure provide stable estimates which (for sample large enough) converge to the point of local minimum of L(τ) which is the closest to the pilot estimator τˆ(0) .

So, this estimator can be considered as an approximate version of a root of likelihood equation estimator (RLE).

The asymptotic behavior of τˆnOS and τˆnEM can be described in terms of Fisher’s information matrix I∗(n,τ)=(Iil∗(n,τ))i,l=1P , where Iil∗(n,τ)= ∑j=1nIil(pj,τ),Iil(p,τ)=E∂L(ξp,p,τ)∂τi∂L(ξp,p,τ)∂τl, where p=(p1,…,pm) , ξp is a random vector with the pdf (10) fp(y,x;τ)= ∑m=1Mpmf(x;ϑm)fε(y−xTbm;ϑm).

Under the regularity assumptions (RR) of theorem 70.5 in [3], (11) (I∗(n,τ))1/2(τˆnMLE−τ)⟶WN(0,E), where E is the R×R unit matrix.

Assumptions (RR) include the assumption of likelihood boundedness, so they do not hold for the normal mixture model. But if the pilot estimate τˆn(0) is n -consistent, one needs only a local version of (RR) to derive asymptotic normality of τˆnOS and τˆnEM , i.e. (RR) must hold in some neighborhood of the true value of estimated τ. These local (RR) hold for the normal mixture model if σm2>0 0$]]> and Σm are nonsingular for all m=1,…,M .

To use all these results for construction of an estimator for bk we need n -consistent pilot estimators for the parameter of interest and nuisance parameters. They can be derived by the nonparametric technique considered in Section 3. To construct confidence ellipsoids we will also need estimators for the dispersion matrix V from (7) in the nonparametric case and estimators for the information matrix I∗(n,τ) in the parametric case. These estimators are discussed in the next section.

Let us start with the estimation of the dispersion matrix V in Theorem 2. In fact, we need to estimate consistently the matrices D and Σ .

Note that D=D(k)=(Dis(k))i,s=1d , where Dis(k)=E[Xi(O)Xs(O)|κ(O)=k]. By theorem 4.2 in [10], (12) Dˆnis(k)=1n ∑j=1najkXjiXjs is a consistent estimate for Dis(k) if E[‖X(O)‖2|κ(O)=m]<∞ for all m=1,…,M and assumption 3 of Theorem 1 holds.

So one can use Dˆn(k)=(Dˆnis(k))i,s=1d as a consistent estimate for D if the assumptions of Theorem 1 hold.

Similarly, Lis(m) can be estimated consistently by (13) Lˆnis(m)=1n ∑j=1najmXjiXjsXjXjT under the assumptions of Theorem 2.

The same idea can be used to estimate σm2 by (14) σˆm;n2(0)=1n ∑j=1najm(Yj−XjTbˆLS(s,n))2. The coefficients αs,q(k) can be approximated by αˆs,q(k)=1n ∑j=1n(aj;nk)2pjspjq. Now replacing true D(m) , Lis(m) , bm , σm2 and αs,q(k) in formula (8) by their estimators Dˆn(m) , Lˆnis(m) , bˆLS(m,n) , σˆm,n2 and αˆs,q(k) , one obtains a consistent estimator Σˆn for Σ.

Then (15) Vˆn=Dˆn−1ΣˆnDˆn−1 is a consistent estimator for V.

To get the pilot estimators for the normal mixture model one can use the same approach. Namely, we define μˆm,n(0)=1n ∑j=1najmXj,Σˆm,n(0)=1n ∑j=1najm(Xj−μˆm,n)(Xj−μˆm,n)T as estimates for μm and Σm .

By theorem 4.3 from [10], μˆm,n(0) , Σˆm,n(0) and σˆm,n2(0) are n -consistent estimators for the corresponding parameters of the normal mixture model. This allows one to use them as pilot estimators for the one-step and EM estimators.

Now let us consider estimation of the Fisher information matrix in the case of normal mixture model. Define (16) Iˆil(n,τ)= ∑j=1n∂L(ξj,pj,τ)∂τi∂L(ξj,pj,τ)∂τl, (17) Iˆ(n,τ)=(Iˆil(n,τ))i,l=1R,Iˆ(n)=Iˆ(n,τˆ) where τˆ can be any consistent estimator for τ (e.g. τˆnOS or τˆnEM ).

In the normal mixture model we will denote τ(l)=(b(l),μl,Σl,σl2) , i.e. the set of all unknown parameters which describe the l-th mixture component. Theorem 3.

Assume that the normal mixture model is taken and 1.

σm2>0 0$]]>, Σm are nonsingular for all m=1,…,M ;

Then 1.

There exist 0<c0<C1<∞ such that con≤‖I∗(n,τ)‖≤C1n for all n=1,2,… .

Let Ξn be any random dataset of size n with distribution dependent of an unknown parameter b∈Rd . Recall that a set Bα=Bα(Ξn)⊂Rd is called an asymptotic confidence set of the significance level α if limn→∞P{b∉Bα(Ξn)}=α.