Research Rutcor Report Large Scale LP Model for Finding Optimal Container Inspection Strategies
E. Borosa L. Fedzhorab P. B. Kantorc K. Saegerd P. Stroude
RRR 26-2006, October, 2006
RUTCOR Rutgers Center for Operations Research Rutgers University 640 Bartholomew Road Piscataway, New Jersey 08854-8003 Telephone: Telefax: Email: 732-445-3804 732-445-5472 [e-mail]
Rutgers University, 640 Bartholomew Road, Piscataway, NJ 08854-8003; [e-mail] b RUTCOR, Rutgers University, 640 Bartholomew Road, Piscataway, NJ 08854-8003; [e-mail] c SCILS, Rutgers University; [e-mail] d LANL; [e-mail] e LANL; [e-mail]
a RUTCOR,
http://rutcor.rutgers.edu/rrr
Rutcor Research Report
RRR 26-2006, October, 2006
Large Scale LP Model for Finding Optimal Container Inspection Strategies
E. Boros
L. Fedzhora
P. B. Kantor
K. Saeger
P. Stroud
Abstract. Cargo ships arriving at USA ports are inspected for unauthorized materials. Ideally, we want to open and check every container they carry, but it is costly and time-consuming. Instead, tests are developed to decide whether a container should be opened. By utilizing a polyhedral description of decision trees, we develop a large scale linear programming model for sequential container inspection that determines an optimal inspection strategy under various limitations, improving on earlier approaches in several ways: (a) we consider mixed strategies and multiple thresholds for each sensor, which provide more eective inspection strategies; (b) our model can accommodate realistic limitations (budget, sensor capacity, time limits, etc.), as well as multiple container types; (c) our model is computationally more tractable allowing us to solve cases that were prohibitive in preceding models, and making it possible to analyze the impact of potentially new sensor technologies.
Acknowledgements: The authors are thankful for the partial support by the National Science Foundation (Grant NSFSES 05-18543) and the Oce of Naval Research (Grant N00014-05-1-0237).
RRR 26-2006
Page 1
1
Introduction
We consider the problem of deploying a set of sensors most eectively when the goal is to detect as many as possible of some extremely rare contraband, which we refer to here as bad contents. The sensors have operating costs (part of which is proportional to the number of containers inspected, and part of which is independent of the inspected containers), and they have individual capacity limitations (which can be extended in discrete steps at certain investment costs). There is also a manual checking procedure, CHK, which is considered 100% eective, and considerably more expensive than the application of other sensor technologies. As is widely known (see e.g., newspaper articles) security at large US ports does not even consider applying sensors to all containers. Rather, most containers are cleared on the basis of other examinations (e.g., some sort of inspection done at the port of origination). Hence we expect the methods presented here to be applied only to the 5-10 percent of all containers that are judged to need further inspection. Whether the specic sensors are hand-held and portable, or xed installations through which the containers must be moved, a container may be routed through the sensors according to the readings it has produced. This can be conceptualized by saying that a particular design of the inspection schedule produces a tree. The rst sensor to be applied is at the root node. According to the results on that sensor, the container may be released, checked manually, or sent to another sensor, etc. Stroud and Saeger [4] have formulated the problem in precisely this way. They postulate that each sensor makes a binary decision, which can be thought of as more suspicious or less suspicious, and they note that the number of possible binary trees is doubly exponential in the number of sensors. They also note that the threshold dividing more suspicious from less suspicious need not be xed in advance, but could be adjusted after the topology of the tree has been set. Any given sensor may appear in dierent parts of the tree, with, in principle, dierent threshold settings. They examined the problem, using the same settings for each sensor, wherever it appears. The problem is still a non-linear optimization problem. In that model, costs were assigned to false negatives and false positives, and the combined cost of sensor operations, plus the expected cost of the outcome is minimized by adjusting the thresholds for each tree. We propose here a number of changes in viewpoint, which, taken together, lead to a more tractable problem, with greater exibility to incorporate additional constraints/objectives, greater potential for generalization to more sensors, and with some attractive operational advantages. The changes in perspective are the following: First, we allow that when a sensor is used in dierent locations in the tree, it may have dierent thresholds. Second, we allow the use of more than one threshold at a given sensor (in a given location in a tree) so that the branching following the application of a sensor may be, in general manifold. We refer to this discretization of sensor readings as the partition of the outcomes. Third, we view the inspection process not with a birds- eye view which sees the entire tree, but with a path- wise view, which considers the possible histories that any container might experience, under a given overall set of partitions of the readings for each of the
Page 2
RRR 26-2006
sensors. Fourth, we reformulate the problem into a comparable form which is more robust in the presence of unknown costs. Specically, we seek to maximize the detection rate, given constraints on operating budget, capacities, time, etc. We feel that this formulation will be more useful to policy-makers, since they can examine the increase in detection as a function of the increase in budget, and use their own cost coecients to complete the decision process. Fifth, we recognize that in the presence of even a single (budgetary, for example) constraint, the optimal detection rate is not necessarily achieved by a pure strategy. When there is a single linear constraint, and a convex problem, the optimum will involve a mixed strategy using at most two pure strategies. As more constraints are added the number of components in the optimal mix rises. Let us remark here that any given decision tree has some vulnerabilities. Namely, assuming intelligent adversaries, the packaging and timing of bad content may be adapted to the particular inspection sequence, so that to minimize detection. A mixed strategy is much less vulnerable in this sense, since our adversaries can never be sure in advance about the type of inspection strategy that will be applied1 . Somewhat surprisingly, complexifying the problem in several ways (moving from trees to paths, considering a larger number of static thresholds rather than a small number of variable ones, seeking mixture solutions directly, etc.) makes the problem of sensor sequencing more tractable, rather than less2 . First of all, the pathwise view of the problem decreases the number of decision variables by an order of magnitude. Secondly, in moving from a problem with adjustable thresholds at each sensor, to a problem with xed partitions at each sensor, we have replaced a really huge number (the number of trees) of small non-linear problems (nding the thresholds) by a single linear optimization problem with a very large (the number of paths) variables, which is still tractable by standard o-the-shelf software packages. We illustrate the approach with examination of the model discussed in [4], and discuss the eects of choosing dierent partitions. To sum up: we formulate a new problem in which each sensor produces a reading which is quantized into one of many bins; we describe the problem not in terms of trees, but in terms of paths, which a container may follow through the tree; we exhibit the need for mixed strategies and develop a tractable linear programming formulation. The mathematical basis for our model is a polyhedral description of all decision trees in the space of paths. We are able to show that this method, which retains and makes use of more specic information gathered at earlier sensors, is able to provide better results at any given budget. We briey examine the eects of rening the partitions on both computational cost and overall quality of the solution found.
Of course, bribing the harbor master remains an ancient method of proven eectiveness. I cant cure your cold, but if you get pneumonia, I can help you. (Joke dating to the discovery of penicillin.)
2 1
RRR 26-2006
Page 3
2
Detailed Discussion and Background
An inspection strategy involves the sequential application of sensors, where the choice of a next sensor depends on the results obtained in the previous steps. Such an inspection sequence terminates with a nal decision about the container. In our study we assume two possible nal decisions, either by declaring a container safe to leave the port and enter the country (OK), or by subjecting it to a manual inspection (CHK). Such inspection strategies customarily are represented as so called decision trees (see Section 4 for details). Even if we allow only two dierent decision to be made after each sensor reading (in other words, if the results of an application of a sensor translated into two classes, e.g., suspicious and not-suspicious), the number of possible (labeled binary) n decision trees is 22 , doubly exponential in the number sensors. Even the number of those decision trees which represent monotone Boolean functions is doubly exponentially large. A complete enumeration of the possible trees, in order to select a best one is a computational challenge even for a small number of sensors (see [4]). The problem of scheduling sensors, and of quantizing the information obtained at the sensors has a long history. Early development was driven by aerospace applications much of it related to the problem of missile detection (see e.g., Varshney [1, 10], Kushner and Pacut [9], Cherikh [2]). Related problems also arise in the use of medical tests (see e.g., Greiner [7]). We have not found any authors who treated the problem using the approach described here. However, in those papers, particularly in Kushner and Pacut [9] attention is drawn to the fact that what is essential for solving problems of this type is to represent each sensor in terms of its ROC performance curve (Radar Operational Curve). Our approach utilizes information about each sensor, which implicitly includes its ROC curve, and so in eect we use the ROC curve, though not in an explicit way. In our study we use several operational characteristics to describe the eectiveness of a particular inspection strategy, including the detection rate (the probability that a truly dangerous container content will be discovered), inspection cost (per container average inspection cost), sensor load (expected fraction of containers inspected by a particular sensor; related to capacity planning/limitations), inspection time (expected time to inspect a container), etc. Ideally, we would like to check manually every container entering the country. However, budget, capacity and time delay limits make this impossible. The application of sensors, which typically can be applied faster and more cheaply than manual inspection, allows for ltering only the truly suspicious containers, which then can feasibly be checked manually. However, since no sensor is perfect, and each provide only a limited amount of information about the content of the inspected containers, such a non-ideal inspection strategy does not guarantee the detection of all dangerous contents. It is therefore important to understand how high a detection rate one can achieve from a given inspection budget, how this achievable maximum detection rate changes if we change the budget, how much sensor capacity we need, how much delay this causes in the container ow, etc. Furthermore, for planning purposes one must nd the best way to build/expand
Page 4
RRR 26-2006
existing inspection stations, study the eects on eciency and cost of the introduction of possible new sensor technologies, etc. To address these and similar questions, we propose a large scale linear programming model, which provides a computationally ecient tool for analysis. The theoretical basis for this model is a polyhedral description of the set of decision trees, in a high dimensional space. Based on this description we show that the above questions can be answered eciently by linear programming techniques. Let us point out two particularly important outcomes, in our view, of our study. The rst is the realization that a best inspection strategy is not simply a decision tree, but rather an ensemble of decision trees, each applied to a certain fraction of the containers (which decision trees and what fraction will be determined by the optimization model). Note that such a mixed strategy is also much more secure, in the sense that it is much harder to adapt to it by packing methods that trick the inspection strategy. The second outcome is perhaps more natural. We show that the more information we retain from a sensor reading, the higher the eciency we can achieve. In other words, using the same set of sensors, the same capacities, and the same inspection budget, one can achieve a higher detection rate (and sometimes substantially higher) by allowing dierent routing decisions based on the actual readings of previous sensors, rather than coding every sensor reading into a binary decision (e.g., suspicious and non-suspicious). Many problems remain open, however. One of these is to decide how to dene the partitions. Since we are, in eect, making a continuous problem discrete, good choice of the partitions will move the polyhedral description closer to the underlying optimum. There is also the possibility that the mixture problem, for a relatively small number of sensors, can be solved directly for a relatively small number of cases. The solution path here is to seek the optimal decision, at later sensors, based on the precise values generated at the earlier sensors. If, as we are assuming in this paper, the variation in sensors is distributed independently, then for N sensors there are not too many dierent sequences. In general, there are N !/k! sequences involving N k of the sensors, and thus in all there are at most eN ! basic sequences. However, the determination of the cost- performance curve for each of these is a subtle problem in non-linear optimization, which will be discussed elsewhere.
2.1
Sensor Model
Currently available sensors provide only partial and noisy information about the contents of a container. Typical output is either a real number, or an image, which then is translated into a real scale by a human operator. The usual physical sensor model therefore involves a distribution function, the distribution of sensor readings on a real scale. Of course, for dierent container types (dierent container content, container material, container size, etc.) we may have dierent distributions corresponding to each sensor. We assume that these distributions are known.
RRR 26-2006
Page 5
For simplicity of presentation, we assume that there are only two types of containers arriving to the port, the good ones containing legal content, and the bad ones containing dangerous illegal content. This assumption is not necessary for our model, but makes the description simpler. Let us denote the set of reals by R (scale of sensor readings), the set of dierent sensors by S (typically |S| is small, say |S| = 4), the set of container types by T (we could consider several dierent types of containers; following [4], in this study we consider only two types T = {good, bad}), and let fs, (x) denote the density function of sensor readings on sensor s S for containers of type T . For a container C and sensor s S let us denote by Rs (C) R the result obtained by applying sensor s to container C. That is Rs (C) is the reading of sensor s on container C. We assume, following [4] that the same reading will be obtained if the container is remeasured. The realistic background for this assumption is that the primary source of randomness in practice is the large variety of container materials, contents, and varying backgrounds, rather than the uncertainties of the actual measurement. We also assume, following [4] that for each sensor and each container type we know the distribution of sensor readings, that is the density functions fs, for s S and T . Given a subset X R of the reals, we denote by (X, s, ) =
X
fs, (x)dx
(1)
the expected fraction of containers of type T which will have a reading in the given range X on sensor s S. In other words, (X, s, ) is the probability that Rs (C) X when C is a randomly selected container of type . Given any partition X1 X2 Xk = R of the set of reals, we clearly have
k
(Xi , s, ) = 1
i=1
for every sensor s S and container type T . We consider the proportions (Xi , s, ) for certain partitions Ps = {X1 , , Xk(s) }, for all sensors s S and container types T as inputs for our model. Let us add that the distribution of readings on a sensor s may also depend on the results of previous sensor readings, since due to those readings, the operator of sensor s may evaluate/interpret the measurement dierently. In other words, the posterior probability that the container is of type may depend on the earlier readings. We might include this in our model, however such a dependence needs more precise data about the sensors. The data available to us does not include such ne detail, and as a consequence, in our computational study we do not consider such dependence. To clarify our terminology, let us next consider, as a small example, a sensor a and the distribution of good and bad container readings as shown in Figure 1.
Page 6
frequency
RRR 26-2006 good
bad
t1
t2
reading on sensor a
Figure 1: An example for a sensor a and the distributions of sensor readings for good and bad containers. Traditionally, sensor readings are discretized by choosing appropriate threshold values. For instance, by choosing two thresholds t1 and t2 , we can partition the set of reals into three sets X1 = {x R|x t1 }, X2 = {x R|t1 < x t2 } and X3 = {x R|t2 < x}. Then, we might get (Xj , a, ) = good = bad j=1 0.50 0.10 0.60 j=2 0.40 j=3 0.10 0.30 To base some decisions on these readings, we can observe that the odds of having a bad container with a reading on sensor a within range X1 are very low, while they are relatively high in X2 and even higher in X3 . This ltering eect can be exploited by applying dierent sensors for further testing in an inspection strategy, depending on the reading on sensor a. Of course, the quality of our decision will strongly depend on which partition we consider and how ne it is.
2.2
Pure Inspection Strategies
Let us consider rst a small example, involving three sensors S = {a, b, c}, and the decision tree D shown in Figure 2. This decision tree may represent an inspection strategy, in which we start testing all containers by sensor a. Those containers C for which Ra (C) is very low are deemed to be safe, and let go without any further testing, while those with medium or very high readings are deemed suspicious, and tested further. For instance, those with very high reading on sensor a are tested further by sensor b, and if they get a high reading again,
RRR 26-2006
Page 7
a
very low medium very high
OK
very low
c
medium very high low
b
high
OK
low
b
high
CHK
low
c
high
CHK CHK
OK
CHK
OK
Figure 2: A decision tree utilizing three sensors S = {a, b, c}.
then checked manually (CHK ), while those with a low reading on sensor b are tested even further by sensor c, etc. Let us note that in this example, sensors c and b are used in two dierent ways. For instance, the readings on sensor c of those containers which have only a medium previous reading on sensor a are categorized as very low , medium or very high, while for those containers which had a very high reading on sensor a and a low reading on sensor b are categorized simply as low or high. Furthermore, though the readings on sensor b are always grouped into low and high, e.g., a low on b may mean two dierent things, depending on the previous readings of the containers, as the topology of D shown in Figure 2 suggests. Summarizing the above example, we may need seven dierent thresholds to realize the above decision tree D, two for sensor a, two possibly dierent ones for the two roles of sensor b, and three for the two dierent roles of sensor c. Denoting by t the vector of these thresholds, we view the pair (D, t) as a pure inspection strategy. In general, once we describe the topology of a decision tree D and select the necessary thresholds t, the pair (D, t) describes precisely how to test containers. Furthermore, using the sensor distributions, as we described in the previous subsection, for any pure strategy (D, t) we can compute the expected fractions of containers of type T which end up in a leaf node labeled CHK , i.e., which are checked manually. Thus, to any pure inspection strategy (D, t) we can associate several operational characteristics, computed from these fractions and other input data. For instance, we denote by (D, t) the detection rate of (D, t), which is the expected fraction of bad containers checked manually, and by (D, t) the rate of false positives, which is the expected fraction of good containers checked manually (unnecessarily). We can also compute the expected per-container-inspection cost C(D, t), the expected per-container-inspection time T (D, t), and expected sensor load s Us (D, t) for s S. Let us remark that we assume the proportion of bad containers to be negligible. Consequently, we can base, without any serious loss of generality, the computation of these
Page 8
RRR 26-2006
operational characteristics only on the ow of good containers (with the exception of the detection rate, which is computed solely from the bad container ow). a
(10, 0) (60, 40) (30, 60)
OK
(30, 4)
c
(25, 6) (5, 30) (24, 30)
b
(6, 30)
OK
(20, 3)
b
(5, 3)
CHK
(20, 10)
c
(4, 20)
CHK CHK
OK
CHK
OK
Figure 3: A pure inspection strategy utilizing three sensors S = {a, b, c}, and a particular threshold selection. Along each edge we indicate the percentages of good and bad containers that receive the corresponding readings. These notions will be dened precisely, later in the paper. Here we illustrate them by using a numerical example. Assume that we choose the thresholds t for the decision tree D shown in Figure 2, such that the percentages of good and bad containers receiving particular readings are as shown in Figure 3. Then we get the following values: (D, t) (D, t) Ua (D, t) Ub (D, t) Uc (D, t) = 0.30 + 0.03 + 0.20 + 0.30 = 0.83 = 0.05 + 0.05 + 0.04 + 0.06 = 0.20 = 1.00 = 0.25 + 0.30 = 0.55 = 0.60 + 0.24 = 0.84
Based on these we can compute C(D, t) and T (D, t) as follows C(D, t) =
sS
Cs Us (D, t) + CCHK (D, t)
= Ca + 0.55Cb + 0.84Cc + 0.20CCHK T (D, t) =
sS
Ts Us (D, t) + TCHK (D, t)
= Ta + 0.55Tb + 0.84Tc + 0.20TCHK where Cs and Ts respectively are the average inspection cost and inspection time of a single container on sensor s S, and CCHK and TCHK respectively are the cost and time of inspecting manually a container.
RRR 26-2006
Page 9
2.3
Finding the Best Pure Inspection Strategy
Using the notation introduced above, we can formulate the problem of nding a best pure inspection strategy as follows: min C(D, t) =
sS
Cs Us (D, t) + CCHK (D, t)
subject to (D, t) 0 where the optimization runs over all possible decision trees D, and all possible threshold selections t for D. This model was considered by [4] in the simplifying case of binary sensor outputs, assuming that 1) if a particular combination of container observations is selected by for CHK, then any combination of observations in which each observation is equally or more suspicious would also be selected for CHK, and 2) whether to CHK or OK depends on all of the sensors. For the case of four sensors observing four independent attributes, there turn out to be 114 feasible selection expressions. These can be implemented by 11,808 distinct decision trees. For each of these trees a nonlinear optimization problem was solved to obtain the best threshold selection by a grid-search method, and the best decision tree was identied. A test case was created assuming standard normal distribution of readings for bad containers, and normal distribution N (Ks , 1) of readings for good containers at sensors s = 1, 2, 3, 4, where the rst sensor can, however, recognize only half of bad containers. Discriminating power (the number of standard deviations separating good and bad containers on that sensor) and cost factors (notionally in dollars per scan) of each of four sensors represent generic sensor types. Sensor 1 has good discriminating power (for those containers that it can recognize at all) and is inexpensive, sensor 2 has low discriminating power but is inexpensive, sensor 3 has medium discriminating power and intermediate cost, and sensor 4 is powerful and expensive. The test case parameters are as follows: K1 K2 K3 K4 = = = = 4.37, 1.53, 2.9, 4.6, C1 C2 C3 C4 = = = = 0.32 0.92 57 176
with CCHK = 600 and 0 = .815. For each of the 11,806 feasible BDTs, a solver is used to nd the set of four threshold settings {t1 , t2 , t3 , t4 } that minimize the combined cost of observations and false positives for a specied minimum detection probability. Ten of the best thirteen BDTs implement the selection logic that suggests a manual inspection of those containers for which either A) sensor 2, sensor 3, and sensor 4 read high; or B) sensor 2, sensor 3, and sensor 1 read high; or C) sensor 1 and sensor 4 read high. There are 148 distinct feasible BDTs that implement this selection logic. The best 7 give
Page 10
RRR 26-2006
a cost of $13.56 per container and the following optimized threshold settings and utilization associated with each sensor: t1 t2 t3 t4 = = = = 2.6, 1.07, 1.16, 2.2, U1 U2 U3 U4 = = = = 1 1 0.14 0.02
Figure 4 describes the optimal pure strategy, where containers with the reading lower than a threshold follow the left branch of a node. The remaining six out of seven best BDTs can be found in the Appendix (see Figures 15-20).
2 1 OK 4 OK CHK OK 3 4 4 1 3 CHK
OK CHK OK CHK Figure 4: A decision tree representing optimal pure strategy for four sensors with one threshold each.
The worst 14 BDT congurations implementing this logic are so poor that they give a cost per container in excess of $148 (these poor BDTs place sensor 3 or 4 at node 1). Many seemingly quite dierent BDTs were found to implement near optimal tness, although at very dierently-tuned threshold levels. From the results of [4] we can conclude also that the above model has several limitations. First of all, it is computationally very expensive, and it is doubtful that one could solve it for n 5 sensors. Secondly, it might be substantially better to allow more than one thresholds at each of the sensors, which makes the formulation even less tractable computationally. Last, but not least, using a pure strategy for container inspection may not be the best thing to do. We will argue in the subsequent sections that the use of mixed strategies (combinations of several pure strategies) not only can improve on detection rates, without increasing costs, times, etc., but is also less vulnerable to adaptation of packaging, etc.
RRR 26-2006
Page 11
2.4
The Advantage of Mixed Inspection Strategies
As we remarked earlier, pure inspection strategies do not necessarily provide us with the best possible solutions, that is with the most eective utilization of our resources. We illustrate in Figure 5 the mixing principle with the simplest possible case: a single sensor. We choose
detection rate 100% 80% 60% 40% 20% budget 1 2 3 4 5
Figure 5: Mixing pure strategies. Random manual checks are plotted as a dotted line, while the single sensor based best strategy is drawn by a dashed curve. The best mixed strategy has three segments. The rst one, shown by a solid line, is a mixture of the single sensor strategy with the strategy letting every container go without any checking. The next segment, shown as a bold dashed line, coincides with the single sensor strategy, while the last segment, shown again as a solid line, is a mixture of the single sensor strategy with the one in which we manually check all containers. a relatively strong sensor for which the distributions are: fb (y) = 3e3y and fg (y) = 10e10y .
The ratio 10/3e7y corresponds to a separation of 3.5 standard deviations. The three pure strategies are: unpack everything manually (until the budget is exhausted; denoted by CHK); use the sensor on all the containers and with the remaining budget, unpack them manually in decreasing order of the odds that they are bad (that is, in decreasing order of the sensor reading y); let the containers go without any checking (denoted by OK). To make a more readable diagram we set the cost (per container) of using the sensor to be 1, and the cost of unpacking to be 5. The performance of the rst two pure strategies is shown in gure 5, while OK clearly has a detection rate equal to zero. The best mixed strategy is also shown in the gure. At very low budgets it is a mixed use of the compound strategy and OK (no examination at all). This can be represented by a diagram as the rst decision tree in Figure 6, where the small square represents the (random) mixing. As the budget increases we encounter a short range of budgets where the optimal solution tracks the compound strategy represented by the second decision tree in Figure 6. Finally the mixed strategy becomes a
Page 12
RRR 26-2006
M
1
S
M
1
OK
S
OK
CHK
S
CHK
OK
CHK
OK
CHK
Figure 6: Three mixed strategies corresponding to Figure 5.
mixture of CHK (unpacking manually) and the compound strategy, represented by the third decision tree in Figure 6. The same structure applies to any case in which there is a single constraint. In general however there will be several curves contributing to the determination of the convex hull. In terms of the polyhedron, the set of vertices determining the optimal point changes from time to time as the budget increases. Note however that as the budget increases, one pure strategy drops out of the optimal mixture, while the other remains active and then begins to mix with a new pure strategy.
3
A Linear Programming Model for Finding the Best Mixed Inspection Strategy
Following the idea of xed-grid search of [4], let us assume that for each sensor s S we x a partition s Ps = {Xj |j = 1, ..., k(s)}
s of the set of reals, and retain only the index j for each container such that Rs (C) Xj , and not the actual sensor reading Rs (C). Let us then denote by D the set of all possible decision trees utilizing the given partitions. For instance, Figure 7 shows a possible decision tree for the case when S = {a, b, c} and |Pa | = 3, |Pb | = 2 and |Pc | = 4. Note that some of the ranges could be merged in a decision tree D D. Note also that D now incorporates the information about the threshold selection, hence each D D represents an inspection strategy, which we denoted earlier by (D, t). Fixing the partitions of the sensor readings incorporated the grid-search over the grid corresponding to the partitions Ps , s S into our notations. This makes it possible to formulate the problem of nding a best mixed inspection as strategy a very large scale linear programming problem. Let us denote by xD the fractions of containers which we plan to inspect by inspection strategy D, for D D. Clearly, we have DD xD = 1 and xD 0 for all D D. Ideally,
RRR 26-2006
Page 13
if we would like to stick to a single decision tree as our inspection strategy, we would like to have only one of these xD variables to be positive. It turns out however, as we noted in the previous section, that it is advantageous to consider mixed strategies, i.e., non-integral values for the variables. If x = (xD |D D) represents a mixed strategy, then we have (x) =
DD
(D)xD
as its corresponding detection rate, C(x) =
DD
C(D)xD
as its corresponding average per-container inspection cost, T (x) =
DD
T (D)xD
as its corresponding average per-container inspection time, (x) =
DD
(D)xD
as its corresponding average rate of false positives, and Us (x) =
DD
Us (D)xD
as its corresponding average load of sensor s, for s S. It is important to note that all these measures are linear functions of the vector x. Consequently, we can write (in terms of x) the problem of nding a (mixed) inspection policy with the highest detection rate, subject to e.g., budget and capacity constraints as a linear programming problem: max
DD
(D)xD
subject to C(D)xD B
DD
(2) for all s S
Us (D)xD K(s)
DD
xD = 1
DD
xD 0
for all D D.
Page 14
RRR 26-2006
Of course, we could also include constraints on the per-container inspection time, or on the rate of false positives, etc. Similarly, the problem of nding the inspection strategy which minimizes the average per-container inspection costs, while still providing a prescribed high detection rate of 1 and not exceeding the available load capacities for each sensor can also be written as a linear programming problem: min
DD
C(D)xD
subject to (D)xD 1
DD
(3) for all s S
Us (D)xD K(s)
DD
xD = 1
DD
xD 0
for all D D.
Many other variants of these problems, involving limits on rates of false positives, inspection times, etc., can be formulated analogously as linear programming problems. The main practical diculty in using the above models is their size. The number of Q |S| variables in these formulations is |D| 2 sS k(s) 22 (double exponential in the number of sensors). As we recalled earlier, for 4 sensors, and binary partitions [4] enumerated over 11000 decision trees. With binary partitions of the sensor readings however we cannot achieve the best performance. To be able to utilize the power of the above models, we want a much ner partition of sensor readings, so that the optimization will be able to select the best decision sequence for the decision trees appearing in the optimal mixed strategy. If we allow e.g., k(s) = 10 for all sensors s S, the number of variables in (2) (or (3)) will jump to well above 100 million, making it almost impossible to generate all those and write down in full detail the linear programming problems (2) and (3). On the positive side, there are only 2 + |S| constraints in the above formulations, and thus the optimal basic solution x involves at most 2 + |S| decision trees with nonzero x D values. This suggests that such a mixed strategy is far from being impractical, since it is not a complicated matter to switch between a small number of decision trees. The small size of a basic solution also implies that column generation might be a good way of solving these very large linear programming problems. This technique was introduced rst by Gilmore and Gomory [5, 6] for the solution of cutting stock problems, but since then column generation has become a standard approach, used widely and eciently in numerous applications. For the sake of completeness we describe briey the main ideas.
RRR 26-2006
Page 15
Let us consider e.g., problem (2) and its linear programming dual, which has only 2 + |S| variables, corresponding to the budget constraint, s , s S corresponding to the load constraints, and corresponding to the last balance equality. Using these we can write the dual of (2) as follows: min B +
sS
K(s) s +
subject to (4) CD +
sS
Us (D) s + (D) , , s 0
for all D D for all s S.
Let us now consider a small subset D D of all decision trees, and suppose that we have solved (2) restricted to this set (i.e., we replace D everywhere in (2) by D). Let us denote by x the optimal solution to this restricted variant of (2), and let , s , s S, and denote the corresponding optimal solution of the dual of the restricted problem (i.e., where we keep constraints in (4) only for D D). The theory of linear programming then implies (see e.g., [3]) that x is an optimal solution for the unrestricted problem (2) if and only if min CD +
sS Us (D) s + (D)
DD
0.
(5)
If this is the case, we can stop, and x is our optimal solution. Otherwise, we choose a decision tree D for which the minimum in (5) is attained, increment D = D {D }, and return to the solution of the restricted (2) problem. Experience shows that on average, this approach generates only a small number of the variables in (2). In fact, as a consequence of the results of Khachiyan [8], the number of times we need to solve (5) is limited by a polynomial of the rank of (2), i.e., a polynomial function of |S| in our case. However, for this to work, one has to be able to solve subproblem (5) eciently. In its current form problem (5) still seem to require the complete enumeration of all decision trees, which is an approach that does scale well, as we noted before. In the following sections we show that problem (5) can in fact be reformulated such that it can be solved eciently.
4
Histories and Decision Trees
Let us note that each container C, subjected to a particular inspection policy D will produce a history (C, D) of sequential sensor readings, and a nal decision. We shall represent such a history as an ordered sequence of pairs (s, j), where s is a sensor on which C was tested and
Page 16
RRR 26-2006
s j is the retained index such that Rs (C) Xj . The order of these pairs in the history (C, D) is the order in which C was tested with the dierent sensors. Finally, the last element of (C, D) is the nal decision made on C, i.e, either OK if C was let go, or CHK if C was manually checked.
a
1 2 3
OK
1 2
c
3 4 1
b
2
OK
1
b
2
CHK
1 2
c
4 3
CHK CHK
OK
CHK
OK
Figure 7: A decision tree utilizing three sensors S = {a, b, c}, and partitions of the sensor readings such that |Pa | = 3, |Pb | = 2 and |Pc | = 4.
Let us consider the decision tree of Figure 7, using sensors S = {a, b, c} and partitions of sensor readings such that |Pa | = 3, |Pb | = 2 and |Pc | = 4. We can see that there are 12 dierent possible histories, which we list here below: 1 2 3 4 5 6 7 8 9 10 11 12 = {(a, 1), OK} = {(a, 2), (c, 1), OK} = {(a, 2), (c, 2), (b, 1), OK} = {(a, 2), (c, 2), (b, 2), CHK} = {(a, 2), (c, 3), (b, 1), OK} = {(a, 2), (c, 3), (b, 2), CHK} = {(a, 2), (c, 4), CHK} = {(a, 3), (b, 1), (c, 1), OK} = {(a, 3), (b, 1), (c, 2), OK} = {(a, 3), (b, 1), (c, 3), CHK} = {(a, 3), (b, 1), (c, 4), CHK} = {(a, 3), (b, 2), CHK}
For a given decision tree D let us denote by (D) the (nite) set of possible container histories, and let denote the set of all those histories, which may appear in some decision trees. Note that our notation assumes that we have xed the set S of sensors and also the partition Ps of sensor readings, for each s S.
RRR 26-2006 To each history (D) let us associate two parameters: g = P rob
(s,j) s Rs (C) Xj
Page 17
C is of type good
denoting the average fraction of good containers that have history , and b = P rob
(s,j) s Rs (C) Xj
C is of type bad
denoting the average fraction of bad containers with this history. Since every container of either type must have some specic path, these parameters satisfy the equalities g =
(D) (D)
b = 1.
(6)
For instance, if we assume that the sensors are stochastically independent we can compute these parameters from the sensor distributions as follows: g =
(s,j) s (Xj , s, good) and b = (s,j) s (Xj , s, bad).
(7)
Let us add that the stochastic independence of the sensors is not an important assumption. What is needed for our model as input are the quantities g and b for all paths . Let us observe next that the performance characteristics of an inspection strategy D introduced in the previous section are in fact linear functions of the parameters g and b , (D). For instance, the detection rate (D) of an inspection strategy D, that is the probability that a bad container will get manually inspected by D, can be written as (D) =
(D) CHK
b .
Let us recall that, we assume that truly dangerous containers are very rare. Consequently, all other operating characteristics of an inspection strategy D will be dominated by the good containers (which are the vast majority of containers actually tested by this policy). The rate of false positives (D) of D, that is the expected fraction of good containers that end up manually inspected (unnecessarily) can be written as (D) =
(D) CHK
g .
Page 18
RRR 26-2006
For a sensor s S the sensor load Us (D) is the fraction of containers that are actually tested on sensor s by the inspection policy D, and it can be written as Us (D) =
(D) (s,j) for some j
g .
The average per-container inspection cost C(D) can be computed as C(D) =
sS
Cs Us (D) + CCHK (D),
where, as before, Cs denotes the average per-container inspection cost of sensor s, and CCHK denotes the average per-container cost of manual inspection. Analogously, the average per-container inspection time T (D) can be written as T (D) =
sS
Ts Us (D) + TCHK (D),
where, as before, Ts denotes the average per-container inspection time of sensor s, and TCHK denotes the average per-container time of manual inspection. There are several further operational characteristics which might be important to consider (such as container delays; separate detection rates for dierent types of dangerous containers when we do not know the relative frequency of those, separate false positive rates for dierent types of good containers, etc.). Since we do not have at this time reliable data for those, and since we want to keep our presentation simple, we do not consider other measures in this study.
5
The Polytope of Decision Trees
Let us consider a given set S of sensors, a xed partition Ps of the possible sensor readings for each s S, and recall that we denote by the set of all possible container histories, while (D) denotes the container histories of containers inspected by strategy D. In other words, we have = (D).
DD
Furthermore, any decision tree D can, equivalently, be represented by the vectors (g | (D)) and (b | (D)). Let us extend these vectors with zeros, and dene in this way g(D), b(D) R , where g(D) = g if (D) 0 if \ (D) and b(D) = b if (D) 0 if \ (D)
RRR 26-2006
Page 19
In fact a mixed strategy x = (xD | D D) can also be represented naturally by the fractions of good and bad containers having a particular history . The corresponding vectors g(x), b(x) R are convex combinations: g(x) =
DD
xD g(D) and b(x) =
DD
xD b(D)
Thus, we can associate naturally two polytopes, Pg and Pb , to the set of sensors and the xed partitions of their readings, corresponding to the two container types we considered. These polytopes are dened as Pg = ConvHull g(D) | D D and Pb = ConvHull b(D) | D D = b(x) xD = 1 xD 0 for all D D
DD
=
g(x)
xD = 1 xD 0 for all D D
DD
.
To characterize these polytopes in terms of a linear system, we need to introduce more notation and terminology. Recall that a history of a container C is an ordered sequence of pairs (s, j), where s is a sensor and j is the index of one of the parts of the partition Ps = {Xis | i = 1, ..., k(s)}, and the last element of is one of the symbols OK and CHK. Let us now consider an arbitrary ordered sequence of pairs (s, j), where s S, and j is an index satisfying 1 j k(s), and such that no sensor appears twice in . Let us call such a sequence a pre-history, and for two such sequences , let us write if is an initial subsequence of . For instance, if we have the history = ((a, 3), (b, 1), (c, 1), OK), then the sequences 1 = ((a, 3)), 2 = ((a, 3), (b, 1)) and 3 = ((a, 3), (b, 1), (c, 1)) are all pre-histories (of ), and we have 1 2 3 . Clearly, if and are pre-histories, and they do not involve the same sensors, then their concatenation (, ) is also a pre-history. The following statement provides a linear characterization of the polytopes Pg and Pb . This result turns out to be very useful for solving problem (2) (and other variants). We state the characterization only for polytope Pg , since the perfectly analogous claim can easily be shown in the same way for Pb , as well. To simplify notation, we assume that the sensors are stochastically independent, and parameters g and b are computed by (7) for all . Let us add however that this is just a technical detail, and analogous linear description for the polytopes Pg and Pb can be derived for stochastically dependent sensors, as well. Theorem 1 Given y R , we have y Pg if and only if y satises the following equalities: y = 1
(8)
Page 20 and
RRR 26-2006
s y = (Xj , s, good) :(,(s,j))
k(s)
y (9)
:(,(s,i))
i=1
for all sensors s S, for all indices 1 j k(s), and for all pre-histories not involving sensor s. Proof. We prove the statement by peeling o vectors g(D) from y, until we arrive to y = 0, each time decreasing the number of histories for which y > 0. To be able to do this peeling o process, we set a real parameter = 1, and claim that if y satises equalities (8) and (9), then there exists a decision tree D D such that y > 0 for all (D). For such a decision tree D, we set y . = min (D) g Due to (6) and (8), we must have 0 1. Let us set xD = . If = 1, then we must have y = g(D), and our proof is complete. Since both y and g(D) satises equalities (8) and (9), by setting y g(D) y= and = (1 ). 1 we obtain a vector y which also satises (8) and (9), and which has at least one less history with a positive y value than vector y had in the previous step. Thus, replacing y by y , by , and repeating the above steps, after a nite number of iterations we can arrive to a convex combination of the vectors g(D), D D which is equal to y, completing the proof of our statement. What is left is to show our claim, namely that if y satises equalities (8) and (9), then there exists a decision tree D D such that y > 0 whenever g > 0 (recall that g(D) = (g | (D)). We prove this claim by induction on the number of sensors. Clearly, if we have only one sensor, then the claim is almost trivial, since we only have to try all possible nal decisions at the leafs of the only possible trivial decision tree (consisting of only one node). Let us now assume that we proved this claim for up to k sensors, and consider the case when |S| = k + 1. Consider an arbitrary history for which y > 0, and let s S be the rst sensor appearing in this history (i.e., for which ((s, j)) for some 1 j k(s)). Since (Xi , s, good) > 0 for all i = 1, ..., k(s), we must also have histories with y > 0 and with ((s, i)) for all i = 1, ..., k(s). Let us denote by the history obtained from by deleting the rst sensor from , and consider the sets i = { | , ((s, i)) }
for i = 1, ..., k(s). Each set i is the set of all possible histories for the smaller set S \ {s} of sensors. Furthermore, for an index 1 i k(s) let us dene i =
,((s,i))
y
RRR 26-2006 and a new vector yi Ri dened by 1 i y = y b i
Page 21
for all i .
Then, yi is a vector satisfying equations (8) and (9), with respect to sensors S \{s}. Thus, we i must have a decision tree Di such that y > 0 for all (Di ) by our inductive hypothesis. b Since this holds for all i = 1, ..., k(s), by gluing these decision trees Di , i = 1, ..., k(s) to sensor s as a root note, we obtain a decision tree D rooted at sensor s, satisfying our claim.
6
A Computationally Tractable LP Model
To provide an ecient linear programming based scheme, we note that by Theorem 1 every decision tree can be viewed as a particular (vertex) y Pg , and conversely, every vector y Pg can be viewed as a mixed strategy. Thus, to solve (5) all we need is to represent D D by the corresponding y vector, and perform the minimization over y Pg instead of D D. If we can do this, such that we obtain an expression which is linear in y, then we can solve (5) as a linear programming problem. To this end, let us associate cost, detection rate, load parameters, etc., to every history . Let us denote by S() the set of sensors appearing in , and dene Cs if OK , sS() C() = CCHK + Cs if CHK ,
sS()
0 () =
if OK ,
b if CHK , where b is dened in (6), and 1 if s S(), Us () = 0 if s S(),
for s S. Furthermore, to a vector y Pg we associate C(y) =
C()y , ()
(y) = Us (y) =
y , and g
Us ()y for s S.
Page 22
RRR 26-2006
In particular, if y = g(D) for a decision tree D D, then we have C(y) = C(D), (y) = (D) and Us (y) = Us (D) for s S, exactly as those parameters were dened in Section 2.2. Thus, we can rewrite the optimization problem in (5) as the following linear programming problem: min Cy +
sS Us (y) s + (y) .
yPg
(11)
Note that in the subsequent iterations only the objective function in (11) changes, and hence we can keep the previously optimal basis, and just re-optimize it with the new objective, which perhaps will save time, and increases eciency.
7
Computational Experiments
We solve model A (problem (4), maximization of detection rate subject to given budget) and model B (problem (5), minimization of operating costs per container subject to given detection rate) for 4 sensors and for up to 7 thresholds, assuming unlimited capacities of each sensor. Using the same data as in [4] we assume that at each sensor bad and good containers have readings distributed normally with variance 1. Readings of good containers have mean 0 for all sensors. For the rst sensor, half of the readings of bad containers have mean K1 = 4.37, another half has the same distribution as good containers (see Figure 8), or in other words sensor 1 can distinguish only one out of two dierent types of bad containers.
good
bad
Figure 8: Distributions of sensor readings of sensor 1 for good and bad sensors, indicated respectively by solid and dotted curves. The cost of inspection of a container by this sensor is C1 = 0.32.
RRR 26-2006
Page 23
good
bad
Figure 9: Distributions of sensor readings of sensor 2 for good and bad sensors, indicated respectively by solid and dotted curves. The cost of inspection of a container by this sensor is C2 = 0.92.
good
bad
Figure 10: Distributions of sensor readings of sensor 3 for good and bad sensors, indicated respectively by solid and dotted curves. The cost of inspection of a container by this sensor is C3 = 57.
For bad containers the mean at sensor s = 2, 3, 4 is Ks (see Figures 9,10 and 11), w...