- Email: [email protected]
Subset selection of best simulated systems
ARTICLE IN PRESS
Journal of the Franklin Institute 344 (2007) 495–506 www.elsevier.com/locate/jfranklin
Subset selection of best simulated systems Mahmoud H. Alrefaei, Mohammad Almomani Jordan University of Science and Technology, Department of Mathematics and Statistics, P.O. Box 3030, Irbid 22110, Jordan Received 8 February 2006; accepted 8 February 2006
Abstract In this paper, we consider the problem of selecting a subset of k systems that is contained in the set of the best s simulated systems when the number of alternative systems is huge. We propose a sequential method that uses the ordinal optimization to select a subset G randomly from the search space that contains the best simulated systems with high probability. To guarantee that this subset contains the best systems it needs to be relatively large. Then methods of ranking and selections will be applied to select a subset of k best systems of the subset G with high probability. The remaining systems of G will be replaced by newly selected alternatives from the search space. This procedure is repeated until the probability of correct selection (a subset of the best k simulated systems is selected) becomes very high. The optimal computing budget allocation is also used to allocate the available computing budget in a way that maximizes the probability of correct selection. Numerical experiments for comparing these algorithms are presented. r 2006 The Franklin Institute. Published by Elsevier Ltd. All rights reserved. Keywords: Statistical selection; Simulation optimization; Ordinal optimization; Ranking and selection
1. Introduction We consider the following optimization problem min JðyÞ, y2Y
(1)
where Y, the feasible solution set, is an arbitrary, has no structure, finite but huge set, and J is the expected performance measure of some complex stochastic system, therefore, J can Corresponding author. Tel.: +962 2 7201000; fax: +962 2 7095014.
E-mail addresses: [email protected] (M.H. Alrefaei), [email protected] (M. Almomani). 0016-0032/$30.00 r 2006 The Franklin Institute. Published by Elsevier Ltd. All rights reserved. doi:10.1016/j.jfranklin.2006.02.020
ARTICLE IN PRESS 496
M.H. Alrefaei, M. Almomani / Journal of the Franklin Institute 344 (2007) 495–506
be represented as JðyÞ ¼ E½Lðy; xÞ, where y is a vector representing the system design parameters and x represents all the random effects of the system and L is a deterministic function that depends on x and y. When the search space Y is small, then ranking and selection (R&S) procedures can be used for ranking the systems and select a subset that contains the best s systems with a prespecified significance level, see for example, Rinott [1] for a two stage R&S procedure. R&S cannot be applied for large-scale problems because it needs a huge computational time. To reduce the computational effort, the idea of ordinal optimization (OO) will be used, of course, the objectives will also be changed, instead of looking for the best design, which may be impossible in the case of our setting where Y is huge, we will focus on finding a good enough design. This will save the amount of time spent in simulating each design, see Ho et al. [2], more details about OO will be given in the next section. In this paper, we propose two sequential algorithms for selecting a subset of k systems that is contained in the set of the top s systems. Each iteration of these algorithms consists of two stages; in the first stage, the OO is used to select a subset of g elements randomly from Y. In the first algorithm, we use R&S procedures in the second stage for selecting a small subset from the subset that is selected in the first stage. In the second stage of the second algorithm, we use the optimal computing budget allocation (OCBA) technique that distribute the available computing budget on the designs in such a way to maximize the probability of correct selection (CS). The rest of the paper is organized as follows, in Section 2, we review the R&S procedures for selecting best systems among a small set of alternative systems and we review the idea of OO. In Section 3, we present the OO with the R&S algorithm. In Section 4, we present the second algorithm that uses the OCBA in the second stage. In Section 5, we give some numerical experiments and we make comparisons between the two proposed algorithms. Finally, in Section 6, we give some concluding remarks. 2. Background 2.1. Ranking and selection procedures R&S procedures are used to make comparisons between several alternative systems whose response is normally distributed when the number of competent systems does not exceed 20, see Bechhofer et al. [3] and Law and Kelton [4]. There are three types of R&S procedures. The first procedure can be used for selecting a design that is indifferent from the actual best design by less than a prespecified value d , called the indifference zone, with a prespecified significance level P , P is the probability that we make the CS. The second procedure can be used for selecting a subset of m systems from Y, so that this selected subset contains the best system with high probability. The third procedure is used for selecting a subset that contains the best k systems from Y with a prespecified probability, for more details on ranking and selections, see Bechhofer et al. [3] and Law and Kelton [4]. In this paper, we are interested in the third procedure that can be used for selecting the best k of n systems. For i ¼ 1; 2; . . . ; n, let xi1 ; xi2 ; . . . ; xini be a sample of ni i.i.d. observations from system i, let mi ¼ EðLðyi ; xij ÞÞ, and let mil be the lth smallest of mi ’s, so
ARTICLE IN PRESS M.H. Alrefaei, M. Almomani / Journal of the Franklin Institute 344 (2007) 495–506
497
that mi1 pmi2 p pmin . The CS, here, is to select k systems with the smallest expected means mi1 ; mi2 ; . . . ; mik . If mik and mikþ1 are close to each others, then we might not care on choosing the system whose mean is mikþ1 as one of the best k systems. We want to avoid making a large amount of simulation time to resolve this unimportant difference. More precisely, we are interested in selecting k systems that are indifferent from the best k systems by less than d ; i.e. fall in the region ½mi1 ; mik þ d , where the indifference d is prespecified by the user. We want the CS to take place with a probability greater than or equal to P , where P is prespecified by the user. Therefore, with probability P , the expected value of each of the k selected systems will not exceed mik þ d . 2.1.1. A two-stage R&S (Rinott’s) procedure The two-stage procedure of Rinott [1] is widely used in R&S. In the first stage, each design is sampled using n0 -independent simulation runs and a first stage sample mean j i ðn0 Þ and sample variance S 2i ðn0 Þ are computed for each design i ¼ 1; . . . ; n. Based on the initial simulation runs n0 , and the sample variance estimate obtained in the first stage, the number of additional simulation runs for each design in the second stage is computed by the formula 2 2 h S ðn0 Þ N i ¼ max n0 þ 1; 1 i 2 , (2) ðd Þ where due is the smallest integer greater than or equal u, d is the indifference zone, and h1 is a constant that solves the Rinott’s integral Z 1 F ðh þ tÞn1 f ðtÞ ¼ P , 1
where f and F denote the probability density and cumulative distribution functions, respectively, of the t distribution with n0 1 degree of freedom. This choice of h guarantee the PðCSÞ ¼ P . h can be obtained from tables of Wilcox [5] or Law and Kelton [4]. Then N i n0 additional simulation runs will be done in the second stage, and the overall sample mean j^i will be computed. The design that has the smallest j^i will be selected as the best one. This two-stage can be used also for selecting the best k of n designs. Chen and Kelton [6] consider reducing the value of N i by replacing the value of h in Eq. (2) by hi for each alternative i, where hi ¼ h=ri and ri ¼ maxfmi mi1 ; d g=d . R&S requires a huge amount of time to guarantee the CS with high probability. This huge amount of computational time restrict the size of the problems that these approaches can solve, in fact they can solve a small size problems, say less than 20. Nelson et al. [7] and Alrefaei and Alawneh [8] propose using subset selection in the first stage to reduce the competent designs. These designs are carried out to the second stage in which the indifference zone approach is used where the information in the first stage are used. 2.2. Ordinal optimization procedure The OO approach has been proposed by Ho et al. [2]. The aim of this approach is to find good, better, or best systems, i.e. ‘‘good enough’’ solution, rather than estimating the performance value of these systems accurately. In many real world problems, one could compromise for the ‘‘good enough’’ rather than the best alternative, especially when the
ARTICLE IN PRESS 498
M.H. Alrefaei, M. Almomani / Journal of the Franklin Institute 344 (2007) 495–506
number of competent alternative systems is very large and it is impossible to simulate all alternatives. The main idea of the OO approach is based on two issues: (1) It is much easier to determine ‘‘order’’ than ‘‘cardinal’’. This is intuitively reasonable, to determine whether A is greater or less than B, is a simpler task than to determine the value of A B, in stochastic situations. (2) Softening the goal of optimization also makes the problem easier. Instead of asking the ‘‘best for sure’’ we settle for the ‘‘good enough with high probability’’. Suppose that the CS is to select a subset G of g systems from Y that contains at least one of the top m% best designs. Since we assume that Y is very huge then the probability of CS is given by m g PðCSÞ 1 1 . 100 Now, suppose that the CS is to select a subset G of g designs that contains at least r of the best s systems. Let S be the subset that contains the actual best s systems, then the probability of CS can be obtained using the hypergeometric distribution as ! ns s g gi X i ! PðCSÞ ¼ PðjG \ SjXrÞ ¼ . n i¼r g Since we assume that jYj ¼ n, and n is very large then PðCSÞ can be approximated by the binomial random variable. Therefore, g X g m i m gi PðCSÞ 1 , 100 100 i i¼r where we assume that s=n 100% ¼ m%. It is also clear that this PðCSÞ increases when the sample size g increases. In the literature, PðCSÞ is denoted as the alignment probability. 3. Ordinal optimization with ranking and selection method In this section, we present the proposed procedure for selecting a subset of k designs that belong to the actual best s designs. The method of selecting one design from the actual top m% designs is presented by Alrefaei and Abdul-Rahman [9]. We first present the method of Alrefaei and Abdul-Rahman for selecting one design from the best m% designs, then we present the proposed algorithm for selecting a subset that contains k designs from the best m% designs. 3.1. A sequential algorithm for selecting one design from the best m% designs This method consists of a sequential procedure. Each iteration of the method is based on a two stage procedure. In the first stage, randomly select a subset G, where jGj ¼ g from
ARTICLE IN PRESS M.H. Alrefaei, M. Almomani / Journal of the Franklin Institute 344 (2007) 495–506
499
the search space Y. In the second stage, R&S is applied for selecting the best design yb in the selected subset G with a specified significance level (the probability of selecting the best design) P X1 a. yb will be carried out to the next iteration and the other g 1 designs will be replaced by newly sampled designs from Y. The algorithm continue in this manner until a specified number of generated samples is exceeded or a desired significance level is achieved. Note that, in the first stage, the probability that the randomly selected subset G contains one of the best m% designs can be evaluated as follows: PðG contains one of the best m% designs) 1 ð1 m=100Þg . Let bi be the selected best design in iteration i of the algorithm, and let G0i be the selected subset in iteration i, where jG 0i j ¼ g 1. Let P ¼ 1 a be the significance level used in the R&S procedure (the second stage of the algorithm, i.e. the probability of selecting the best design in the second stage is P ¼ 1 a). Then a lower bound of the probability of selecting one of the top m% best designs in the first iteration is m g P1 ¼ ð1 aÞ 1 1 100 and a lower bound of the probability of correct selection in the jth iteration is given by Pj ¼ Pðbj1 is in the top m% best designs or G 0j contains one of the top m% best designs)ð1 aÞ. Which equals to Pj ¼ [Pðbj1 is in the top m% best designs) + PðG 0j contains one of the top m% best designs) Pðbj1 is in the top m% best designs) PðG 0j contains one of the top m% best designs)]ð1 aÞ. As an example, consider m% ¼ 10% and 1 a ¼ 0:95 then P1 ¼ 0:619, P2 ¼ 0:808, P3 ¼ 0:925. This means that in three iterations, we can select a design that belongs to the top 10% best designs with probability greater than or equal to 0.925 regardless of the size of the optimization problem given in Eq. (1). 3.2. Subset selection approach In each iteration of this approach, we also use two stages, in the first stage, we select a random subset G of g systems from Y, in the second stage, we use the R&S procedure to select a subset of the kp20 best designs of G with indifference zone d . We repeat this procedure until the probability of CS becomes large enough. The algorithm is described as follows: Algorithm 2 (OO þ R&S). Step 0: Select a subset G randomly from Y, such that jGj ¼ g, gp20. Step 1: Generate n0 i.i.d. samples xij , j ¼ 1; 2; . . . ; n0 , for all i ¼ 1; . . . ; g. Step 2: For i ¼ 1; . . . ; g, compute the first stage sample means J ð1Þ i ðn0 Þ ¼
n0 1X Lðyi ; xij Þ, n0 j¼1
and the sample variances S2i ðn0 Þ ¼
n0 1 X 2 ½Lðyi ; xij Þ J ð1Þ i ðn0 Þ n0 1 j¼1
ARTICLE IN PRESS M.H. Alrefaei, M. Almomani / Journal of the Franklin Institute 344 (2007) 495–506
500
and let
h2 S 2 ðn0 Þ N i ¼ max n0 þ 1; 3 i 2 ðd Þ
,
where due is the smallest integer greater than or equal u, d is the indifference zone, and h3 depends on g, k, n0 , and P , and can be obtained from Table 10.13 of Law and Kelton [4]. Step 3: For i ¼ 1; . . . ; g, generate N i n0 i.i.d. samples xij . Compute the overall sample means J^ i ðN i Þ as follows: Ni 1 X J^ i ðN i Þ ¼ Lðyi ; xij Þ. N i j¼1
Step 4: Select the best k systems from the set G, call this subset as G 0 . Step 5: If a stopping rule is not satisfied; randomly select a subset G 00 of g k alternatives from YnG 0 , let G ¼ G 0 [ G 00 , and Go to Step 1. Remark. Note that the R&S procedures can be applied only when the number of alternatives is small, say less than 20, this is why in Step 0, we select gp20. Let S be the subset that contains the actual s best systems, then the calculation of the probability of CS is very complicated. As an example if k ¼ 2, then a lower bound of the probability of CS in the first iteration, can be found as follows: !3 2 ns s 6X gi 7 i 6 g 7 7 P . ! P1 ðCSÞ ¼ PðjG \ SjX2Þ P ¼ 6 6 7 n 4 i¼2 5 g In the ith iteration, the probability of CS can be calculated as follows: Pi ðCSÞ ¼ PðjG \ SjX2Þ P ¼ P½jðG0 [ G 00 Þ \ SjX2 P ¼ P½ðjG0 \ SjX2 or jG 00 \ SjX2 or ðjG 0 \ SjX1 and jG 00 \ SjX1ÞÞ P . Let A be the event A ¼ jG 0 \ SjX2, B be the event B ¼ jG00 \ SjX2, and C be the event C ¼ ðjG 0 \ SjX1 and jG00 \ SjX1Þ then by assuming that the events A, B, and C are mutually independent, we get Pi ðCSÞ ¼ PðA [ B [ CÞ P ¼ ½PðAÞ þ PðBÞ þ PðCÞ PðAÞPðBÞ PðAÞPðCÞ PðBÞPðCÞ þ PðAÞPðBÞPðCÞ P . It is clear that PðAÞ ¼ Pi1 ðCSÞ,
ARTICLE IN PRESS M.H. Alrefaei, M. Almomani / Journal of the Franklin Institute 344 (2007) 495–506
501
! n2s s g2 g2i X i ! . PðBÞ ¼ n2 i¼2 g2 Note also that PðjG 0 \ SjX1ÞXPi1 ðjG0 \ SjX2Þ ¼ Pi1 ðCSÞ and PðjG 00 \ SjX1Þ ¼
1 1
s g2 , n2
so that PðCÞ ¼ ðPi1 ðCSÞÞ 1 1
s g2 . n2
It is clear that PðCSÞ increases as the number of iterations increases. As an example, if k ¼ 2, jYj ¼ n ¼ 1000 systems, s ¼ 100, and g ¼ 10. Then one can evaluate PðCSÞ in each iteration as follows. In the first iteration, note that ! ! ns s PðjG \ SjX2Þ ¼
g X
gi ! n
i
i¼2
100 ¼
10 X
g !
i
1000 100
!
10 i ! 1000
i¼2
10 ¼ 0:263703353. If we use the binomial approximation, this probability becomes PðjG \ SjX2Þ ¼ 0:26390. Let Pi ðCSÞ be the probability of selecting at least two systems that belong to the top 10% systems in iteration i, let P ¼ 0:95 in Step 2, then P1 ðCSÞ is given by ! !1 0 100 1000 100 BX C B 10 C i 10 i C ð0:95Þ ! P1 ðCSÞ ¼ B B C 1000 @ i¼2 A 10 ¼ ð0:263703353Þ ð0:95Þ, P1 ðCSÞ ¼ 0:25052. By simple calculations, we can find a lower bound of the probability of CS in the second iteration, P2 ðCSÞ ¼ 0:45387 and P3 ðCSÞ ¼ 0:63739. For a lower bound of Pi ðCSÞ,
ARTICLE IN PRESS M.H. Alrefaei, M. Almomani / Journal of the Franklin Institute 344 (2007) 495–506
502
Table 1 The lower bound of the probability of correct selection Pk ðCSÞ, k ¼ 1; . . . ; 10, when a subset of the best 10% systems is desired k
1
2
3
4
5
6
7
8
9
10
pk
0.251
0.454
0.637
0.772
0.851
0.891
0.909
0.916
0.919
0.920
i ¼ 1; . . . ; 10, see Table 1. The maximum value of the probability of CS for this choice of parameters does not exceed 0:9211, which can be reached at iteration 13. To increase this value, one needs to increase the value of P . 4. Ordinal optimization with optimal computing budget allocation method The OCBA approach is highly efficient, and significantly reduces the total simulation cost. This approach has been proposed by Chen [10], Chen et al. [11], and Chen et al. [12]. The goal of this approach is to allocate the total simulated samples from all systems in a way that maximizes the probability of selecting the best system within a given computing budget. Suppose that we are allowed to take a total of T samples from all systems for solving the optimization problem given in Eq. (1). The question is where to allocate these computing simulated samples to maximize the CS. Where the CS is to select a system b that has the best performance. Furthermore, instead of equally simulating all systems, we will further improve the performance of OO by determining the best numbers of simulated samples for each system. Suppose that Y ¼ fy1 ; . . . ; yn g, then we wish to choose T 1 ; T 2 ; . . . ; T n such that PðCSÞ is maximized, subject to a limited computing budget T. max
T 1 ;...;T n
s:t:
PðCSÞ n X
T i ¼ T,
i¼1
Pn where T i 2 N for i ¼ 1; 2; . . . ; n and N is the set of non-negative integers and i¼1 T i denotes the total computational samples, assuming the simulation times for different systems are roughly the same. The proof of the following theorem has been given by Chen et al. [12]. Theorem 1. Given a total number of simulated samples T to be allocated to n competing systems whose performance is depicted by random variables with means Jðy1 Þ; Jðy2 Þ; . . . ; Jðyn Þ, and finite variances s21 ; s22 ; . . . ; s2n , respectively. As T!1, the approximate probability of CS can be asymptotically maximized when (1) T i =T j ¼ ððsi =db;i Þ=ðsj =db;j ÞÞ2 , i, j 2 f1; 2; . . . ; ng, and iajab. qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi Pn 2 2 (2) T b ¼ sb i¼1;iab T i =si , where T i is the number of samples allocated to system i, db;i the estimated difference between the performance P i of the two designs ðdb;i ¼ J b J i Þ, and J b pmini J i for all i. Here J i ¼ ð1=T i Þ Tj¼1 Lðyi ; xij Þ, where xij is a sample from xi for j ¼ 1; . . . ; T i .
ARTICLE IN PRESS M.H. Alrefaei, M. Almomani / Journal of the Franklin Institute 344 (2007) 495–506
503
Now, we present the second algorithm in which we use the OCBA in the second stage for selecting a subset of the k best systems in G instead of R&S. Algorithm (OO þ OCBA). Step 0: Determine the size of the subset G, gXk. Choose the number of initial simulations samples t1 X2 for each system. Let l ¼ 0 and let T l1 ¼ T l2 ¼ ¼ T lg ¼ t1 . Step 1: Select a subset G of g alternatives randomly from Y. Take random samples of t1 observations xij ð1pjpt1 Þ for each system i in G, where i ¼ 1; . . . ; g. Step 2: Calculate the first stage sample mean J ð1Þ i for all i ¼ 1; . . . ; g Tl
l J ð1Þ i ðT i Þ
i 1 X ¼ l Lðyi ; xij Þ. T i j¼1
Step 3: Order the systems in G according to their sample averages ð1Þ l l ð1Þ l J ð1Þ 1 ðT 1 ÞpJ 2 ðT 2 Þp pJ g ðT g Þ. 0 Step 4: Select Pg thel best k systems from the set G, call this subset as G . Step 5: If i¼1 T i XT, stop. Step 6: Randomly select a subset G 00 of g k alternatives from Y G 0 , let G ¼ G 0 [ G 00 . Step 7: Increase the computing budget by D and compute the new budget allocation, lþ1 lþ1 T lþ1 1 ; T 2 ; . . . ; T g , according to Theorem 1. Step 8: Perform additional maxf0; T lþ1 T li g samples for each system i, i ¼ 1; . . . ; g, let i l l þ 1. Go to Step 2.
5. Numerical example Suppose we have 10,000 different systems Y ¼ y1 ; . . . ; y10;000 , where Jðyi Þ ¼ i=100 and Lðyi ; xÞ is the normal distribution NðJðyi Þ; 1Þ, with mean m ¼ Jðyi Þ and variance s2 ¼ 1. The optimization problem is min
i¼1;...;10;000
Jðyi Þ.
We want to select a subset that is contained in the best 10% systems, it is obvious that y1 ; y2 ; . . . ; y1000 are the actual best 10% systems. We apply our algorithms where we use R&S in each iteration with significance level P ¼ 0:95, the initial sample size t1 ¼ 40, and the indifference zone d ¼ 0:3. The maximum computing budget in the second algorithm is T ¼ 1500, D ¼ 100, and t1 ¼ 40. In both algorithms, we randomly select a set that consists of 10 systems in the first iteration, and select a subset G 0 that contains the best two systems among these 10 systems using the two proposed algorithms. Then we replace the other eight systems by eight new systems randomly selected from the Y G0 . We repeat this Table 2 The probability of selecting two alternatives from the top 10% of n ¼ 10; 000 systems using g ¼ 10, t1 ¼ 40 Iteration
1
2
3
4
5
6
7
8
9
10
OO þ R&S OO þ OCBA
0.23 0.26
0.55 0.6
0.74 0.8
0.88 0.88
0.95 0.94
0.96 0.99
0.99 0.99
0.99 0.99
1.0 0.99
1.0 1.0
ARTICLE IN PRESS 504
M.H. Alrefaei, M. Almomani / Journal of the Franklin Institute 344 (2007) 495–506
Table 3 The probability of selecting five alternatives from the top 10% of n ¼ 10; 000 systems using g ¼ 10, t1 ¼ 40 Iteration
2
4
6
8
10
12
14
15
OO þ R&S OO þ OCBA
0.0 0.0
0.04 0.13
0.21 0.37
0.44 0.55
0.56 0.72
0.72 0.8
0.88 0.88
0.92 0.92
Table 4 The probability of selecting two alternatives from the top 1% of n ¼ 10; 000 systems using g ¼ 10, t1 ¼ 40 Iteration
1
5
10
15
20
25
30
40
50
OO þ R&S OO þ OCBA
0.01 0.0
0.06 0.07
0.2 0.25
0.32 0.43
0.48 0.52
0.66 0.59
0.78 0.69
0.87 0.82
0.92 0.93
experiment 100 replications and in each replication we use only 10 iterations. Table 2 contains the probability of CS obtained after 10 iterations over the 100 replications by using the R&S and OCBA. It is clear that the OCBA gives better results in almost all iterations. In all cases, we use almost the same computational time for both algorithms. Now, suppose that we want to select a subset that contains five alternatives from the top 10% alternative systems using the same parameters in the previous setting. Table 3 includes the results for this experiment. Again the OO þ OCBA algorithm gives better results. The OO þ R&S also behaves well, it starts slow, but it moves faster in later iterations. We repeat the two algorithms with the objective of selecting a set that contains two systems among the best 1% systems, it is obvious that y1 ; y2 ; . . . ; y100 are the actual best 1% systems. For the ðOO þ R&SÞ algorithm, we use a significance level P ¼ 0:95, the initial sample size t1 ¼ 40, and the indifference zone d ¼ 0:6. The maximum computing budget for the ðOO þ OCBAÞ is T ¼ 1500, D ¼ 100, and t1 ¼ 40. In both algorithms, we select a subset that consists of 10 systems in the first iteration and select a subset G0 that contains the best two systems among these 10 systems, and then replace the other eight systems by eight new systems randomly selected from the search space Y G0 . We run this experiment for 100 replications and in each replication we use 50 iterations. Table 4 contains the probability of CS obtained after 50 iterations over the 100 replications. It is clear that both of the algorithms quickly locate the desired set and that the OO þ OCBA starts faster and then the OO þ R&S accelerates later. To see the robustness of the proposed algorithms, we have applied these algorithms on the following example. Suppose, we have nðM=M=1Þ queuing systems and we want to select a subset of k, where k ¼ 1; . . . ; 5 systems that is contained in the best 10% systems, where the best system is the system that has the minimum average waiting time per customer. The arrival rate is assumed to be a fixed number, l ¼ 1, and the service rate is m 2 Y ¼ f4:0; 4:0001; 4:0002; . . . ; 5:0g, so we have 10,000 different ðM=M=1Þ queuing systems. In the first stage, we use the OO to select a set G of 10 systems and in the second stage, we use the R&S procedures to select the best k systems from the set G with indifference zone d ¼ 0:05. This experiment is repeated 1000 replications, each replication
ARTICLE IN PRESS M.H. Alrefaei, M. Almomani / Journal of the Franklin Institute 344 (2007) 495–506
505
Fig. 1. The performance of the OO þ R&S algorithm when they are applied to select a subset of k, k ¼ 1; . . . ; 5 systems from the top 10% M/M/1 queuing systems.
Fig. 2. A comparison between the OO þ R&S and the OO þ OCBA algorithms when they are applied to select a subset of 2 of the best 10% M/M/1 queuing systems.
ARTICLE IN PRESS 506
M.H. Alrefaei, M. Almomani / Journal of the Franklin Institute 344 (2007) 495–506
consists of 20 iterations; the probability of CS is shown in Fig. 1. It is clear from the figure that the probability of CS increases very rapidly. To compare the OO þ R&S algorithm with the OO þ OCBA algorithm, we have implemented the OO þ OCBA algorithm to select a subset of two systems that belongs to the best 10% systems. We use initial sample size of t1 ¼ 10 samples only in this experiment. Fig. 2 shows the relation between the probability of CS when we use the OO þ R&S algorithm and the OO þ OCBA for 1000 replications for the first 20 iteration. It is clear that the OCBA gives better results with computational time approximately the same as in R&S. 6. Conclusions We have presented two sequential algorithms for selecting a subset of k systems from the top s systems when the number of alternatives is very large. Each iteration of the algorithm consists of two stages. In the first stage, we randomly select a subset G of g alternatives. OO is used for evaluating the probability that this set is contained in the best s alternatives set. In the second stage, we use one of the two procedures, in the first procedure we use the R&S procedure for ranking the alternative systems in the selected subset in the first stage and then select a subset that contains the best k systems of this subset. The other g k alternatives are replaced by newly g k alternatives selected randomly from the search space and the algorithm is repeated until a stopping criterion is reached. In the second procedure, we use the OCBA, in which the number of samples are distributed over the competent alternatives in a way that maximizes the probability of CS. We apply these algorithms on two numerical examples. From our numerical results, we note that these algorithms achieve the CS with high probability, quickly. References [1] Y. Rinott, On two-stage selection procedures and related probability–inequalities, Commun. Statist.: Theory Methods A 7 (1978) 799–811. [2] Y.C. Ho, R.S. Sreenivas, P. Vakili, Ordinal optimization of DEDS, J. Discrete Event Dynamic Syst. 2 (1992) 61–88. [3] R.E. Bechhofer, T.J. Santner, D.M. Goldsman, Design and Analysis of Experiments for Statistical Selection, Screening, and Multiple Comparisons, Wiley, New York, 1995. [4] A.M. Law, W.D. Kelton, Simulation Modeling and Analysis, third ed., McGraw-Hill, New York, 2000. [5] R.R. Wilcox, A table for Rinott’s selection procedure, J. Quality Technol. 16 (2) (1984) 97–100. [6] E.J. Chen, W.D. Kelton, An enhanced two-stage selection procedure, in: J.A. Joines, R.R. Barton, K. Kang, P.A. Fishwick (Eds.), Proceedings of the 2000 Winter Simulation Conference, Institute of Electrical and Electronics Engineers, Piscataway, New Jersey, 2000, pp. 727–735. [7] B.L. Nelson, J. Swann, D. Goldsman, W. Song, Simple procedures for selecting the best simulated system when the number of alternatives is large, Oper. Res. 49 (2001) 950–963. [8] M.H. Alrefaei, A.J. Alawneh, Selecting the best stochastic system for large scale problems in DEDS, Math. Comput. Simulation 64 (2004) 237–245. [9] M.H. Alrefaei, H. Abdul-Rahman, Two sequential algorithms for selecting one of the best simulated systems, WSEAS Trans. Syst. 3 (2004) 2517–2522. [10] C.H. Chen, An effective approach to smartly allocated computing budget for discrete event simulation, IEEE Conference on Decision and Control, vol. 34, 1995, pp. 2598–2605. [11] C.H. Chen, C.D. Wu, L. Dai, Ordinal comparison of heuristic algorithms using stochastic optimization, IEEE Trans. Robotics Autom. 15 (1) (1999) 44–56. [12] C.H. Chen, E. Yucesan, S.E. Chick, Simulation budget allocation for further enhancing the efficiency of ordinal optimization, J. Discrete Event Dynamic Syst.: Theory Appl. 10 (2000) 251–270.
Journal of the Franklin Institute 344 (2007) 495–506 www.elsevier.com/locate/jfranklin
Subset selection of best simulated systems Mahmoud H. Alrefaei, Mohammad Almomani Jordan University of Science and Technology, Department of Mathematics and Statistics, P.O. Box 3030, Irbid 22110, Jordan Received 8 February 2006; accepted 8 February 2006
Abstract In this paper, we consider the problem of selecting a subset of k systems that is contained in the set of the best s simulated systems when the number of alternative systems is huge. We propose a sequential method that uses the ordinal optimization to select a subset G randomly from the search space that contains the best simulated systems with high probability. To guarantee that this subset contains the best systems it needs to be relatively large. Then methods of ranking and selections will be applied to select a subset of k best systems of the subset G with high probability. The remaining systems of G will be replaced by newly selected alternatives from the search space. This procedure is repeated until the probability of correct selection (a subset of the best k simulated systems is selected) becomes very high. The optimal computing budget allocation is also used to allocate the available computing budget in a way that maximizes the probability of correct selection. Numerical experiments for comparing these algorithms are presented. r 2006 The Franklin Institute. Published by Elsevier Ltd. All rights reserved. Keywords: Statistical selection; Simulation optimization; Ordinal optimization; Ranking and selection
1. Introduction We consider the following optimization problem min JðyÞ, y2Y
(1)
where Y, the feasible solution set, is an arbitrary, has no structure, finite but huge set, and J is the expected performance measure of some complex stochastic system, therefore, J can Corresponding author. Tel.: +962 2 7201000; fax: +962 2 7095014.
E-mail addresses: [email protected] (M.H. Alrefaei), [email protected] (M. Almomani). 0016-0032/$30.00 r 2006 The Franklin Institute. Published by Elsevier Ltd. All rights reserved. doi:10.1016/j.jfranklin.2006.02.020
ARTICLE IN PRESS 496
M.H. Alrefaei, M. Almomani / Journal of the Franklin Institute 344 (2007) 495–506
be represented as JðyÞ ¼ E½Lðy; xÞ, where y is a vector representing the system design parameters and x represents all the random effects of the system and L is a deterministic function that depends on x and y. When the search space Y is small, then ranking and selection (R&S) procedures can be used for ranking the systems and select a subset that contains the best s systems with a prespecified significance level, see for example, Rinott [1] for a two stage R&S procedure. R&S cannot be applied for large-scale problems because it needs a huge computational time. To reduce the computational effort, the idea of ordinal optimization (OO) will be used, of course, the objectives will also be changed, instead of looking for the best design, which may be impossible in the case of our setting where Y is huge, we will focus on finding a good enough design. This will save the amount of time spent in simulating each design, see Ho et al. [2], more details about OO will be given in the next section. In this paper, we propose two sequential algorithms for selecting a subset of k systems that is contained in the set of the top s systems. Each iteration of these algorithms consists of two stages; in the first stage, the OO is used to select a subset of g elements randomly from Y. In the first algorithm, we use R&S procedures in the second stage for selecting a small subset from the subset that is selected in the first stage. In the second stage of the second algorithm, we use the optimal computing budget allocation (OCBA) technique that distribute the available computing budget on the designs in such a way to maximize the probability of correct selection (CS). The rest of the paper is organized as follows, in Section 2, we review the R&S procedures for selecting best systems among a small set of alternative systems and we review the idea of OO. In Section 3, we present the OO with the R&S algorithm. In Section 4, we present the second algorithm that uses the OCBA in the second stage. In Section 5, we give some numerical experiments and we make comparisons between the two proposed algorithms. Finally, in Section 6, we give some concluding remarks. 2. Background 2.1. Ranking and selection procedures R&S procedures are used to make comparisons between several alternative systems whose response is normally distributed when the number of competent systems does not exceed 20, see Bechhofer et al. [3] and Law and Kelton [4]. There are three types of R&S procedures. The first procedure can be used for selecting a design that is indifferent from the actual best design by less than a prespecified value d , called the indifference zone, with a prespecified significance level P , P is the probability that we make the CS. The second procedure can be used for selecting a subset of m systems from Y, so that this selected subset contains the best system with high probability. The third procedure is used for selecting a subset that contains the best k systems from Y with a prespecified probability, for more details on ranking and selections, see Bechhofer et al. [3] and Law and Kelton [4]. In this paper, we are interested in the third procedure that can be used for selecting the best k of n systems. For i ¼ 1; 2; . . . ; n, let xi1 ; xi2 ; . . . ; xini be a sample of ni i.i.d. observations from system i, let mi ¼ EðLðyi ; xij ÞÞ, and let mil be the lth smallest of mi ’s, so
ARTICLE IN PRESS M.H. Alrefaei, M. Almomani / Journal of the Franklin Institute 344 (2007) 495–506
497
that mi1 pmi2 p pmin . The CS, here, is to select k systems with the smallest expected means mi1 ; mi2 ; . . . ; mik . If mik and mikþ1 are close to each others, then we might not care on choosing the system whose mean is mikþ1 as one of the best k systems. We want to avoid making a large amount of simulation time to resolve this unimportant difference. More precisely, we are interested in selecting k systems that are indifferent from the best k systems by less than d ; i.e. fall in the region ½mi1 ; mik þ d , where the indifference d is prespecified by the user. We want the CS to take place with a probability greater than or equal to P , where P is prespecified by the user. Therefore, with probability P , the expected value of each of the k selected systems will not exceed mik þ d . 2.1.1. A two-stage R&S (Rinott’s) procedure The two-stage procedure of Rinott [1] is widely used in R&S. In the first stage, each design is sampled using n0 -independent simulation runs and a first stage sample mean j i ðn0 Þ and sample variance S 2i ðn0 Þ are computed for each design i ¼ 1; . . . ; n. Based on the initial simulation runs n0 , and the sample variance estimate obtained in the first stage, the number of additional simulation runs for each design in the second stage is computed by the formula 2 2 h S ðn0 Þ N i ¼ max n0 þ 1; 1 i 2 , (2) ðd Þ where due is the smallest integer greater than or equal u, d is the indifference zone, and h1 is a constant that solves the Rinott’s integral Z 1 F ðh þ tÞn1 f ðtÞ ¼ P , 1
where f and F denote the probability density and cumulative distribution functions, respectively, of the t distribution with n0 1 degree of freedom. This choice of h guarantee the PðCSÞ ¼ P . h can be obtained from tables of Wilcox [5] or Law and Kelton [4]. Then N i n0 additional simulation runs will be done in the second stage, and the overall sample mean j^i will be computed. The design that has the smallest j^i will be selected as the best one. This two-stage can be used also for selecting the best k of n designs. Chen and Kelton [6] consider reducing the value of N i by replacing the value of h in Eq. (2) by hi for each alternative i, where hi ¼ h=ri and ri ¼ maxfmi mi1 ; d g=d . R&S requires a huge amount of time to guarantee the CS with high probability. This huge amount of computational time restrict the size of the problems that these approaches can solve, in fact they can solve a small size problems, say less than 20. Nelson et al. [7] and Alrefaei and Alawneh [8] propose using subset selection in the first stage to reduce the competent designs. These designs are carried out to the second stage in which the indifference zone approach is used where the information in the first stage are used. 2.2. Ordinal optimization procedure The OO approach has been proposed by Ho et al. [2]. The aim of this approach is to find good, better, or best systems, i.e. ‘‘good enough’’ solution, rather than estimating the performance value of these systems accurately. In many real world problems, one could compromise for the ‘‘good enough’’ rather than the best alternative, especially when the
ARTICLE IN PRESS 498
M.H. Alrefaei, M. Almomani / Journal of the Franklin Institute 344 (2007) 495–506
number of competent alternative systems is very large and it is impossible to simulate all alternatives. The main idea of the OO approach is based on two issues: (1) It is much easier to determine ‘‘order’’ than ‘‘cardinal’’. This is intuitively reasonable, to determine whether A is greater or less than B, is a simpler task than to determine the value of A B, in stochastic situations. (2) Softening the goal of optimization also makes the problem easier. Instead of asking the ‘‘best for sure’’ we settle for the ‘‘good enough with high probability’’. Suppose that the CS is to select a subset G of g systems from Y that contains at least one of the top m% best designs. Since we assume that Y is very huge then the probability of CS is given by m g PðCSÞ 1 1 . 100 Now, suppose that the CS is to select a subset G of g designs that contains at least r of the best s systems. Let S be the subset that contains the actual best s systems, then the probability of CS can be obtained using the hypergeometric distribution as ! ns s g gi X i ! PðCSÞ ¼ PðjG \ SjXrÞ ¼ . n i¼r g Since we assume that jYj ¼ n, and n is very large then PðCSÞ can be approximated by the binomial random variable. Therefore, g X g m i m gi PðCSÞ 1 , 100 100 i i¼r where we assume that s=n 100% ¼ m%. It is also clear that this PðCSÞ increases when the sample size g increases. In the literature, PðCSÞ is denoted as the alignment probability. 3. Ordinal optimization with ranking and selection method In this section, we present the proposed procedure for selecting a subset of k designs that belong to the actual best s designs. The method of selecting one design from the actual top m% designs is presented by Alrefaei and Abdul-Rahman [9]. We first present the method of Alrefaei and Abdul-Rahman for selecting one design from the best m% designs, then we present the proposed algorithm for selecting a subset that contains k designs from the best m% designs. 3.1. A sequential algorithm for selecting one design from the best m% designs This method consists of a sequential procedure. Each iteration of the method is based on a two stage procedure. In the first stage, randomly select a subset G, where jGj ¼ g from
ARTICLE IN PRESS M.H. Alrefaei, M. Almomani / Journal of the Franklin Institute 344 (2007) 495–506
499
the search space Y. In the second stage, R&S is applied for selecting the best design yb in the selected subset G with a specified significance level (the probability of selecting the best design) P X1 a. yb will be carried out to the next iteration and the other g 1 designs will be replaced by newly sampled designs from Y. The algorithm continue in this manner until a specified number of generated samples is exceeded or a desired significance level is achieved. Note that, in the first stage, the probability that the randomly selected subset G contains one of the best m% designs can be evaluated as follows: PðG contains one of the best m% designs) 1 ð1 m=100Þg . Let bi be the selected best design in iteration i of the algorithm, and let G0i be the selected subset in iteration i, where jG 0i j ¼ g 1. Let P ¼ 1 a be the significance level used in the R&S procedure (the second stage of the algorithm, i.e. the probability of selecting the best design in the second stage is P ¼ 1 a). Then a lower bound of the probability of selecting one of the top m% best designs in the first iteration is m g P1 ¼ ð1 aÞ 1 1 100 and a lower bound of the probability of correct selection in the jth iteration is given by Pj ¼ Pðbj1 is in the top m% best designs or G 0j contains one of the top m% best designs)ð1 aÞ. Which equals to Pj ¼ [Pðbj1 is in the top m% best designs) + PðG 0j contains one of the top m% best designs) Pðbj1 is in the top m% best designs) PðG 0j contains one of the top m% best designs)]ð1 aÞ. As an example, consider m% ¼ 10% and 1 a ¼ 0:95 then P1 ¼ 0:619, P2 ¼ 0:808, P3 ¼ 0:925. This means that in three iterations, we can select a design that belongs to the top 10% best designs with probability greater than or equal to 0.925 regardless of the size of the optimization problem given in Eq. (1). 3.2. Subset selection approach In each iteration of this approach, we also use two stages, in the first stage, we select a random subset G of g systems from Y, in the second stage, we use the R&S procedure to select a subset of the kp20 best designs of G with indifference zone d . We repeat this procedure until the probability of CS becomes large enough. The algorithm is described as follows: Algorithm 2 (OO þ R&S). Step 0: Select a subset G randomly from Y, such that jGj ¼ g, gp20. Step 1: Generate n0 i.i.d. samples xij , j ¼ 1; 2; . . . ; n0 , for all i ¼ 1; . . . ; g. Step 2: For i ¼ 1; . . . ; g, compute the first stage sample means J ð1Þ i ðn0 Þ ¼
n0 1X Lðyi ; xij Þ, n0 j¼1
and the sample variances S2i ðn0 Þ ¼
n0 1 X 2 ½Lðyi ; xij Þ J ð1Þ i ðn0 Þ n0 1 j¼1
ARTICLE IN PRESS M.H. Alrefaei, M. Almomani / Journal of the Franklin Institute 344 (2007) 495–506
500
and let
h2 S 2 ðn0 Þ N i ¼ max n0 þ 1; 3 i 2 ðd Þ
,
where due is the smallest integer greater than or equal u, d is the indifference zone, and h3 depends on g, k, n0 , and P , and can be obtained from Table 10.13 of Law and Kelton [4]. Step 3: For i ¼ 1; . . . ; g, generate N i n0 i.i.d. samples xij . Compute the overall sample means J^ i ðN i Þ as follows: Ni 1 X J^ i ðN i Þ ¼ Lðyi ; xij Þ. N i j¼1
Step 4: Select the best k systems from the set G, call this subset as G 0 . Step 5: If a stopping rule is not satisfied; randomly select a subset G 00 of g k alternatives from YnG 0 , let G ¼ G 0 [ G 00 , and Go to Step 1. Remark. Note that the R&S procedures can be applied only when the number of alternatives is small, say less than 20, this is why in Step 0, we select gp20. Let S be the subset that contains the actual s best systems, then the calculation of the probability of CS is very complicated. As an example if k ¼ 2, then a lower bound of the probability of CS in the first iteration, can be found as follows: !3 2 ns s 6X gi 7 i 6 g 7 7 P . ! P1 ðCSÞ ¼ PðjG \ SjX2Þ P ¼ 6 6 7 n 4 i¼2 5 g In the ith iteration, the probability of CS can be calculated as follows: Pi ðCSÞ ¼ PðjG \ SjX2Þ P ¼ P½jðG0 [ G 00 Þ \ SjX2 P ¼ P½ðjG0 \ SjX2 or jG 00 \ SjX2 or ðjG 0 \ SjX1 and jG 00 \ SjX1ÞÞ P . Let A be the event A ¼ jG 0 \ SjX2, B be the event B ¼ jG00 \ SjX2, and C be the event C ¼ ðjG 0 \ SjX1 and jG00 \ SjX1Þ then by assuming that the events A, B, and C are mutually independent, we get Pi ðCSÞ ¼ PðA [ B [ CÞ P ¼ ½PðAÞ þ PðBÞ þ PðCÞ PðAÞPðBÞ PðAÞPðCÞ PðBÞPðCÞ þ PðAÞPðBÞPðCÞ P . It is clear that PðAÞ ¼ Pi1 ðCSÞ,
ARTICLE IN PRESS M.H. Alrefaei, M. Almomani / Journal of the Franklin Institute 344 (2007) 495–506
501
! n2s s g2 g2i X i ! . PðBÞ ¼ n2 i¼2 g2 Note also that PðjG 0 \ SjX1ÞXPi1 ðjG0 \ SjX2Þ ¼ Pi1 ðCSÞ and PðjG 00 \ SjX1Þ ¼
1 1
s g2 , n2
so that PðCÞ ¼ ðPi1 ðCSÞÞ 1 1
s g2 . n2
It is clear that PðCSÞ increases as the number of iterations increases. As an example, if k ¼ 2, jYj ¼ n ¼ 1000 systems, s ¼ 100, and g ¼ 10. Then one can evaluate PðCSÞ in each iteration as follows. In the first iteration, note that ! ! ns s PðjG \ SjX2Þ ¼
g X
gi ! n
i
i¼2
100 ¼
10 X
g !
i
1000 100
!
10 i ! 1000
i¼2
10 ¼ 0:263703353. If we use the binomial approximation, this probability becomes PðjG \ SjX2Þ ¼ 0:26390. Let Pi ðCSÞ be the probability of selecting at least two systems that belong to the top 10% systems in iteration i, let P ¼ 0:95 in Step 2, then P1 ðCSÞ is given by ! !1 0 100 1000 100 BX C B 10 C i 10 i C ð0:95Þ ! P1 ðCSÞ ¼ B B C 1000 @ i¼2 A 10 ¼ ð0:263703353Þ ð0:95Þ, P1 ðCSÞ ¼ 0:25052. By simple calculations, we can find a lower bound of the probability of CS in the second iteration, P2 ðCSÞ ¼ 0:45387 and P3 ðCSÞ ¼ 0:63739. For a lower bound of Pi ðCSÞ,
ARTICLE IN PRESS M.H. Alrefaei, M. Almomani / Journal of the Franklin Institute 344 (2007) 495–506
502
Table 1 The lower bound of the probability of correct selection Pk ðCSÞ, k ¼ 1; . . . ; 10, when a subset of the best 10% systems is desired k
1
2
3
4
5
6
7
8
9
10
pk
0.251
0.454
0.637
0.772
0.851
0.891
0.909
0.916
0.919
0.920
i ¼ 1; . . . ; 10, see Table 1. The maximum value of the probability of CS for this choice of parameters does not exceed 0:9211, which can be reached at iteration 13. To increase this value, one needs to increase the value of P . 4. Ordinal optimization with optimal computing budget allocation method The OCBA approach is highly efficient, and significantly reduces the total simulation cost. This approach has been proposed by Chen [10], Chen et al. [11], and Chen et al. [12]. The goal of this approach is to allocate the total simulated samples from all systems in a way that maximizes the probability of selecting the best system within a given computing budget. Suppose that we are allowed to take a total of T samples from all systems for solving the optimization problem given in Eq. (1). The question is where to allocate these computing simulated samples to maximize the CS. Where the CS is to select a system b that has the best performance. Furthermore, instead of equally simulating all systems, we will further improve the performance of OO by determining the best numbers of simulated samples for each system. Suppose that Y ¼ fy1 ; . . . ; yn g, then we wish to choose T 1 ; T 2 ; . . . ; T n such that PðCSÞ is maximized, subject to a limited computing budget T. max
T 1 ;...;T n
s:t:
PðCSÞ n X
T i ¼ T,
i¼1
Pn where T i 2 N for i ¼ 1; 2; . . . ; n and N is the set of non-negative integers and i¼1 T i denotes the total computational samples, assuming the simulation times for different systems are roughly the same. The proof of the following theorem has been given by Chen et al. [12]. Theorem 1. Given a total number of simulated samples T to be allocated to n competing systems whose performance is depicted by random variables with means Jðy1 Þ; Jðy2 Þ; . . . ; Jðyn Þ, and finite variances s21 ; s22 ; . . . ; s2n , respectively. As T!1, the approximate probability of CS can be asymptotically maximized when (1) T i =T j ¼ ððsi =db;i Þ=ðsj =db;j ÞÞ2 , i, j 2 f1; 2; . . . ; ng, and iajab. qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi Pn 2 2 (2) T b ¼ sb i¼1;iab T i =si , where T i is the number of samples allocated to system i, db;i the estimated difference between the performance P i of the two designs ðdb;i ¼ J b J i Þ, and J b pmini J i for all i. Here J i ¼ ð1=T i Þ Tj¼1 Lðyi ; xij Þ, where xij is a sample from xi for j ¼ 1; . . . ; T i .
ARTICLE IN PRESS M.H. Alrefaei, M. Almomani / Journal of the Franklin Institute 344 (2007) 495–506
503
Now, we present the second algorithm in which we use the OCBA in the second stage for selecting a subset of the k best systems in G instead of R&S. Algorithm (OO þ OCBA). Step 0: Determine the size of the subset G, gXk. Choose the number of initial simulations samples t1 X2 for each system. Let l ¼ 0 and let T l1 ¼ T l2 ¼ ¼ T lg ¼ t1 . Step 1: Select a subset G of g alternatives randomly from Y. Take random samples of t1 observations xij ð1pjpt1 Þ for each system i in G, where i ¼ 1; . . . ; g. Step 2: Calculate the first stage sample mean J ð1Þ i for all i ¼ 1; . . . ; g Tl
l J ð1Þ i ðT i Þ
i 1 X ¼ l Lðyi ; xij Þ. T i j¼1
Step 3: Order the systems in G according to their sample averages ð1Þ l l ð1Þ l J ð1Þ 1 ðT 1 ÞpJ 2 ðT 2 Þp pJ g ðT g Þ. 0 Step 4: Select Pg thel best k systems from the set G, call this subset as G . Step 5: If i¼1 T i XT, stop. Step 6: Randomly select a subset G 00 of g k alternatives from Y G 0 , let G ¼ G 0 [ G 00 . Step 7: Increase the computing budget by D and compute the new budget allocation, lþ1 lþ1 T lþ1 1 ; T 2 ; . . . ; T g , according to Theorem 1. Step 8: Perform additional maxf0; T lþ1 T li g samples for each system i, i ¼ 1; . . . ; g, let i l l þ 1. Go to Step 2.
5. Numerical example Suppose we have 10,000 different systems Y ¼ y1 ; . . . ; y10;000 , where Jðyi Þ ¼ i=100 and Lðyi ; xÞ is the normal distribution NðJðyi Þ; 1Þ, with mean m ¼ Jðyi Þ and variance s2 ¼ 1. The optimization problem is min
i¼1;...;10;000
Jðyi Þ.
We want to select a subset that is contained in the best 10% systems, it is obvious that y1 ; y2 ; . . . ; y1000 are the actual best 10% systems. We apply our algorithms where we use R&S in each iteration with significance level P ¼ 0:95, the initial sample size t1 ¼ 40, and the indifference zone d ¼ 0:3. The maximum computing budget in the second algorithm is T ¼ 1500, D ¼ 100, and t1 ¼ 40. In both algorithms, we randomly select a set that consists of 10 systems in the first iteration, and select a subset G 0 that contains the best two systems among these 10 systems using the two proposed algorithms. Then we replace the other eight systems by eight new systems randomly selected from the Y G0 . We repeat this Table 2 The probability of selecting two alternatives from the top 10% of n ¼ 10; 000 systems using g ¼ 10, t1 ¼ 40 Iteration
1
2
3
4
5
6
7
8
9
10
OO þ R&S OO þ OCBA
0.23 0.26
0.55 0.6
0.74 0.8
0.88 0.88
0.95 0.94
0.96 0.99
0.99 0.99
0.99 0.99
1.0 0.99
1.0 1.0
ARTICLE IN PRESS 504
M.H. Alrefaei, M. Almomani / Journal of the Franklin Institute 344 (2007) 495–506
Table 3 The probability of selecting five alternatives from the top 10% of n ¼ 10; 000 systems using g ¼ 10, t1 ¼ 40 Iteration
2
4
6
8
10
12
14
15
OO þ R&S OO þ OCBA
0.0 0.0
0.04 0.13
0.21 0.37
0.44 0.55
0.56 0.72
0.72 0.8
0.88 0.88
0.92 0.92
Table 4 The probability of selecting two alternatives from the top 1% of n ¼ 10; 000 systems using g ¼ 10, t1 ¼ 40 Iteration
1
5
10
15
20
25
30
40
50
OO þ R&S OO þ OCBA
0.01 0.0
0.06 0.07
0.2 0.25
0.32 0.43
0.48 0.52
0.66 0.59
0.78 0.69
0.87 0.82
0.92 0.93
experiment 100 replications and in each replication we use only 10 iterations. Table 2 contains the probability of CS obtained after 10 iterations over the 100 replications by using the R&S and OCBA. It is clear that the OCBA gives better results in almost all iterations. In all cases, we use almost the same computational time for both algorithms. Now, suppose that we want to select a subset that contains five alternatives from the top 10% alternative systems using the same parameters in the previous setting. Table 3 includes the results for this experiment. Again the OO þ OCBA algorithm gives better results. The OO þ R&S also behaves well, it starts slow, but it moves faster in later iterations. We repeat the two algorithms with the objective of selecting a set that contains two systems among the best 1% systems, it is obvious that y1 ; y2 ; . . . ; y100 are the actual best 1% systems. For the ðOO þ R&SÞ algorithm, we use a significance level P ¼ 0:95, the initial sample size t1 ¼ 40, and the indifference zone d ¼ 0:6. The maximum computing budget for the ðOO þ OCBAÞ is T ¼ 1500, D ¼ 100, and t1 ¼ 40. In both algorithms, we select a subset that consists of 10 systems in the first iteration and select a subset G0 that contains the best two systems among these 10 systems, and then replace the other eight systems by eight new systems randomly selected from the search space Y G0 . We run this experiment for 100 replications and in each replication we use 50 iterations. Table 4 contains the probability of CS obtained after 50 iterations over the 100 replications. It is clear that both of the algorithms quickly locate the desired set and that the OO þ OCBA starts faster and then the OO þ R&S accelerates later. To see the robustness of the proposed algorithms, we have applied these algorithms on the following example. Suppose, we have nðM=M=1Þ queuing systems and we want to select a subset of k, where k ¼ 1; . . . ; 5 systems that is contained in the best 10% systems, where the best system is the system that has the minimum average waiting time per customer. The arrival rate is assumed to be a fixed number, l ¼ 1, and the service rate is m 2 Y ¼ f4:0; 4:0001; 4:0002; . . . ; 5:0g, so we have 10,000 different ðM=M=1Þ queuing systems. In the first stage, we use the OO to select a set G of 10 systems and in the second stage, we use the R&S procedures to select the best k systems from the set G with indifference zone d ¼ 0:05. This experiment is repeated 1000 replications, each replication
ARTICLE IN PRESS M.H. Alrefaei, M. Almomani / Journal of the Franklin Institute 344 (2007) 495–506
505
Fig. 1. The performance of the OO þ R&S algorithm when they are applied to select a subset of k, k ¼ 1; . . . ; 5 systems from the top 10% M/M/1 queuing systems.
Fig. 2. A comparison between the OO þ R&S and the OO þ OCBA algorithms when they are applied to select a subset of 2 of the best 10% M/M/1 queuing systems.
ARTICLE IN PRESS 506
M.H. Alrefaei, M. Almomani / Journal of the Franklin Institute 344 (2007) 495–506
consists of 20 iterations; the probability of CS is shown in Fig. 1. It is clear from the figure that the probability of CS increases very rapidly. To compare the OO þ R&S algorithm with the OO þ OCBA algorithm, we have implemented the OO þ OCBA algorithm to select a subset of two systems that belongs to the best 10% systems. We use initial sample size of t1 ¼ 10 samples only in this experiment. Fig. 2 shows the relation between the probability of CS when we use the OO þ R&S algorithm and the OO þ OCBA for 1000 replications for the first 20 iteration. It is clear that the OCBA gives better results with computational time approximately the same as in R&S. 6. Conclusions We have presented two sequential algorithms for selecting a subset of k systems from the top s systems when the number of alternatives is very large. Each iteration of the algorithm consists of two stages. In the first stage, we randomly select a subset G of g alternatives. OO is used for evaluating the probability that this set is contained in the best s alternatives set. In the second stage, we use one of the two procedures, in the first procedure we use the R&S procedure for ranking the alternative systems in the selected subset in the first stage and then select a subset that contains the best k systems of this subset. The other g k alternatives are replaced by newly g k alternatives selected randomly from the search space and the algorithm is repeated until a stopping criterion is reached. In the second procedure, we use the OCBA, in which the number of samples are distributed over the competent alternatives in a way that maximizes the probability of CS. We apply these algorithms on two numerical examples. From our numerical results, we note that these algorithms achieve the CS with high probability, quickly. References [1] Y. Rinott, On two-stage selection procedures and related probability–inequalities, Commun. Statist.: Theory Methods A 7 (1978) 799–811. [2] Y.C. Ho, R.S. Sreenivas, P. Vakili, Ordinal optimization of DEDS, J. Discrete Event Dynamic Syst. 2 (1992) 61–88. [3] R.E. Bechhofer, T.J. Santner, D.M. Goldsman, Design and Analysis of Experiments for Statistical Selection, Screening, and Multiple Comparisons, Wiley, New York, 1995. [4] A.M. Law, W.D. Kelton, Simulation Modeling and Analysis, third ed., McGraw-Hill, New York, 2000. [5] R.R. Wilcox, A table for Rinott’s selection procedure, J. Quality Technol. 16 (2) (1984) 97–100. [6] E.J. Chen, W.D. Kelton, An enhanced two-stage selection procedure, in: J.A. Joines, R.R. Barton, K. Kang, P.A. Fishwick (Eds.), Proceedings of the 2000 Winter Simulation Conference, Institute of Electrical and Electronics Engineers, Piscataway, New Jersey, 2000, pp. 727–735. [7] B.L. Nelson, J. Swann, D. Goldsman, W. Song, Simple procedures for selecting the best simulated system when the number of alternatives is large, Oper. Res. 49 (2001) 950–963. [8] M.H. Alrefaei, A.J. Alawneh, Selecting the best stochastic system for large scale problems in DEDS, Math. Comput. Simulation 64 (2004) 237–245. [9] M.H. Alrefaei, H. Abdul-Rahman, Two sequential algorithms for selecting one of the best simulated systems, WSEAS Trans. Syst. 3 (2004) 2517–2522. [10] C.H. Chen, An effective approach to smartly allocated computing budget for discrete event simulation, IEEE Conference on Decision and Control, vol. 34, 1995, pp. 2598–2605. [11] C.H. Chen, C.D. Wu, L. Dai, Ordinal comparison of heuristic algorithms using stochastic optimization, IEEE Trans. Robotics Autom. 15 (1) (1999) 44–56. [12] C.H. Chen, E. Yucesan, S.E. Chick, Simulation budget allocation for further enhancing the efficiency of ordinal optimization, J. Discrete Event Dynamic Syst.: Theory Appl. 10 (2000) 251–270.