Feasible criterion for designs based on fixed effect ANOVA model

Statistics and Probability Letters 87 (2014) 134–142

Contents lists available at ScienceDirect

Statistics and Probability Letters journal homepage: www.elsevier.com/locate/stapro

Feasible criterion for designs based on fixed effect ANOVA model Xiaodi Wang a,∗ , Yingshan Zhang b , Yincai Tang b a

School of Statistics and Mathematics, Central University of Finance and Economics, Beijing, China

b

School of Finance and Statistics, East China Normal University, Shanghai, China

article

info

Article history: Received 2 April 2013 Received in revised form 20 January 2014 Accepted 20 January 2014 Available online 27 January 2014

abstract A design is called feasible if all the ANOVA parameters in it are estimable. This paper provides a simple method to judge the feasibility of a design using matrix image. Furthermore, a broad family of feasible designs is obtained. © 2014 Elsevier B.V. All rights reserved.

Keywords: ANOVA Estimable Feasible design Matrix image

1. Introduction Analysis of variance (ANOVA) model with fixed effects has received a lot of attention in the past decades. As an exploratory tool to explain observations, ANOVA model has been widely used in the study of effects of multiple factors for designed experiments. In broad terms, when the m-way ANOVA model is used, there are m factors A1 , A2 , . . . , Am , each having pj different states or levels. Some possible level combinations of the m factors (can be repeated) consist of an experimental design H = (aij )n×m with aij ∈ {1, 2, . . . , pj } being the level of Aj at the ith run. Each combination is applied to each of the experimental units and one value of response (y) is taken on each experimental unit, then the ANOVA model, yi = µ + A1 (ai1 ) + A2 (ai2 ) + · · · + Am (aim ) + IA1 A2 (ai1 , ai2 ) + · · · + IA1 ...Am (ai1 , . . . , aim ) + εi ,

(1)

where µ is the grand mean, Aj (aij ) is the main effect of the aij th level of factor Aj , IA1 A2 (ai1 , ai2 ) is the interaction term of the ai1 th level of A1 and the ai2 th level of A2 , IA1 ···Am (ai1 , . . . , aim ) is the interaction term of A1 , . . . , Am at the level combination (ai1 · · · aim ), εi is the unobserved error of the ith observation with mean = 0 and some variance σ 2 , and i = 1 to n. Consider a two-way ANOVA model with factors A and B, each having 2 levels, for example. If the experimental design with 4 runs is



1 1

1 2

2 1

′

2 2

, then the model can be expressed as:

y1 µ + A(1) + B(1) + IAB (1, 1) ε1 y2  µ + A(1) + B(2) + IAB (1, 2) ε2  y  = µ + A(2) + B(1) + I (2, 1) + ε  . 3 AB 3 y4 µ + A(2) + B(2) + IAB (2, 2) ε4

 

∗





 

Corresponding author. E-mail address: [email protected] (X. Wang).

0167-7152/$ – see front matter © 2014 Elsevier B.V. All rights reserved. http://dx.doi.org/10.1016/j.spl.2014.01.020

(2)

X. Wang et al. / Statistics and Probability Letters 87 (2014) 134–142

135

In most cases, the investigated model is a part of (1) since some interactions are assumed to be nonexistent based on the researchers knowledge of the experimental situation. Therefore, in this paper, Model (1) with one interaction IAi ···Ais (1 ≤ i1 ≤ · · · ≤ is ≤ m) is studied. It is among the simplest that reveals the general structure of the ANOVA model. 1 The principles we present hold, for more complex models involving more interactions. Let Ω = {0, 1, . . . , m, {i1 , i2 , . . . , is }} indicate the summands of the considered main effects and interactions, then the studied model can be expressed in vector notation as: Y = Θ0 +

m 



Θj + Θi1 ···is + ε =

ΘM + ε,

(3)

M ∈Ω

j=1

where

Θ0 = (µ, . . . , µ)′ , Θj = (Aj (a1j ), . . . , Aj (anj ))′ , j = 1, . . . , m, Θi1 ···is = (IAi1 ···Ais (a1i1 , . . . , a1is ), . . . , IAi1 ···Ais (ani1 , . . . , anis ))′ , Y = (y1 , . . . , yn )′ ,

ε = (ε1 , . . . , εn )′ .

The parameters in Model (3) are more than that can be estimated from the observations, and the model is overparameterized (Searle, 1971, 1987, Milliken and Johnson, 1984, Chapter 6, Shaw, 1993). To solve this problem, restrictions are suggested to be imposed on the parameters. Placing restrictions on the parameters can be accomplished in many ways, a generalized sum-to-zero one is considered in this study: p1 

A1 (k) ·

k=1 pi 1 

nk1 n

= 0, . . . ,

pm 

Am (k) ·

nkm

k=1

IAi ···Ais (k, i′2 , . . . , i′s ) · 1

k=1

nki′ ···i′s 2

ni′ ···i′s

n

= 0,

(4)

= 0,

(5)

2

··· pis 

IAi ···Ais (i′1 , i′2 , . . . , k) · 1

k=1

ni′ i′ ···k 1 2

ni′ ···i′ 1

= 0,

(6)

s−1

where nkj is the occurrence number of level k in the jth column of the design, i′1 · · · i′t −1 i′t +1 · · · i′s (t = 1, . . . , s) is any level combination of Ai1 , . . . , Ait −1 , Ait +1 , . . . , Ais appears in the design and ni′ ···i′ i′ ···i′s is the corresponding occurrence number, 1

t −1 t +1

ni′ ···k···i′s is the occurrence number of the level combination i′1 · · · k · · · i′s of Ai1 , . . . , Ais in the design. Note that 1

pit 

ni′ ···k···i′s = ni′ ···i′ i′ ···i′s , 1 1 t −1 t +1

t = 1, 2, . . . , s.

k=1

In particular, when H is balanced (i.e., equal subclass number), these restrictions are the common used sum-to-zero restrictions (Milliken and Johnson, 1984) as: pj 

Aj (k) = 0,

j = 1, . . . , m,

k=1

pi 1  k=1

IAi ···Ais (k, i′2 , . . . , i′s ) = 0, . . . , 1

pis 

IAi ···Ais (i′1 , i′2 , . . . , k) = 0. 1

k=1

To answer the following questions, it is necessary to get good estimates for ΘM through Y : 1. Which is the best level combination of A1 , . . . , Am for producing y at the highest value? 2. Do the levels of factor Aj differ in their effects on the response variable? However, differing designs produce distinct estimates of ΘM . Bad designs may result in the inestimability of the parameters. Herein, the term estimability refers to the existence of linear unbiased estimators of the parameters (i.e. there is a matrix B such that E (BY ) = ΘM ). A design is said to be feasible if all the ΘM in it are estimable, and the ANOVA procedure should be conducted on the feasible designs. The goal of this paper is to solve the problem of selecting feasible designs for ANOVA models. An essential feasibility criterion for designs is the linear-based criterion (Searle, 1971, 1987; Wang et al., 2012), which is proposed by a linear transformation of Model (3) into Y =

 M ∈Ω

CM βM + ε = C β + ε,

(7)

136

X. Wang et al. / Statistics and Probability Letters 87 (2014) 134–142

where βM = (θ1 , . . . , θlM )′ is the vector of lM (< n) independent parameters of ΘM , CM is the n×lM matrix such that CM βM = ΘM , which is obtained by expressing each element of ΘM in terms of βM with the restrictions, C = (C0 , C1 , . . . , Cm , Ci1 ···is ) ′ and β = (β0′ , β1′ , . . . , βm , βi′1 ···is )′ . The superscript prime here indicates the transport of a vector or matrix. Since the estimability of ΘM is equivalent to that of βM , H is feasible if and only if β is estimable. From the theory of linear model, β is estimable if and only if C is full column rank. Hence, the feasibility criterion for H is to verify whether the columns of C have full rank. In spite of its effectiveness in deciding feasible designs, the linear-based criterion may suffer from a cumbersome linear transformation task when the design is complex (unbalanced with large number of levels). Because the restrictions of the ANOVA parameters in this case are complex, for example, the number of restrictions for a three-factor interaction with each factor having four levels might be 48 according to Eqs. (5) and (6). Moreover, it is difficult to see what types of design structures satisfy this criterion. The only result can be found for orthogonal arrays (Cheng, 1978, 1980). This poses a problem for choosing or constructing a feasible design before the process of data analysis, especially when orthogonal arrays are not available. Matrix image (MI) is a powerful tool for studying designs and has been widely employed in the construction of orthogonal arrays (Zhang, 1993; Zhang et al., 1999; Pang et al., 2004). In this paper, we present and prove an interesting equivalent relationship between feasible designs and their MIs. As a result, a broad family of designs, including orthogonal and nonorthogonal types, are found feasible. Meanwhile, this equivalent relationship can be used as an alternative criterion to determine the feasibility of any given design. Compared to the traditional criterion, the criterion is much simpler to use without the need for linear transformation of the ANOVA model. In Section 2, we present the definition and notions for matrix image. In Section 3, we present our main result in terms of MI. This is followed by an example to illustrate the procedure of feasibility judgment using this result. Some kinds of newfound feasible designs are presented in Section 4.

2. Matrix image For the design H = (aij )n×m and D = {l1 , . . . , lh } ⊆ {1, 2, . . . , m}, we denote (ail1 · · · ailh ) as aiD . Suppose G = (g1 , . . . , gn )′ is an n × 1 vector with each element corresponding to each row of H, then an ID is defined as the n × n matrix such that

 1   g  |J | ∈J s  1D     1D s g1   1  gs   g2   |J |   2D s∈J2D  , ID G = ID    ..  =  .. .     . gn      1 gs |JnD | s∈J nD

where JiD = {s : asD = aiD } and |JiD | is the number of elements in JiD (i = 1, . . . , n). Note that |JiD | = naiD . In particular, Ji∅ = {s : s = 1, 2, . . . , n}. Example 2.1. Consider





1 1 H = 2 2

1 2 , 1 2

  y1

y  G = Y =  2 . y 3

y4

For D = ∅, {1}, {2} and {1, 2}, the corresponding I∅ , I{1} , I{2} and I{1,2} are the matrices such that

1

 (y1 + y2 + y3 + y4 )  4  1  (y1 + y2 + y3 + y4 )  4 , I∅ Y =  1   (y1 + y2 + y3 + y4 ) 4    1 (y1 + y2 + y3 + y4 ) 4

1

 (y1 + y2 ) 2  1   (y1 + y2 ) 2  , I{1} Y =  1   (y3 + y4 ) 2    1 (y3 + y4 ) 2

X. Wang et al. / Statistics and Probability Letters 87 (2014) 134–142

1 2 1 

(y1 + y3 )

137



    2 (y2 + y4 ) , I{2} Y =  1   (y1 + y3 ) 2    1 (y2 + y4 )

y  1

y2   I{1,2} Y =  y3  , y4

2

which yields



1 1 1 I∅ =  4 1 1

1 1 1 1



1 1 0 I{2} =  2 1 0



1 1 1 1

0 1 0 1

1 1 , 1 1

1 0 1 0



1 1 0 0

0 0 1 1

0 0 , 1 1



0 1 0 0

0 0 1 0

0 0 . 0 1

1 1 1 I{1} =  2 0 0



0 1 , 0 1

1 0 I{1,2} =  0 0





It is not hard to see that ID is similar to an indicator matrix and can be directly calculated by the design H since it only depends on the rows of H. Actually, the (i, j)th element of ID is

  1 , ID (i, j) = |JiD | 0,

if ajD = aiD , else.

(8)

Then we present the definition of matrix image. Definition 2.1. For the design H = (aij )n×m and D = {l1 , . . . , lh } ⊆ {1, . . . , m}, AD =



(−1)|D|−|N | IN

(9)

N ⊆D

is called the matrix image of columns l1 , . . . , lh of H . Consider Example 2.1, the matrix images for D = ∅, {1}, {2} and {1, 2} are respectively A∅ = I∅ , A{1} = I{1} − I∅ , A{2} = I{2} − I∅ and A{1,2} = I{1,2} − I{1} − I{2} + I∅ . 3. Feasible criterion for H In this section, we give our feasible criterion for H . To prove the criterion, two steps will be done. The first is to establish a relationship between matrix image and the linear model (7), and the second is to obtain the main result in terms of the ranks of matrix images. 3.1. Main theorems First we give and prove a lemma. Lemma 3.1. In Model (7), AM CM = CM for ∀M ∈ Ω . Proof. The restrictions (4)–(6) are equivalent to IMk− ΘM = 0

for M ∈ Ω ,

where Mk− = M \ {ik }, k = 1, . . . , s, if M = {i1 , . . . , is } and Mk− = ∅ if M = {1}, . . . , {m}, which yields IN Θ M =



ΘM

if N = M if N ⊂ M .

0

Thus AM ΘM =



(−1)|M |−|N | IN ΘM = ΘM .

N ⊆M

From ΘM = CM βM , we have AM CM βM = CM βM for any βM . Hence A M CM = CM .  Now we show the main result.

138

X. Wang et al. / Statistics and Probability Letters 87 (2014) 134–142

Theorem 3.1. Under the condition that rank(AM ) ≤ lM , H is a feasible design if and only if



 

rank



=

AM

M ∈Ω

rank(AM ).

(10)

M ∈Ω

Proof. From Lemma 3.1, we have R(CM ) ⊂ R(AM ), and lM = rank(CM ) ≤ rank(AM ). In addition, rank(AM ) ≤ lM according to the theorem condition. Thus we have rank(CM ) = rank(AM ) and R(CM ) = R(AM ). Since PCM and AM are the projection matrices of R(CM ), from the uniqueness of projection matrices, we obtain PCM = AM . Now we prove the result that ‘‘if ’’: since H is a feasible design, that is, ΘM s are estimable, β in Model (7) is estimable. According to Lemma 2.1 of Wang et al. (2012), C has full column rank, i.e. rank(C ) =



lM .

M

′ Thus CM CM are invertible matrices, and



PCM =

M



′ ′ CM (CM CM )−1 CM

M

= (C0 , C1 , . . . , Ci1 ···is )diag ((C0′ C0 )−1 , . . . , (Ci′1 ···is Ci1 ···is )−1 )(C0′ , C1′ , . . . , Ci′1 ···is )′ = CKC ′ . Notice that K = diag ((C0′ C0 )−1 , . . . , (Ci′1 ···is Ci1 ···is )−1 ) is a positive definite matrix and there exists a invertible matrix L such that K = LL′ . Then, we have



PCM = CLL′ C ′ = (CL)(CL)′ .

M

Thus,



 

rank

= rank

AM



 

M

= rank((CL)(CL)′ ) = rank(CL) = rank(C ) =

PCM

 M

M

Furthermore, from rank(CM ) = rank(AM ) = lM , we have



 

rank

=

AM



M

M

‘‘only if ’’: from rank

 rank

rank(AM ) =



lM .

M



M



AM =



M

rank(AM ) and rank(CM ) = rank(AM ) = lM , we have

 

=

AM

M



rank(AM ) =

M



lM .

M

Furthermore,



AM =



=



M

PCM =

M



′ ′ CM (CM CM )−1 CM

M

CM DM = (C0 , C1 , . . . , Ci1 ···is )(D′0 , D′1 , . . . , D′i1 ···is )′

M

= CD, where C = (C0 , C1 , . . . , Ci1 ···is ) and D = (D′0 , D′1 , . . . , D′i1 ···is )′ . Thus

 rank(C ) ≥ rank

  M

AM

=

 M

lM .

lM .

X. Wang et al. / Statistics and Probability Letters 87 (2014) 134–142

On the other hand, the number of the columns of C is

 rank(C ) = rank



M lM ,

and thus rank(C ) ≤



M lM .

139

Hence,

 

=

AM



M

lM ,

M

that is, C has full column rank. Then β is estimable and βM is estimable. Thus ΘM is estimable. That is, H is a feasible design.  From the theorem above, the feasibility of a given design H can be judged only by its matrix images. Steps are as follows: 1. Compute AM for M ∈ Ω by Eqs. (8) and (9). 2. Verify rank(AM) ≤ lM for  M ∈Ω . 3. Compare rank M rank(AM ). If they are equal, then H is feasible. Else, H is unfeasible. M AM and Since AM and their ranks can be obtained by a computer program (such as SAS) following the expression of (8), the method is rapid and easy to apply. In addition, the matrix images are got directly by H itself, and the proposed method does not need the linear transformation for the ANOVA model. Remark 1. For MIs of main effects, the conditions rank(Aj ) ≤ lj hold naturally. Proof can be seen in Corollary 4.1. For MIs of interactions AM , Zhang (1993) showed that rank(AM ) = lM when the design columns belonging to M compose an orthogonal array with strength |M | or a pseudo level orthogonal array by merging some levels in one or more columns of an |M |-strength orthogonal array. Therefore, step 2 above can be skipped when no interaction is considered, or the columns corresponding to each interaction are chosen as an orthogonal array or a pseudo level orthogonal array (levels in other columns and level combinations of the columns from different interactions can be arbitrary). Remark 2. When step 2 is necessary, lM can be got by sequentially deleting the restrained parameters of ΘM according to the restrains (4)–(6). Especially, lM = j∈M (pj − 1) if all possible level combinations of the columns belonging to M are appeared. Remark 3. If the condition rank(AM ) ≤ lM is not satisfied, the conclusion of Theorem 3.1 does not hold. Because we select several designs which do not meet the condition, and find some of them are feasible even (10) is not tenable. 3.2. Example In this subsection, we present an example to illustrate and verify the new feasibility criterion. Example 3.1. Consider a 3-way ANOVA model with factors A, B, C and one interaction IAB in the following design,

1 1 1  1  H = 1 1  2  2 2

1 2 3 1 2 3 1 2 3

1 1 2  1  1 . 2  2  2 3

According to (4)–(6), the restrictions on parameters are: 2A(1) + A(2) = 0,

B(1) + B(2) + B(3) = 0,

IAB (1, 1) + IAB (1, 2) + IAB (1, 3) = 0, 2IAB (1, 1) + IAB (2, 1) = 0,

4C (1) + 4C (2) + C (3) = 0,

IAB (2, 1) + IAB (2, 2) + IAB (2, 3) = 0,

2IAB (1, 2) + IAB (2, 2) = 0,

2IAB (1, 3) + IAB (2, 3) = 0.

Now we judge the feasibility of H by the proposed method. Step 1: Notice that columns 1 and 2, corresponding to IAB , consist of a pseudo level orthogonal array, thus the condition of Theorem 3.1 holds. Step 2: From (9), we obtain

1 1 1   1 1 A ø = 1 9 1  1  1 1

1 1 1 1 1 1 1 1 1

1 1 1 1 1 1 1 1 1

1 1 1 1 1 1 1 1 1

1 1 1 1 1 1 1 1 1

1 1 1 1 1 1 1 1 1

1 1 1 1 1 1 1 1 1

1 1 1 1 1 1 1 1 1

1 1 1  1  1 , 1  1  1 1

140

X. Wang et al. / Statistics and Probability Letters 87 (2014) 134–142

A{1}

A{2}

1 1 1   1  1 = 1 18  1   −2  −2 −2 2 −1 −1   1 2 = −1 9 −1  2  −1 −1 5

1 1 1 1 1 1 −2 −2 −2

−1

−1 −1

2

−1 −1 2

2

−1 5 5 −4 5 5 −4 −4 −4 −4

2  −1  −1   1  2 =  −1 18 −1   −4  2 2

1 1 1 1 1 1 −2 −2 −2

−4 −4 5

−4 −4 5 5 5 −4

−1 2

−1 −1 2

−1 2

−4 2

1 1 1 1 1 1 −2 −2 −2

−1

2 −1 −1 2 −1 −1 2 −1 −1

2 −1 −1 2 −1 −1 2

−1 −1

5  −4   1  5 A{3} = 5 36 −4   −4  −4 −4

A{1,2}

1 1 1 1 1 1 −2 −2 −2

2

2

2

−1 5 5 −4 5 5 −4 −4 −4 −4

2 −1 −1 2 −1 −1 −4 2 2

2 −1 −1 2 2 2 −4

2 −1 −1 2 −1 −1 2

−1 −1

5 5 −4 5 5 −4 −4 −4 −4

−1 −1

−1 −1

−1 −1

−2 −2 −2 −2 −2 −2

−2 −2 −2 −2 −2 −2

4 4 4

4 4 4

2 −1 −1 2 −1 −1 2 −1 −1

−1

1 1 1 1 1 1 −2 −2 −2

−4 −4

4 4

−1  −1  2   −1   −1  , 2   −1   −1

2

−1 −1 2

−1 −1 2

−1

−4 −4

2

−4 −4

5

5

5

−4 −4

−4 −4

−4 −4

5 5 5 −4

5 5 5 −4

5 5 5 −4

−1 2

−1 −1 2

−1 2

−4 2

−1 −1 2 −1 −1 2 2 2 −4

−2 −2 −2  −2  −2 , −2  4  

−4 2 2 −4 2 2 8 −4 −4

2 −4 2 2 −4 2 −4 8 −4

−4 −4 −4  −4  −4 , −4  −4  −4 32 2  2  −4  2   2 . −4  −4  −4 8

Step 3: Calculate

 rank

AM





 



= 6,

M

rank(AM ) = 1 + 1 + 2 + 2 + 2 = 8 ̸= rank



AM

.

M

M

Hence, H is not a feasible design. To verify the result, we transform Example 3.1 into a linear model. From the restrictions, the independent parameters in Θ∅ , Θ1 , Θ2 , Θ3 and Θ12 are respectively µ, A(1), B(1) and B(2), C (1) and C (2), and IAB (1, 1) and IAB (1, 2). By expressing the parameters in the ANOVA model in terms of µ, A(1), B(1), B(2), C (1), C (2), IAB (1, 1) and IAB (1, 2), we obtain

1 1 1  1  Y = C β + ε = 1 1  1  1 1

1 1 1 1 1 1 −2 −2 −2

1 0 −1 1 0 −1 1 0 −1

0 1 −1 0 1 −1 0 1 −1

1 1 0 1 1 0 0 0 −4

0 0 1 0 0 1 1 1 −4

1 0 −1 1 0 −1 −2 0 2

0   µ 1   A(1)  −1    B(1)   0    B(2)   1   C (1)  + ε.  −1    C (2)  0    IAB (1, 1) −2 IAB (1, 2) 2

Since rank(C ) = 6, C is not full column rank. From linear model theory, β is unestimable. Thus H is not a feasible design, which supports the conclusion of the proposed method. 4. Further results In this section, a class of designs is found feasible following Theorem 3.1.

X. Wang et al. / Statistics and Probability Letters 87 (2014) 134–142

141

Theorem 4.1. Under the condition that rank(AM ) ≤ lM , H is feasible if AM AN = 0,

for M ̸= N , M , N ⊆ Ω .

(11)

Proof. From the proof of Theorem 3.1, AM is a projection matrix when rank(AM ) ≤ lM . Given AM AN = 0 for M ̸= N , M , N ⊆

Ω,

′

 

=

AM

M

That is,



′

AM =



M

M

2

 

AM



and

M

=

AM

M



A2M =

M



AM .

M

AM is a projection matrix. Therefore,

 rank

 

AM



 

= tr

M

AM

=

M



tr (AM ) =

M

According to Theorem 3.1, H is feasible.



rank(AM ).

M



Corollary 4.1. If the investigated model is a main-effect model, that is Ω = {0, 1, 2, . . . , m}, then designs with the following structure are feasible, naij ,atk

=

naij natk

1

(12)

n

for i, t = 1, . . . , n and 1 ≤ j ̸= k ≤ m. Proof. From the definition of Aj and the proof of Theorem 3.1, we have pj − 1 ≤ rank(Aj ) ≤ pj . Furthermore, notice that the sum of all the rows of Aj equals 0, thus rank(Aj ) = pj − 1 = lj . While Aj Ak = 0 is equivalent to Ij Ik = (Aj + A∅ )(Ak + A∅ ) = A∅ A∅ = A∅ =

1 ′ 1n 1n , n

which is Eq. (12). Here 1n is the n × 1 vector of 1’s. Following Theorem 4.1, H is a feasible design.



We call this kind of design generalized orthogonal design (GOD) with strength 2. Such designs include orthogonal arrays, pseudo level orthogonal arrays, fractional additional designs, independent contingency tables, orthogonal balanced block designs (Zhang et al., 2009) and so on. Moreover, by using the existing designs, many other such designs can be constructed. For example, they can be constructed by merging orthogonal arrays with suitable fractional additional designs. Corollary 4.2. If the investigated model involving interactions IM1 , . . . , IMk , that is Ω = {0, 1, . . . , m, M1 , . . . , Mk }, then designs with the following structure are feasible, naiM ,atN naiM natN

=

1

(13)

n

for i, t = 1, . . . , n, M ∩ N = ∅, |M + N | ≤ h, where h = maxM ,N ⊆Ω |M ∪ N |. Proof. In this case, (11) holds when AM AN = 0 for all M ̸= N , |M ∪ N | ≤ h, which is equivalent to

 IM IN =

  K ⊆M

AK

 

T ⊆N

AT

=



AK = IM ∩N .

K ⊆M ∩N

It can be verified that this condition is equivalent to Eq. (13).  Moreover, (13) includes the equations for M ∩ N = ∅, M , N ⊆ Mi , i = 1, . . . , k. According to Zhang (1993), rank(AMi ) = j∈Mi (pj − 1) = lMi . Following Theorem 4.1, H is feasible.  We call this kind of design a generalized orthogonal design with strength h. Such designs include orthogonal arrays with strength h, pseudo level orthogonal arrays obtained from h-strength orthogonal arrays. Moreover, using GOD with strength 2, many other such designs can be constructed. For example, let Lsn−1 be a (s − 1)-strength GOD with s − 1 columns, then it is easy to verify that H = ((p) ⊗ 1n , 1p ⊗ Lsn−1 ) is a GOD with strength s, where (p) = (1, . . . , p)′ and ⊗ is the Kronecker product. Acknowledgments The authors are grateful to Professor H. Koul and the anonymous referees for their careful reading of this paper and valuable comments.

142

X. Wang et al. / Statistics and Probability Letters 87 (2014) 134–142

References Cheng, C., 1978. Optimality of certain asymmetrical experimental designs. Ann. Statist. 6, 1239–1261. Cheng, C., 1980. Orthogonal arrays with variable numbers of symbols. Ann. Statist. 8, 447–453. Milliken, G.A., Johnson, D.E., 1984. Analysis of Messy Data. Volume 1: Designed Experiments. Van Nostrand Reinhold, New York, USA. Pang, S.Q., Zhang, Y.S., Liu, S.Y., 2004. Further results on the orthogonal arrays obtained by generalized Hadamard product. Statist. Probab. Lett. 68, 17–25. Searle, S.R., 1971. Linear Models. Wiley, New York, USA. Searle, S.R., 1987. Linear Models for Unbalanced Data. Wiley, New York, USA. Shaw, R.G., 1993. ANOVA for unbalanced data: an overview. Ecology 74, 1638–1645. Wang, X.D., Tang, Y.C., Zhang, Y.S., 2012. Orthogonal arrays for estimating global sensitivity indices of non-parametric models based on ANOVA highdimensional model representation. J. Statist. Plann. Inference 142, 1801–1810. Zhang, Y.S., 1993. Theory of Mulitilateral Matrices. Chinese Statistical Press, Beijing, China. Zhang, Y.S., Lu, Y.Q., Pang, S.Q., 1999. Orthogonal arrays obtained by orthogonal decomposition of projection matrices. Statist. Sinica 9, 595–604. Zhang, Y.S., Wang, X.D., Zhang, X.Q., Pan, C.Y., Chen, X.P., Wang, H., Tian, J.T., Mao, S.S., 2009. Essential theory of orthogonal balanced block designs. in: The 6th Statistics and Probability Seminar across the Taiwan Strait. Nanjing University.