An Algorithm for Constructing Orthogonal
and Nearly-Orthogonal Arrays With
Mixed Levels and Small Runs
Hongquan Xu
Department of Statistics
University of California
Los Angeles, CA 90095
( hqxu@stat.ucla.edu)
Orthogonal arrays are used widely in manufacturing and high-technology industries for quality and
productivity improvement experiments. For reasons of run size economy or exibility, nearly-orthogonal
arrays are also used. The construction of orthogonal or nearly-orthogonal arrays can be quite challeng-
ing. Most existing methods are complex and produce limited types of arrays. This article describes
a simple and effective algorithm for constructing mixed-level orthogonal and nearly-orthogonal arrays
that can construct a variety of small-run designs with good statistical properties ef ciently.
KEY WORDS:
D
-optimality; Exchange algorithm; Interchange algorithm;
J
2
-optimality.
1. INTRODUCTION
Consider an experiment to screen factors that may in uence
the blood glucose readings of a clinical laboratory testing
device. One two-level factor and eight three-level factors
are included in the experiments. The nine factors (Wu
and Hamada 2000, table 7.3) are (A) wash (no or yes),
(B) microvial volume (2.0, 2.5, or 3.0 mL), (C) caras H
2
O
level (20, 28, or 35 mL), (D) centrifuge speed (2,100, 2,300,
or 2,500 rpm), (E) centrifuge time (1.75, 3, or 4.5 minutes),
(F) sensitivity (.1 0, .25, or .50), (G) temperature (25, 30,
or 37
C), (H) dilution ratio (1:51, 1:101, or 1:151), and (I)
absorption (2.5, 2, or 1.5). To ensure t hat all the main effects
are estimated clearly from one another, i t is desirable to
use an
orthogonal array
(OA). The smallest OA found for
one two-level factor and eight three-level factors requires
36 runs. However, the scientist wants to reduce the cost of
this experiment and plans to use an 18-run design. A good
solution then is to use an 18-run
nearly-orthogonal array
(NOA).
The concept of OA dates back to Rao (1947). OAs have
been used widely in manufacturing and high-technology
industries for quality and productivity improvement experi-
ments, as evidenced by many industrial case studies and recent
design textbooks (Myers and Montgomery 1995; Wu and
Hamada 2000). Applications of NOAs have been described
by Wang and Wu (1992), Nguyen (1996b), and the references
cited therein.
Formally, an OA of strength two, denoted by
OA4N 1 s
1
¢ ¢ ¢
s
n
5
, is an
N
n
matrix of which the
i
th column has
s
i
levels
and for any two columns all of their level combinations appear
equally often. An OA is
mixed
if not all
s
i
’s are equal. An
NOA, denoted by
OA
0
4N 1 s
1
¢ ¢ ¢ s
n
5
, is optimal under the
J
2
cri-
terion (de ned in Sec. 2.1). From an estimation standpoint, all
of the main effects of an OA are estimable and orthogonal to
each other, whereas all of the main effects of an NOA are still
estimable, but some are partially aliased with others. Because
balance is a desired and important property in practice, in
this a rticle only balanced
OA
0
4N 1 s
1
¢ ¢ ¢ s
n
5
are considered, in
which all levels appear equally oft en for any column. When an
array is used as a factorial design, each column is assigned to
a factor and each row corresponds to a run. Here the terms
of “array” and “ design,” “row” and “run,” and “column” and
“factor” are freely exchanged.
The purpos e of this article is to present a simple and
effective algorithm for constructing OAs and NOAs with
mixed levels and small runs. The algorithm can ef ciently
construct various designs with good statistical properties.
Section 2 introduces the concept of
J
2
-optimality and other
optimality criteria. Section 3 describes an algorithm for
constructing mixed-level OAs and NOAs. Section 4 gives the
performance and compares the algorith m with others in terms
of speed and ef ciency. Section 5 revisits the blood glucose
experiment, and Section 6 gives concluding remarks.
2. OPTIMALITY CRITERIA
A combinatorial criterion,
J
2
-optimality, is introduced in
Section 2.1. This criterion has the advantages of convenience
for programming and ef ciency for computation. The statisti-
cal justi cation of
J
2
-optimality an d other optimality criteria
is given in Section 2. 2.
2.1 The Concept of J
2
-Optimality
For an
N
n
matrix d
D
6x
ik
7
, weight
w
k
>
0 is assign ed
for column
k
, which has
s
k
levels. For 1
µ i1 j µ N
, let
„
i1 j
4
d
5
D
n
X
k
D
1
w
k
„4x
ik
1 x
jk
51
(1)
where
„4x1 y5
D
1 if
x
D
y
and 0 otherwise. The
„
i1 j
4
d
5
value
measures the similarity between the
i
th and
j
th rows of d.
© 2002 American Statistical Association and
the American Society for Quality
TECHNOMETRICS, NOVEMBER 2002, VOL. 44, NO. 4
DOI 10.1198/004017002188618554
356
ORTHOGONAL AND NEARLY-ORTHOGONAL ARRAYS 357
In particular, if
w
k
D
1 is chosen for all
k
, then
„
i1 j
4
d
5
is the
number of coincidences between the
i
th and
j
th rows. De ne
J
2
4
d
5
D
X
1
µi<jµN
6„
i1 j
4
d
57
2
0
A design is
J
2
-optimal
if it minimizes
J
2
. Obviously, by min-
imizing
J
2
4
d
5
, it is desired that the rows of d be as dissimilar
as possible. The following lemma shows an important lower
bound of
J
2
.
Lemma 1.
For an
N
n
matrix d whose
k
th column has
s
k
levels and weight
w
k
,
J
2
4
d
5 ¶ L4n5
D
2
ƒ
1
"
³
n
X
k
D
1
Ns
ƒ
1
k
w
k
´
2
C
³
n
X
k
D
1
4s
k
ƒ
1
54Ns
ƒ
1
k
w
k
5
2
´
ƒ
N
³
n
X
k
D
1
w
k
´
2
#
1
(2)
and th e equality holds if and only if d is an OA.
The proof is given in Appendix A. From Lemma 1, an OA
is
J
2
-optimal with any choice of weights if it exists, whereas
an NOA under
J
2
-optimality may depend on the choice of
weights.
Example 1.
Consider the 12
10 matrix given in Table 1.
The rst column has three levels, and the other nine columns
have two levels each. For illustration,
w
k
D
1 is chosen for all
k
. First, consider a design comprising the rst ve columns.
Table 1. OA
0
(
12
,
3
1
2
9
)
Run 1 2 3 4 5 6 7 8 9 10
1 0 0 1 0 1 1 1 0 0 0
2 0 1 0 0 1 1 0 0 1 0
3 0 0 1 1 0 1 0 1 1 1
4 0 1 0 1 0 0 1 1 0 0
5 1 0 0 0 0 0 0 0 0 1
6 1 1 1 0 0 0 1 0 1 1
7 1 0 1 1 1 0 0 1 1 0
8 1 1 0 1 1 1 1 1 0 1
9 2 0 0 1 0 1 1 0 1 0
10 2 1 1 0 0 1 0 1 0 0
11 2 0 0 0 1 0 1 1 1 1
12 2 1 1 1 1 0 0 0 0 1
NOTE: The ’ rst ’ ve columns form an OA
(
121 3
1
2
4
)
. The pairs of columns (1, 6) and (1, 10)
are nonorthogonal and have an A
2
value of .167; the pairs of col umns (2, 9), (3, 7), (4, 8),
and (6, 10) are nonorthogonal and have an A
2
value of .111; and all other pairs of columns
are orthogonal.
The coincidence matrix
4„
i1 j
4
d
55
of the 12 rows is
5 3 3 1 2 2 3 1 1 2 3 2
3 5 1 3 2 2 1 3 1 2 3 2
3 1 5 3 2 2 3 1 3 2 1 2
1 3 3 5 2 2 1 3 3 2 1 2
2 2 2 2 5 3 2 2 3 2 3 0
2 2 2 2 3 5 2 2 1 4 1 2
3 1 3 1 2 2 5 3 2 1 2 3
1 3 1 3 2 2 3 5 2 1 2 3
1 1 3 3 3 1 2 2 5 2 3 2
2 2 2 2 2 4 1 1 2 5 2 3
3 3 1 1 3 1 2 2 3 2 5 2
2 2 2 2 0 2 3 3 2 3 2 5
1
and
J
2
is the sum of squares of the elements above the d iag-
onal. It is easy to verify that
J
2
D
330 and that the lower
bound in (2) is also 330 for one three-level and four two-level
columns with
w
k
D
1. Therefore, the rst ve columns form
an
OA4
12
1
3
1
2
4
5
, because the
J
2
value equals the lower bound.
Next, consider the whole array, comprising all 10 columns.
Simple calculation shows that
J
2
D
1
1
284 and that the lower
bound in (2) is 1,260. Therefore, the whole array is not an
OA, because the
J
2
value is greater than the lower b ound.
Now consider the change in the
J
2
value if a column is
added to d or if two symbols are switched in a column. If
a column c
D
4c
1
1 : : : 1 c
N
5
0
is added to d and d
C
is the new
N
4n
C
1
5
design, and if c has
s
k
levels and weight
w
k
, then
for 1
µ i1 j µ N
,
„
i1 j
4
d
C
5
D
„
i1 j
4
d
5
C
w
k
„
i1 j
4
c
51
(3)
where
„
i1 j
4
c
5
D
„4c
i
1 c
j
5
. In addition, if the added column c
is balanced, then it is easy to show that
J
2
4
d
C
5
D
J
2
4
d
5
C
2
w
k
X
1
µi<jµN
„
i1 j
4
d
5„
i1 j
4
c
5
C
2
ƒ
1
Nw
2
k
4Ns
ƒ
1
k
ƒ
1
50
(4)
The summation in the second term on the right side of (4) does
not involve any multiplication, because
„
i1 j
4
c
5
is either 0 or 1.
Therefore, calculating
J
2
4
d
C
5
as in (4) is much faster than by
taking the sum of squares of all
„
i1 j
4
d
C
5
in (3). Now suppose
that two distinct symbols,
c
a
6D
c
b
, in rows
a
and
b
in the added
column are switched. Then all
„
i1 j
4
c
5
are unchanged, except
that
„
a1 j
4
c
5
D
„
j1 a
4
c
5
and
„
b1 j
4
c
5
D
„
j1 b
4
c
5
are switched for
j
6D
a1 b
. Hence
J
2
4
d
C
5
is reduced by 2
w
k
ã4a1 b5
, where
ã4a1b5
D
X
1
µj
6D
a1bµN
6„
a1j
4
d
5
ƒ
„
b1j
4
d
576„
a1j
4
c
5
ƒ
„
b1j
4
c
570
(5)
Calculation of
ã4a1 b5
involves no multiplication, because
both
„
a1 j
4
c
5
and
„
b1 j
4
c
5
are either 0 or 1 . These formulas
provide a n ef cient way to update the
J
2
value and are used
in the algorithm.
TECHNOMETRICS, NOVEMBER 20 02, VOL. 44, NO. 4
358 HONGQUAN XU
2.2 Other Optimality Criteria and
Statistical Justi’ cation of J
2
-Optimality
To gain an understanding of the statistical justi cation of
J
2
-optimality, recall other optimality criteria. It is wel l known
that an
s
-level factor has
s
ƒ
1 degrees of freedom. Commonly
used contrasts are from orthogonal polynomials, especially for
quantitative factors. For example, the orthogonal polynomials
corresponding to levels 0 and 1 of a two-level factor are
ƒ
1
and
C
1, and the orthogonal polynomials co rresponding to lev-
els 0, 1, and 2 of a three-level factor are
ƒ
1, 0, and
C
1 for
linear effects and 1 ,
ƒ
2, and 1 for quadratic effects.
For an
N
n
matrix d
D
6x
ik
7
, whose
k
th column has
s
k
levels, consider the main-effects model
Y
D
‚
0
1
C
X
1
‚
1
C
˜ 1
where
Y
is the vector of
N
observations,
‚
0
is the general
mean,
‚
1
is the vector of all the main effects, 1 is the vector of
1s,
X
1
is the matrix of contrast coef cients for
‚
1
, and
˜
is the
vector of independent random errors. Let
X
1
D
4
x
1
1 : : : 1
x
m
5
and
X
D
4
x
1
=
˜
x
1
˜
1 : : : 1
x
m
=
˜
x
m
˜
5
, where
m
D
P
4s
i
ƒ
1
5
.
In the literature, d is known as the design matrix and X
1
is the model matrix (of the main-effects model). A design
is
D-optimal
if it maximizes
—
X
0
X
—
. It is well known that
—
X
0
X
—
µ
1 for any d esign and that
—
X
0
X
— D
1 if and only if the
original design d is an OA. Wang and Wu (1992) proposed
the
D
criterion
D
D —
X
0
X
—
1
=m
(6)
to measure the overall ef ciency of an NOA. Note that
R
D
X
0
X
is the correlation mat rix of
m
columns of
X
1
.
A good surrogate for the
D
criterion is the
4M1 S5
criterion
(Eccleston and Hedayat 1974). A design is
4M 1 S5-optimal
if it maximizes tr
4X
0
X5
and minimizes tr
64X
0
X5
2
7
among
those designs t hat maximize tr
4X
0
X5
. The
4M 1 S5
criterion is
cheaper to compute th an the
D
criterion and has been used
in the construction of computer-aided designs (see, e.g., Lin
1993; Nguyen 2001). Because all diagonal elements of
X
0
X
are 1s, the
4M 1 S5
criterion reduces to the minimization of
tr
64X
0
X5
2
7
, which is the sum of squares of elements of
X
0
X
,
or, equivalently, to the minimization of the sum of squares of
off-diagonal elements of
X
0
X
. This minimization leads to the
following concept of
A
2
-optimality.
Formally, if
X
0
X
D
6r
ij
7
, let
A
2
D
X
i<j
r
2
ij
1
which measures th e overall aliasing (or nonorthogonality)
between all possible pairs of columns. In particular, for a
two-level design,
A
2
equals the sum of squares of correlation
between all possible pairs of columns, and therefore it is
equivalent to the popular
ave4s
2
5
criterion in the context
of two-level supersaturated designs. A design is
A
2
-optimal
if it minimizes
A
2
. This is a good optimality criterion for
NOAs because
A
2
D
0 if and only if d is an OA. Further,
A
2
-optimality is a special case of the
generalized minimum
aberration
criterion proposed by Xu and Wu (2001) for
assessing nonregular designs.
The statistical justi cation for
J
2
-optimality arises from the
following lemma, which shows an important identity relating
the
J
2
and
A
2
criteria.
Lemma 2.
For a balanced design d of
N
runs and
n
fac-
tors, if the weight equals the number of levels for each factor
4
i.e.,
w
k
D
s
k
5
, then
J
2
4
d
5
D
N
2
A
2
4
d
5
C
2
ƒ
1
N
h
Nn4n
ƒ
1
5
C
N
X
s
k
ƒ
X
s
k
¢
2
i
0
The proof is given in Appendix A. For convenience, the
choice of
w
k
D
s
k
is called
natural weights
.
J
2
-optimality with
natural weights is equivalent to
A
2
-optimality and thus is a
good surrogate for
D
-optimality.
Advantages of the use of
J
2
over
D
,
4M 1 S5
,
A
2
and
ave4s
2
5
as an objective function include the following:
1. It is simple to program.
J
2
works wi th the desig n matrix,
whereas all other criteria work with the model matrix.
2. It is cheap to compute. Neither the calculation of
„
i1 j
4
d
5
in (1) nor that of
„
i1 j
4
d
C
5
in (3) involves any multiplication,
because both
„4x
ik
1 x
jk
5
and
„
i1 j
4
c
5
are either 0 or 1, and this
speeds up the algorithm.
3. It works with columns of more than two levels. Note
that the NOA algorithm of Nguyen (1996b) works only with
two-level columns. To construct an
OA
0
4
18
1
2
1
3
8
5
, for exam-
ple, Nguyen has to use a separate blocking algorithm (see
Nguyen 2001) to divide an
OA4
18
1
2
1
3
7
5
into three blocks.
4. It works wi th any choice of weights. By choosing
proper weights, one can construct d ifferent types of NOAs
with a single algorithm. Note that to construct two types of
OA
0
4
12
1
3
1
2
9
5
’s, Nguyen (1996b) has to code the three-level
column differently in his NOA algorithm. This advantage is
discussed in more detail at the end of Section 4.2.
5. It is very ef cient when the number of runs is less than
the number of parameters, as in the case of supersaturated
designs.
3. AN ALGORITHM
The basic idea of the algorithm is to add columns sequen-
tially to an existing design. The sequential operation is
adopted for speed and simplicity. This operation avoids
an exhaustive search of columns for improvement, which
could be complex and inef cient in computation. The two
operations when adding a column are
interchange
and
exchange
. The interchange procedure, also called the pairwise
switch, switches a pair of distinct symbols in a column. For
each candidate column, the algorithm searches all possible
interchanges and makes an interchange th at reduces
J
2
the
most. The interchange procedure is repeated until a lowe r
bound is achieved or until no further improvement is possible.
The exchange procedure replaces the candidate column by a
randomly selected column. This procedure is allowed to repeat
at most
T
times if no lower bound is achieved. The value
of
T
depends on the orthogonality of the previous design.
If the previous design is an OA, then
T
D
T
1
; otherwise,
T
D
T
2
, where
T
1
and
T
2
are two constants controlled by the
user. With any speci ed weights
w
1
1 : : : 1 w
n
, the algorithm
constructs an
OA
0
4N 1 s
1
¢ ¢ ¢ s
n
5
, in which the rst
n
0
columns
form an
OA4N 1 s
1
: : : s
n
0
5
.
The algorithm proceeds as follows:
1. For
k
D
1
1 : : : 1n
, compute the lower bound
L4k5
accord-
ing to (2).
TECHNOMETRICS, NOVEMBER 2002, VOL. 44, NO. 4
ORTHOGONAL AND NEARLY-ORTHOGONAL ARRAYS 359
2. Specify an initial design d with two columns,
4
0
1 : : : 1
0
1
1
1 : : : 1
1
1 : : : 1s
1
ƒ
1
1 : : : 1s
1
ƒ
1
5
and
4
0
1 : : : 1s
2
ƒ
1
1
0
1 : : :
,
s
2
ƒ
1
1 : : : 1
0
1 : : : 1s
2
ƒ
1
5
. Compute
„
i1 j
4
d
5
and
J
2
4
d
5
by de -
nition. If
J
2
4
d
5
D
L4
2
5
, then let
n
0
D
2 and
T
D
T
1
; otherwise,
let
n
0
D
0 and
T
D
T
2
.
3. For
k
D
3
1 : : : 1n
, do the following:
a. Randomly generate a balanced
s
k
-level column c.
Compute
J
2
4
d
C
5
by (4). If
J
2
4
d
C
5
D
L4k5
, go to (d).
b. For all pairs of rows
a
and
b
with distinct sym-
bols, compute
ã4a1 b5
as in (5). Choose a pair of rows
with largest
ã4a1 b5
and exchange the symbols in rows
a
and
b
of column c. Reduce
J
2
4
d
C
5
by 2
w
k
ã4a1 b5
. If
J
2
4
d
C
5
D
L4k5
, then go to (d); otherwise, repeat (b) until
no further improvement is made.
c. Repeat (a) and (b)
T
times and choose a column c
that produces the smallest
J
2
4
d
C
5
.
d. Add column c as the
k
th column of d, let
J
2
4
d
5
D
J
2
4
d
C
5
, and update
„
i1 j
4
d
5
by (3). If
J
2
4
d
5
D
L4k5
, then let
n
0
D
k
; otherwise, let
T
D
T
2
.
4. Return the nal
N
n
design d, of which the rst
n
0
columns form an OA.
This is an example of a columnwise algorithm. As noted
by Li and Wu (1997), the ad vantage of columnwise instead of
rowwise operation is that the balance property of a design is
retained at each iteration. A simple way of adding an
s
-level
column is to choose a best column from all possible candi-
date columns. However, it is computationally impossible to
enumerate all possible candidate columns if the run size
N
is
not small:
N
N =
2
¢
balanced columns for
s
D
2 and
N
N =
3
¢
2
N =
3
N =
3
¢
balanced columns for
s
D
3. The n umbers grow exponentially
with
N
; for example,
24
12
¢
D
2
1
704
1
156;
32
16
¢
D
601
1
080
1
390;
18
12
¢
12
6
¢
D
17
1
153
1
136; and
27
18
¢
18
9
¢
2
0
3
10
11
. The inter-
change and exchange procedures used in the algorithm yield
a feasible approach to minimizing
J
2
for computational ef -
ciency. An interchange operation searches
N
2
=
4 columns for
s
D
2,
N
2
=
3 columns for
s
D
3, and fewer than
N
2
=
2 columns
for any
s
. The interchange procedure usually involves a few
(typically less than six) iterations. Compared with the size
of all candidate columns, the interchange operation searches
only a rather small portion of the whole space. Thus it is an
ef cient local learning procedure, but often en ds up with a
local minimum. For this reason, global exchange procedures
are also incorporated into the algorithm to allow the search to
move around the whole space and not be limited to a small
neighborhood. As discussed later, the global exchange proce-
dure with moderate
T
1
and
T
2
improves the performance of
the algorithm tremendously.
The values of
T
1
and
T
2
determines the speed and perfor-
mance of the al gorithm. A large
T
i
value allows the algorithm
to spend more effort in searching fo r a good column, which
takes more time. The choice of
T
1
and
T
2
depends on the type
of design to be constructed. For constructing OAs, a moderate
T
1
, say 100, is recommended and
T
2
can be 0; for constructing
NOAs, moderate
T
1
and
T
2
are recommended. More details
are given in the next section.
Remark 1.
Both interchange and exchange algorithms
have been proposed by a number of authors for various pur-
poses. (See Nguyen 1996a and Li and Wu 1997 in the context
of constructing supersaturated designs.)
Remark 2.
The performance of the algorithm may depend
on t he order of levels. Experience suggests that it is most
effective if all levels are a rranged in a de creasing order (i.e.,
s
1
¶ s
2
¶ ¢ ¢ ¢ ¶ s
n
), because the number of balanced columns
is much larger for a higher level than a lower level.
Remark 3.
The speed of the algorithm is maximized
because only integer operations are required if integral
weights are used. For ef ciency and exibility, the alg orithm
is implemented as a function in C and can be called from S.
Both C and S source code s are available from the author on
request.
4. PERFORMANCE AND COMPARISON
This section reports the performance and comparison of the
algorithm with others for the construction of OAs and NOAs.
4.1 Orthogonal Arrays
In the construction of OAs, the weights can be xed at
w
i
D
1, and
T
2
should be 0 because it is unnecessary to con-
tinue adding columns if the current design is not orthogonal.
Here the choice of
T
1
is studied in more detail, because it
determines the speed and performance of the algorithm.
The algorithm is tested with four choices of
T
1
: 1
1
10
1
100
and 1
1
000. For each OA and
T
1
, the algorithm is repeated
1
1
000 times with different random seeds on a Sun SPARC
400-MHz workstation. It either succeeds or fails in construct-
ing an OA each time. Table 2 shows the success rate and
the average time in seconds over 1,000 repetitions. In the
construction of a mixed-level OA, as stated in Remark 2,
the levels are arranged in a decreasing order. Table 2 shows
clearly the trade-off between the success rate and speed, which
depends on the choice of
T
1
. The success rate increases and the
speed decreases as
T
1
increases. A good measure is the num-
ber of OAs constructed per CPU time. The algorithm is least
ef cient for
T
1
D
1 and is more ef cient for
T
1
D
10 or 100
than
T
1
D
1
1
000. Overall, the choice of
T
1
D
100 balances
success rate and speed and so is g enerally recommended.
The construction of OAs continues to be an active re search
topic since Rao (1947) introduced the concept. Construction
methods include combinatorial, geometrical, algebraic, cod-
ing theoretic, and algorithmic approaches. State-of-the-art con-
struction of OAs h as been described by Hedayat, Sloane, and
Stufken (1999). Here the focus is on algorithmic construction
and c omparison wit h existing algorithms.
Many exchange algorithms have been proposed for con-
structing exact
D
-optimal designs. (For reviews, see Cook
and Nachtsheim 1980 and Nguyen and Miller 199 2.) These
algorithms can be used to construct OAs; however, they are
inef cient, and the largest OA constructed and published so
far is
OA4
12
1
2
11
5
(Galil and Kiefer 1980). By modifying the
exchange procedure, Nguyen (1996a) proposed an interchange
algorithm fo r constructing supersaturated designs. His pro-
gram can be used to construct two-level OAs; the largest OA
constructed and published is
OA4
20
1
2
19
5
.
Global optimization algorithms, including simulated anneal-
ing (Kirkpatrick, Gelatt, and Vecchi 1983), thresholding
accepting (Dueck and Scheuer 1990), and genetic algorithms
TECHNOMETRICS, NOVEMBER 20 02, VOL. 44, NO. 4
360 HONGQUAN XU
Table 2. Performance in the Construction of OAs
T
1
D
1 T
1
D
10 T
1
D
100 T
1
D
1
1
000
Array Rate Time Rate Time Rate Time Rate Time
OA(9
1
3
4
) 10000% .001 10000% .001 10000% 0001 10000% 0001
OA(12
1
2
11
) 9306% .002 9305% .002 9509% 0002 9406% 0010
OA(16
1
8
1
2
8
) 202% .002 9608% .006 10000% 0007 10000% 0007
OA(16
1
2
15
) 702% .002 9909% .009 10000% 0008 10000% 0008
OA(16
1
4
5
) 506% .003 1504% .013 1507% 0107 1705% 0889
OA(18
1
3
7
2
1
) 03% .003 4209% .022 8207% 0051 8108% 0260
OA(18
1
6
1
3
6
) 09% .004 1604% .017 1806% 0115 1602% 10110
OA(20
1
2
19
) 0% .005 5409% .029 6304% 0090 6308% 0444
OA(20
1
5
1
2
8
) 0% .003 404% .020 3202% 0107 3305% 0782
OA(24
1
2
23
) 0% .009 901% .055 3004% 0370 2804% 10903
OA(24
1
4
1
2
20
) 0% .008 1605% .056 4505% 0309 4304% 10609
OA(24
1
3
1
2
16
) 0% .008 04% .049 305% 0398 303% 20568
OA(24
1
12
1
2
12
) 0% .006 2801% .059 9808% 0104 9805% 0129
OA(24
1
4
1
3
1
2
13
) 0% .007 04% .045 506% 0383 503% 20507
OA(24
1
6
1
4
1
2
11
) 0% .006 102% .045 1001% 0327 809% 20339
OA(25
1
5
6
) 05% .009 903% .077 1200% 0608 1007% 50989
OA(27
1
9
1
3
9
) 0% .012 1004% .106 9700% 0433 10000% 0450
OA(27
1
3
13
) 0% .013 0% .091 02% 0878 03% 70782
OA(28
1
2
27
) 0% .015 0% .078 104% 0764 08% 50544
OA(32
1
16
1
2
16
) 0% .017 0% .144 8801% 10079 10000% 10158
OA(32
1
8
1
4
2
2
18
) 0% .018 103% .134 3801% 10364 4000% 50644
OA(40
1
20
1
2
20
) 0% .037 0% .258 801% 30146 6809% 130972
NOTE: The entries in the colu mns are the success rate of constructing an OA and the average time in seconds per repetition.
(Goldberg 1989), may be used to construct OAs. These algo-
rithms often involve a large number of iterations and are very
slow to converge. These algorithms have been applied to many
hard problems and are documented to be powerful. However,
they are not effective in the construction of OAs (Hedayat
et al. 1999, p. 3 37); for example, t hresholding accepting
algorithms of Fang, Lin, Winker, and Zhang (2000) and Ma,
Fang, an d Liski (2000) failed to produce any
OA4
27
1
3
13
5
or
OA4
28
1
2
27
5
.
DeCock and Stufken (2000) proposed an algorithm for con-
structing mixed-level OAs via searching some existing two-
level OAs. Their purpose is to c onstruct mixed-level OAs with
as many two-level columns as possible, and their algorithm
succeeded i n constructing several new large mixed-level OAs.
In contrast, the purpose in the present article is to construct
as many nonisomorphic mixed-level OAs (with small runs)
as possible, for which the p roposed algorithm is more ex-
ible and effective. For example, the proposed algorithm is
quite effective in constructing an
OA4
20
1
5
1
2
8
5
that is known
to have maximal two-level columns whereas the algorithm
of DeCock and Stufken fails to produce any
OA4
20
1
5
1
2
7
5
.
Furthermore, the proposed algorithm successfully constructs
several new 36-run OAs not constructed by their algorithm.
Appendix B lists nine new 36-run OAs.
It is interesting to have some head-to-head timing compar-
isons between this and ot her algorithms. A Fedorov exchange
algorithm and Nguyen’s NOA a lgorithm are chosen for
comparison. Cook and Nachtsheim (1980) reported that the
Fedorov exchange algorithm produces the best result but takes
the longest CPU time among several
D
-optimal exchange
algorithms. The Fortran so urce code due to Miller and Nguyen
(1994) is used for the Fedorov algorithm, downloaded from
StatLib (
http://www.lib.stat.cmu.edu
). The Nguye n algorithm
is implemented by replacing
ave4 s
2
5
with
J
2
for convenience
because no source code is available. This modi cation will
affect the speed, but not the ef ciency in constructing OAs.
Table 3 shows the comparisons of the algorithms in terms of
speed and ef ciency. All algorithms were compiled and run
on an iMac PowerPC G4 computer for a fair comparison.
The iMac computer has a 867-MHz CPU, about three times
faster than the Sun workstation d escribed earlier. In the
simulation, the Fedorov algorithm was repeated 1,000 times
for
OA4
12
1
2
11
5
and
OA4
16
1
2
15
5
because it is very slow,
and other algorithms were repeated 10,000 times for all
arrays. Table 3 clearly shows that the proposed algorithm
performs t he best and the Fedorov algorithm performs the
worst in terms of both speed and ef ciency. The Fedorov
algorithm is slow because it uses an exhaustive search of
points for improvement and uses
D
-optimality as the objective
function, which i nvolves real-valued matrix operations. The
(modi ed) Nguyen algorithm is slower than the proposed
algorithm, because the former uses a nonsequential approach
and the latter uses a sequential approach. The nonsequential
approach stops only when no swap i s made on an y column
for consecutive
n
times (where
n
is the number of columns),
whereas the sequential approach stops when it fails to add
an orthogonal column for consecutive
T
1
times. When
T
1
is less than
n
, the sequential approach stops earlier in the
case of failure. This explains why the sequential approach
is faster. Furthermore, the high success rate of the proposed
algorithm shows that the sequential approach is more ef cient
than the nonsequential approach. Note that with the increased
computer power, the Fedorov algorithm succeeds in generating
some
OA4
16
1
2
15
5
’s, whereas the Nguyen algorithm still fails
to generate any
OA4
24
1
2
23
5
in 10,000 repetitions.
4.2 Nearly-Orthogonal Arrays
Wang and Wu (1992) sys tematically studied NOAs and
proposed some general combinatorial construction methods.
Nguyen (1996b) proposed an algorithm for constructing
TECHNOMETRICS, NOVEMBER 2002, VOL. 44, NO. 4
