Advances in Systems Science and Applications (2012) Vol.12 No.4 388-398 Optimal Designs in Random Intercept Model with Heteroscedastic Errors Jing Cheng1 and Rongxian Yue2 1Department of Mathematics, Chaohu College, Anhui 238000, China 2Department of Mathematics, Shanghai Normal University, Shanghai 200234, China Scientific Computing Key Laboratory of Shanghai Universities, and Division of Scientific Computation of E-Institute of Shanghai Universities Abstract This paper considers optimal designs based on the D-, G-, A-, I- and Ds-optimality criteria for a random intercept model with heteroscedastic errors. It is shown that the search of optimal approximate designs can be confined at extreme settings of the design region if heteroscedastic structure satisfies specified conditions. Closed expressions for the optimal proportions are given. Keywords Optimal design, Random intercept model, Heteroscedastic errors, Ide- ntical design 1 Introduction Random coefficient models have been widely used for the researching in the area of biosciences, psychology and population pharmacokinetics, where repeated mea- surements are available from different individuals. These models have been in- troduced by Longford[1], for recent researching we refer to Pena and Yohai[2] and Yu[3]. In recent years, the problem of optimal designs for random coefficient models has attracted growing interest. Schmelter[4-5] showed that optimal de- signs in the linear mixed models could be restricted to the class of group-wise identical designs, and optimal designs in the class of single-group designs were also optimal designs in the larger class of more group designs when the design criteria satisfied some assumptions. Schwabe and Schmelter[6], Schmelter et al[7] and Luoma et al[8] . investigated optimal designs in random intercept model, random slope model and random coefficient cubic regression model, respective- ly. Entholzner et al[9] obtained optimal and efficient designs in mixed models. Debusho and Haines[10] provided V-optimal and D-optimal designs with longi- tudinal data in linear regression models with a random intercept. There are many other results of optimal designs are obtained, such as Wang et al[11], Yu[12] and Wen et al[13]. In this article, we investigate the problem of op- timal designs based on some common optimality criteria for a random intercept model with heteroscedastic errors. In Section 2, we introduce the model with necessary notations. Section 3 provides a lemma which makes it sure that we can confine the search of optimal designs at extreme settings of the design region if the optimality criteria satisfy an assumption. Simple expressions of these optimal Advances in Systems Science and Applications (2012) Vol.12 No.4 389 designs are given in this section. Section 4 introduces some examples. Proof of Lemma 1 is given in Appendix. 2 The Random Intercept Model with Heteroscedastic Errors We investigate a linear regression model on the unit interval with a random intercept and heteroscedastic errors. It is assumed that there are individuals with observations each, and the jth observation of ith individual is described by yij = µi + xijβ + e(xij), i = 1, ..., n; j = 1, ...,mi. (1) Where, xij ∈ [0, 1] is the experimental setting; µi denotes the ith individual effect with unknown mean µ and known variance σ2 µ; β is the unknown slope parameter; observational errors e(xij) are assumed to be heteroscedastic with zero mean and variance σ2/λ(xij), here σ2 is known and λ(xij) is a positive real-valued contin- uous function defined on [0,1]. We assume that COV (µi, µi′) = 0, i ̸= i′ COV (µi, e(xi′j)) = 0, ∀i, i′; COV (e(xij), e(xi′j′)) = 0, (i, j) ̸= (i′, j′). For the ith individual, denote Yi =  yi1 ... yimi  , Xi =  xi1 ... ximi  , e(xi) =  e(xi1) ... e(ximi)  , Fi = (1mi , Xi) Here 1mi is a vector of length mi with all entries equal one. Then the model (1) can be expressed by Yi = Fi ( µ β ) + 1mi(µi − µ) + e(xi) , Fiθ + 1mi(µi − µ) + e(xi), i = 1, . . . , n. By the assumptions we have (µi − µ) ∼ (0, σ2 µ) and Vi , COV (Yi) = σ2diag{1/λ(xi1), . . . 1/λ(ximi)}+σ2 µ1mi1 T mi , σ2(Di+d1mi1 T mi ). Here Di = diag{1/λ(xi1), . . . 1/λ(ximi)} and d = σ2 µ/σ 2. For all n individuals, the vector of all observations can be expressed by Y =  Y1 ... Yn  =  F1 ... Fn  θ +  1m1 0 . . . 0 1mn   µ1 − µ ... µn − µ +  e(x1) ... e(xn)  (2) 390 Jing Cheng:Optimal Designs in Random Intercept Model with Heteroscedastic Errors The design matrix for random intercepts is block diagonal, e.g., 1m1 0 . . . 0 1mn  . Consequently, the covariance matrix of Y is COV (Y ) = diag{V1, . . . , Vn}. The best linear unbiased estimate of θ is given by θ̂ = ( n∑ i=1 F T i V −1 i Fi )−1 n∑ i=1 F T i V −1 i Yi. (3) And we can get COV (θ̂) = ( n∑ i=1 F T i V −1 i Fi )−1 . 3 Optimal Designs In this section, we investigate the optimal designs based on D-, G-, A-, I- and Ds- optimality criteria for the models described in previous section. The D-optimal design minimizes the generalized variance of parameter estimates, the G-optimal design minimizes the maximum variance of the predicted value of the response over the design region, the A-optimal design minimizes the total variance of the parameter estimates, the I-optimal design minimizes the integrated mean squared error and the interest of Ds-optimal design is in estimating the slope. In some practical situations like human or animal pharmaceutics studies or medical diagnostics there are often restrictions, e.g., technical implementations, which force the experiment to be performed with identical regimes for all individ- uals. This means that for each individual the number mi of repeated measure- ments equals m and experimental settings xij = xj are identical across all the individuals. So we only consider identical designs in the following, i.e., mi = m, xi = x1 and hence, Fi = F1, Vi = V for all i. Then the best linear unbiased estimate of θ can be written as θ̂ = ( nF T 1 V −1 1 F1 )−1 F T 1 V −1 1 n∑ i=1 Yi. Here V −1 1 = 1 σ2 ( D1 + d1m1Tm )−1 = 1 σ2 ( D−1 1 − dD−1 1 1m1TmD−1 1 1 + d1TmD−1 1 1m ) = 1 σ2 ( diag{λ(xj)} − dD−1 1 1m1TmD−1 1 1 + d ∑n j=1 λ(xj) ) . Advances in Systems Science and Applications (2012) Vol.12 No.4 391 Without loss of generality, we assume σ2 = 1 in the followings. Furthermore, we will consider approximate designs. For any approximate design ξ of the following form ξ = ( x1, . . . , xp ω1, . . . , ωp ) , 2 < p < m, p∑ j=1 ωj = 1. (4) Denote νk = ∫ 1 0 xkλ(x)dξ(x) = p∑ j=1 ωjx k jλ(xj), k = 0, 1, 2. Then the information matrix corresponding to the design ξ of the form (4) can be expressed by M(ξ) = mn 1 + γν0 ( ν0 ν1 ν1 ν2 + γ(ν0ν2 − ν21) ) . (5) Here we note γ = md. For regression model without any random effects, optimal designs are obtained at extreme settings of the design region and Schwabe et al[6] discussed optimal designs of random intercept models. We can’t use the conclusions in Schwabe et al[6] directly in the random intercept model with heteroscedastic errors, but we have the following lemma. Lemma 1 In the model (2), assume that mi = m, (i = 1, . . . , n) and λ(x) satis- fies the following condition 1 λ(x) ≥ 1− x λ0 + x λ(1) , x ∈ [0, 1]. (6) Where λ0 = λ(0) and λ1 = λ(1). Then for any approximate design ξ of the form (4), there exists an approximate design of the form ξ∗ = ( 0, 1 1− ω, ω ) , 0 < ω < 1. Such that M(ξ∗) ≥ M(ξ). The proof of Lemma 1 can be found in the Appendix. The criteria, Φ(·), considered in this paper are functions of the information matrices which are required to satisfy the following assumptions: A1 Φ(·), is a real-valued function defined on the whole set M of 2×2 symmetric non-negative definite matrices, Φ : M → (−∞,∞]; A2 Φ(·) is monotone (the Loewner order (e.g., Pukelsheim[14], p.101)) on M in the sense that M1,M2 ∈ M,M∞ ≥ M∈ ⇒ ⊕(M∞) ≤ ⊕(M∈). 392 Jing Cheng:Optimal Designs in Random Intercept Model with Heteroscedastic Errors These assumptions are satisfied for most of the common criteria including the D-, G-, A-, I- and Ds-optimality. So, by majorization we can confine the search of optimal designs at extreme settings x = 0 and x = 1 if λ(x) satisfies the condition (6) in Lemma 1. Therefore, in what follows we only consider approximate designs ξ of the form ξ = ( 0 1 1− ω ω ) . (7) For approximate designs of the form (7), we have ν0 = λ1ω + λ0(1− ω), ν1 = ν2 = λ1ω. First, we consider the D-optimality Φ(M(ξ)) = ∣∣M−1(ξ) ∣∣. Note that |M(ξ)| , |M(ω)| = (mn)2λ0λ1ω(1− ω) 1 + γ[ωλ1 + (1− ω)λ0] . It is easy to verify that |M(ω)| is maximized at ω = √ 1 + γλ0/( √ 1 + γλ0 + √ 1 + γλ1) and hence ∣∣M−1(ω) ∣∣ is minimized. Therefore we have Theorem 1 For the model (2) with mi = m (i = 1, . . . , n) and λ(x) satisfying (6) the D-optimal design is ξ∗D = ( 0, 1 1− ωD, ωD ) , ωD = √ 1 + γλ0√ 1 + γλ0 + √ 1 + γλ1 For the Ds-optimality, note that the covariance matrix of θ̂ can be calculated by M−1(ω) = 1 mn(ν0 − ν2) ( 1 + γ(ν0 − ν1) −1 −1 ν0 ν1 ) = 1 mn ( 1 λ0(1−ω) + γ − 1 λ0(1−ω) − 1 λ0(1−ω) 1 λ0(1−ω) + 1 λ1(ω) ) The variance of the estimate for β is given by COV (β̂) = [M−1(ω)]22 = 1 mn [ 1 λ0(1− ω) + 1 λ1ω ]. The variance of β̂ is minimized at ω = √ λ0/( √ λ0 + √ λ1). Therefore we have Theorem 2 For the model (2) with mi = m (i = 1, . . . , n), and λ(x) satisfying (6) the Ds-optimal design is ξ∗Ds = ( 0, 1 1− ωDs , ωDs ) , ωDs = √ λ0√ λ0 + √ λ1 . Advances in Systems Science and Applications (2012) Vol.12 No.4 393 Consider the G-optimality, then Φ(M(ω)) = max x∈[0,1] d(ω, x). Here d(ω, x) is the variance of the predicted value of the response, which is given by d(x, ω) = ( 1 x ) M−1(ω) ( 1 x ) = 1 mn {[ 1 λ0(1− ω) + 1 λ1ω ] x2 − 2x λ0(1− ω) + 1 λ0(1− ω) + γ } . As d(ω, x) is a polynomial of degree 2 with positive leading term, its maximum is attained either x = 0 or x = 1 or both, i.e.,max d(ω, x) = max x∈[0,1] {d(0, ω), d(1, ω)} . Note that d(0, ω) = 1 mn { 1 λ0(1− ω) + γ } is strictly increasing in ω, d(1, ω) = 1 mn { 1 λ1(ω) + γ } is strictly decreasing in ω. Thus min ω∈[0,1] max x∈[0,1] d(x, ω) is attained when d(0, ω) = d(1, ω), i.e., 1 λ0(1− ω) = 1 λ1(ω) . So we have Theorem 3 For the model (2) with mi = m (i = 1, . . . , n), and λ(x) satisfying (6) the G-optimal design is ξG ∗ = ( 0, 1 1− ωG, ωG ) , ωG = λ0 λ0 + λ1 . For the I-optimality, Φ(M(ω)) = ∫ 1 0 d(x, ω)dx = 1 mn [ γ + 1 3λ0(1− ω) + 1 3λ1ω ] . Which is minimized at ω = √ λ0/( √ λ0 + √ λ1) = ωDs . So we have Theorem 4 For the model (2) with mi = m (i = 1, . . . , n), and λ(x) satisfying (6) the I-optimal design is ξ∗I = ( 0, 1 1− ωI , ωI ) , ωI = √ λ0√ λ0 + √ λ1 . 394 Jing Cheng:Optimal Designs in Random Intercept Model with Heteroscedastic Errors For the A-optimality Φ(M(ω)) = tr ( M−1(ω) ) = 1 mn [ 1 λ1ω + 2 λ0(1− ω) + γ ] . It is easy to verify that tr ( M−1(ω) ) is minimized at ω = √ λ0√ λ0+ √ 2λ1 .Therefore we have Theorem 5 For the model (2) with mi = m (i = 1, . . . , n), and λ(x) satisfying (6) the A-optimal design is ξA ∗ = ( 0, 1 1− ωA, ωA ) , ωI = √ λ0√ λ0 + √ 2λ1 . From above discussion, we observe that the G-, Ds-, I- and A-optimal designs only depend on the variances at extreme settings; the D-optimal design depends repeated times m and variance proportion d and error variances at the extreme settings. Note that the particular shape of λ(x) is immaterial for the results, but only its values at 0 and 1, as long as condition (6) is satisfied. Specially, when heteroscedastic structure satisfies λ0 = λ1 = max x∈[0,1] λ(x) , the optimal designs discussed above are independent of the variance ratio d. These optimal designs are the same as the corresponding optimal designs in the linear regression model without any random effects, i.e., ωD = ωG = ωI = ωDs = 1 2 , ωA = √ 2− 1. 4 Examples In this section, we consider three random intercept models with the following heteroscedastic errors λ(x) = x2 + 1, λ(x) = 1 1 + x , λ(x) = x4 + 1 x2 + 1 . It is easy to verify that these three λ(x) satisfy the condition (6). These het- eroscedastic structures are also considered in Chang[15] for D-optimal designs in weighted polynomial regression models. We will give the optimal designs for the three models in terms for the result- s given in Section 3. We also compare the D- and G-optimal designs for the three models with the equireplicated design ω0 = 0.5 which is simultaneously D- and G-optimal for the fixed effects only model (d = 0) in terms of the D- and G-efficiency which are defined as following EffD(ω0) = ( |M(ω0)| |M(ω)| ) 1 2 , EffD(ω0) = max x∈[0,1] d(x, ωG) max x∈[0,1] d(x, ω0) (8) Advances in Systems Science and Applications (2012) Vol.12 No.4 395 Example 1 For the model (2) with mi = m(i = 1, . . . , n), and λ(x) = x2 + 1 , from the theorems in Section 3, we obtain the D-, G-, A-, I- and Ds-optimal proportions as follows: ωD = √ 1 + γ√ 1 + γ + √ 1 + 2γ , ωG = ωA = 1 3 , ωI = ωDs = √ 2− 1. The D- and G-efficiencies defined by (8) of the equireplicated design ω0 = 0.5 are as follows: EffD(ω0) = √ 1 + 2γ + √ 1 + γ√ 4 + 6γ , EffD(ω0) = γ + 1.5 γ + 2 . It is clear that EffD(ω0) decreases strictly in γ and ultimately tends to ( √ 2 + 1)/ √ 6 , and EffD(ω0) increases strictly in γ and ultimately tends to one. Fig.1 shows the plots of these two efficiencies. Fig.1 The efficiencies of EffD(ω0) and EffG(ω0) with different λ Example 2 For the model (2) with mi = m(i = 1, . . . , n),and λ(x) = 1 1+x , the D-, G-, A-, I- and Ds-optimal proportions as follows: ωD = √ 2 + 2γ√ +γ + √ 2 + 2γ , ωG = 2 3 , ωA = 1 2 , ωI = ωDs = 2− √ 2. The D- and G-efficiencies defined by (8) of the equireplicated design ω0 = 0.5 are as follows: EffD(ω0) = √ 2 + 2γ + √ 2 + γ√ 8 + 6γ , EffD(ω0) = γ + 3 γ + 4 . It is clear that EffD(ω0) decreases strictly in γ and ultimately tends to ( √ 2 + 1)/ √ 6 = 0.9856, and EffD(ω0) increases strictly in γ and ultimately tends to one. Fig.2 shows the plots of these two efficiencies. 396 Jing Cheng:Optimal Designs in Random Intercept Model with Heteroscedastic Errors Fig.2 The efficiencies of EffD(ω0) and EffG(ω0) with different λ Example 3 For the model (2) with mi = m (i = 1, . . . , n), and λ(x) = x4+1 x2+1 , the D-, G-, A-, I- and Ds-optimal proportions as follows: ωD = ωG = ωI = ωDs = 1 2 , ωA = √ 2− 1. That is, the D-, G-, I- and Ds-optimal designs are all the equireplicated designs. Appendix Proof of Lemma 1 From Liski et al[16], we get M−1(ξ) = M−1 0 (ξ) + ( d n 0 0 0 ) Here M0(ξ) is the corresponding generalized information matrix when there are no individual intercepts, i.e., M0(ξ) = mn ( ν0 ν1 ν1 ν2 ) . Let the proportion ω in ξ∗ be of the form ω = ν1/λ1. It follows that M0(ξ ∗) = mn ( ν∗0 ν∗1 ν∗1 ν∗2 ) , and M0(ξ ∗)−M0(ξ) = mn ( ν∗0 − ν0 0 0 ν∗2 − ν2 ) . Here ν∗0 = λ1ω + λ0(1− ω) and ν∗2 = ν∗1 = λ1ω . Since 1 λ(x) ≥ 1− x λ0 + x λ1 ≥ x λ1 , Advances in Systems Science and Applications (2012) Vol.12 No.4 397 so λ1 ≥ xλ(x). It implies 0 < ω < 1. [ M0(ξ ∗)−M0(ξ) ] 11 = mn [ λ1ω + λ0(1− ω)− p∑ j=1 ωjλ(xj) ] = mn p∑ j=1 ωj [ λ(xj)xj(λ1 − λ0)− λ1λ(xj) + λ1λ0 ] Condition (6) implies λ(x)x(λ1 − λ0)− λ1λ(x) + λ1λ0 ≥ 0. So we obtain [ M0(ξ ∗)−M0(ξ) ] 11 ≥ 0 By ν∗2 = ν∗1 = ν1 ≥ ν2 , we have[ M0(ξ ∗)−M0(ξ) ] 22 ≥ 0 So we get M0(ξ ∗) ≥ M0(ξ) and hence M(ξ∗) ≥ M0(ξ). Acknowledgements This work was partially supported by a NSFC grant (11071168), Special Funds for Doctoral Authorities of Education Ministry (20103127110002), E-Institutes of Shanghai Municipal Education Commission (E03004), Shanghai Leading Aca- demic Discipline Project (S30405), the Innovation Program of Shanghai Munici- pal Education Commission (11zz116), and the Scientific Research Foundation of Chaohu College. References [1] N.T. Longford (1993), Random Coefficient Regression Models, Clarendon Press, Oxford. [2] D. Pena, V. Yohai (2006), “A dirichlet random coefficient regression mod- el for quality indicators”, Journal of Statistical Planning and Inference, Vol.136, pp.942-961. [3] F. Yu (2007), “Quadratic design criterion for nonlinear models”, Advances in Systems Science and Applications, Vol.7, No.2, pp.155-160. [4] T. Schmelter (2007), “The optimality of single-group designs for certain mixed models”, Metrika, Vol.65, pp.183-193. 398 Jing Cheng:Optimal Designs in Random Intercept Model with Heteroscedastic Errors [5] T. Schmelter (2007), “Consideration on group-wise identical designs for lin- ear mixed models”, Journal of Statistical Planning and Inference, Vol.137, pp.4003-4010. [6] R. Schwabe, T. Schmelter. (2008), “On optimal designs in random intercept models”, Tatar Mountains Mathematical Publications, Vol.39, pp.189-195. [7] T. Schmelter, et al. (2007), “Some curiosities in optimal designs for ran- dom slopes”, In: Contributions to statistics. mODa 8 - Advances in Model- Oriented Design and Analysis,Physica-Verlag Heidelberg, pp.189-195. [8] A. Luoma, et al. (2007), “Optimal designs in random coefficient cubic re- gression models”, Journal of Statistical Planning and Inference, Vol.137, pp.3611-3617. [9] M. Entholzner, et al. (2005), “A note on designs for estimating population parameter”, Listy Biometryczne-Biometrical Letters, Vol.42, pp.25-41. [10] L.K. Debusho, L.M. Haines (2007), “V- and D-optimal population designs for the simple linear regression model with a random intercept term”, Journal of Statistical Planning and Inference, Vol.138, pp.1116-1130. [11] G. Z. Wang, J. Zhao, J. B. Chen (2006), “The uniform design modeling and analysis”, Advances in Systems Science and Applications, Vol.6, No.4, pp.699-703. [12] S.H. Yu (2007), “The linear minimax estimator of stochastic regression coef- ficients and parameters under quadratic loss function”, Statistics and Prob- ability Letters, Vol.77, pp.54-62. [13] G.Y. Wen, et al. (2008), “Orthogonal experimental design in matrix form and its application to the introduction selection of thymus genus”, Advances in Systems Science and Applications, Vol.8, No.3, pp.437-446. [14] F. Pukelsheim (1993), Optimal Design of Experiment, JohnWiley, New York. [15] F.C. Chang (2005), “D-optimal designs for weighted polynomial regression– a functional approach”, Annals of the Institute of Statistical Mathematics, Vol.57, No.4, pp.833-844. [16] E.P. Liski, N.K. Mandal, K.R. Shah, B.K. Sinha (2002), Topics in Optimal Designs, Springer, New York. Corresponding author Nicholas Nechval can be contacted at:yue2@shnu.edu.cn