FUZZY AND POSSIBILISTIC PROGRAMMING TECHNIQUES IN THE
RADIATION THERAPY PROBLEM: AN IMPLEMENTATIONBASED
ANALYSIS
by
Elizabeth A. Untiedt
B.A., Princeton, 1997
A thesis submitted to the
University of Colorado at Denver
and Health Sciences Center
in partial fulfillment
of the requirements for the degree of
Master of Science
Mathematical Sciences
July 5, 2006
This thesis for the Master of Science
degree by
Elizabeth A. Untiedt
has been approved
by
K. David Jamison
Date
DEDICATION
For C.J., who always says Im a professional at math.
ACKNOWLEDGMENT
This thesis would not have been possible without the support of my loving
husband, Bill, who afforded me many unplanned (or underplanned) respites
from my familial obligations as various deadlines loomed. In addition to logis
tical support, he provided emotional and financial support for this monumental
undertaking.
I also owe a debt of gratitude to my adviser, Weldon Lodwick, for providing
seasoned and knowledgeable mathematical guidance, tempered with an under
standing of my need to juggle academic pursuits and family needs. His input
during the last critical weeks transformed this thesis from a formless collection
of facts and findings into a focused result.
None of this work would have resulted in a diploma without the help of
Marcia Kelly, who kept me advised of required documents and deadlines, and
pulled off a few eleventh hour miracles to keep me on track.
Untiedt, Elizabeth A. (M.S., Mathematical Sciences)
Fuzzy and Possibilistic Programming Techniques in the Radiation Therapy
Problem: An ImplementationBased Analysis
Thesis directed by Professor Weldon Lodwick
ABSTRACT
Fuzzy and possibility theory have made contributions to the field of opti
mization under uncertainty for over 30 years. A plethora of formulations have
been proposed in the literature, but few have been implemented in practice.
This paper investigates the theory, semantics, and efficacy of a selection of sig
nificant fuzzy and possibilistic optimization algorithms via their application to
a wellknown largescale problem: the radiation therapy planning problem. The
algorithms are compared, critiqued, and organized with the following objective
in mind: to guide a decision maker in the selection and implementation of a
fuzzy or possibilistic optimization algorithm.
This abstract accurately represents the content of the candidates thesis. I
recommend its publication.
Signed ____________________________
Weldon Lodwick
m
CONTENTS
Figures .............................................................. vii
1. Introduction....................................................... 1
2. Concepts and Semantics of Fuzzy and Possibilistic Optimization ... 3
2.1 Fuzzy Number Theory................................................... 3
2.2 Fuzzy Relations..........;....................................... 9
2.3 Fuzzy Optimization.................................................. 13
3. Review of Selected Possibilistic and Fuzzy Programming Formulations 17
3.1 Fuzzy Programming ................................................... 17
3.1.1 Bellman and Zadeh................................................. 18
3.1.2 Tanaka, Okuda, and Asai........................................... 19
3.1.3 Zimmer mann...................................................... 21
3.1.4 Ver degay ........................................................ 23
3.1.5 Surprise.......................................................... 26
3.2 Possibilistic Programming........................................... 27
3.2.1 Buckley........................................................... 29
3.2.2 Tanaka, et al..................................................... 33
3.2.3 Fuzzy Max......................................................... 37
3.2.4 Jamison and Lodwick .............................................. 40
3.2.5 Luhandjula........................................................ 42
3.3 Programming With Fuzzy and Possibilistic Components ................ 44
IV
3.3.1 Robust Programming .............................................. 44
3.3.2 Delgado, et al................................................... 47
3.3.3 Comments on Fuzzy/Possibilistic Programming...................... 48
4. Application of Possibilistic and Fuzzy Programming Techniques to the
Radiation Therapy Problem............................................ 50
4.1 The Radiation Therapy Planning Problem.......................... 50
4.2 Fuzzy RTP Problem.................................................. 53
4.2.1 Flexible Programming............................................. 53
4.2.2 Verdegay ........................................................ 55
4.2.3 Surprise......................................................... 57
4.3 Possibilistic RTP Problem...................................... 60
4.3.1 Buckley.......................................................... 61
4.3.2 hLevel Technique ............................................... 66
4.3.3 Fuzzy Max Technique ............................................. 69
4.3.4 JamisonLodwick Penalization Method.............................. 72
4.3.5 Robust Programming .............................................. 76
4.4 The Crisp RTP Problem.............................................. 81
5. Comparison of Optimization of Uncertainty and Techniques ......... 83
5.1 Formulation Selection ............................................. 83
5.1.1 Hard or Soft Constraints......................................... 83
5.1.2 Fuzzy or Possibilistic Parameters................................ 83
5.1.3 Feasibility or Optimality........................................ 84
5.1.4 Fuzzy of Crisp Solution.......................................... 85
5.1.5 Complexity and Speed............................................. 86
v
5.2 Additional Thoughts............................................ 86
6. Avenues for Further Research..................................... 88
Appendix
A. Symbols and Terms................................................ 90
References.......................................................... 91
vi
FIGURES
Figure
2.1 Symmetric triangular fuzzy number, x = (cc, ft)........................ 8
2.2 Triangular fuzzy number, x = (a, /?,7)..................................... 9
2.3 Trapezoidal fuzzy interval x = (a, (3,6,7)................................ 10
2.4 An illustration of a LR fuzzy interval A (aR, aL, fiR, (3l)lr 11
2.5 Ranking indexes for fuzzy intervals. Heavy lines represent fiM and
Hn, and light lines represent 1 hm and 1 yjy. Point a is
Nec[a; > y], Point b is Poss[x+ > y~], point c is Poss[(a:+ > y+],
and Nec[x_ > y+] = 0 is not shown...................................... 12
2.6 Illustration of a fuzzy threshold. The membership function of M is
represented by line (a). A fuzzy threshold of type [M,+ inf), the
fuzzy set of numbers that is possibly greater than y, is represented
by line (b). A fuzzy threshold of type ]M,+ inf), the fuzzy set of
numbers that is necessarily greater than y, is represented by line (c). 15
2.7 Illustration of Weak Feasibility (a), and Strong Feasibility II (b), and
Strong Feasibility I (c). Heavy lines represent and /j,Ni and light
lines represent 1 hm and 1 ................................ 16
3.1 An illustration of a fuzzy decision...................................... 19
3.2 An illustration of sup^a A maxcQ I^g(^)] sup.,, y^..................... 21
vii
3.3 Derivation of surprise function for (Tx)i < bramor The membership
function, HbTumor{{Tx)i) is on the left, and the surprise function,
s(x\(Tx)i < bTumor) is on the right................................
3.4 Example with d (2,2) and e = (3,2)...............................
3.5 An illustration of Yi>0............................................
3.6 Illustration of fuzzy max, c = max(a, b). Fuzzy max might be
thought of as the rightmost of the two membership functions. . .
3.7 Illustration of A >h B.............................................
4.1 Dose volume histogram for RTP problem solved via Zimmermanns
fuzzy technique. Optimal a* = .22..................................
4.2 Dose volume histogram for RTP problem solved via Verdegays fuzzy
technique. Snapshot taken at a = 0.................................
4.3 Dose volume histogram for RTP problem solved via Verdegays fuzzy
technique. Snapshot at a = .22.....................................
4.4 Dose Volume Histogram for tumor dosage. The heavy line repre
sents ideal dosage distribution, dotted line represents typical dosage
distribution, and hatched areas represent areas the surprise function
should penalize....................................................
4.5 Dose volume histogram obtained with fuzzy surprise methodology. .
4.6 Dose volume histogram for RTP problem with possibilistic right hand
side, Buckleys method, a 0.22...................................
4.7 Dose volume histogram for RTP problem with possibilistic right hand
side, Buckleys method a = 0.......................................
27
32
36
38
39
55
57
58
59
60
62
63
4.8 Dose volume histogram for RTP problem with possibilistic right hand
side and A matrix, a = 0.70.......................................... 64
4.9 Dose volume histogram for RTP problem with possibilistic right hand
side and A matrix, a = 0.50.......................................... 65
4.10 Dose volume histogram for RTP problem with possibilistic right hand
side and A matrix, a = 0.20.......................................... 66
4.11 Dose volume histogram for RTP problem with possibilistic right hand
side and A matrix, a = 0.00......................................... 67
4.12 Dose volume histogram for RTP problem solved via Tanakas hlevel
technique, h* = .21..........................:....................... 68
4.13 Dose volume histogram for RTP problem solved via the fuzzy max
technique, h* = .05.................................................... 72
4.14 In the figure on the left, A
4.15 Dose volume histogram for Jamison Lodwick method....................... 77
4.16 Dose volume histogram for robust formulation, with a resolution of 3. 81
4.17 Dose volume histogram for the crisp formulation........................ 82
IX
1. Introduction
The goal of optimization, in theory as well as in practice, is to achieve an
outcome which best satisfies a set of goals and constraints set forth by the deci
sion maker. Since much of the decisionmaking in the real world depends upon
constraints, goals, and consequences of potential actions which are not explicitly
known, the field of optimization is faced with the challenge of quantifying impre
cise data in a meaningful way. The most welldeveloped branch of optimization
under uncertainty is stochastic programming, which utilizes probability distri
butions for all uncertain data. Not all uncertainty, however, is inherently ran
dom. Some uncertainty stems from vagueness about the quantitative meaning
of a constraint or goal (as in, We wish for a return on investment significantly
higher than the prime rate.) Some uncertainty stems from a lack of information
about, or ambiguity concerning, the value of a goal or constraint (as in, the
body temperature which constitutes a fever is around 100 degrees.) Known
respectively as fuzziness and possibility, these two forms of uncertainty require
quantitative representation which differs from their probabilistic counterparts.
Fuzzy and possibility theory have made contributions to the field of op
timization under uncertainty for over 30 years. The literature contains many
proposed formulations, but few have been implemented in practice. The desired
outcome of this paper is to guide a decision maker in the selection and implemen
tation of a fuzzy or possibilistic optimization algorithm. Chapter Two provides
necessary background on fuzzy and possibility theory, and on optimization under
1
uncertainty. Chapter Three investigates the theory and semantics of a selection
of significant fuzzy and possibilistic optimization algorithms. Chapter Four in
vestigates the efficacy of a selection of these algorithms via their application to
a wellknown largescale problem: the radiation therapy problem. In Chapter
Five, the algorithms are compared, critiqued, and organized for the ease of use
for potential decisionmakers.
2
2. Concepts and Semantics of Fuzzy and Possibilistic Optimization
In this chapter, we define fuzzy and possibilistic numbers. We delve into
their semantics and the relationship between the two. An examination of fuzzy
and possibilistic binary operations will set the stage for fuzzy and possibilistic
optimization.
2.1 Fuzzy Number Theory
Fuzzy entities are sets without sharp boundaries in which there is a transition
between elements that are members and elements that are nonmembers. Grades
of membership in the set are defined by a membership function.
Definition 2.1 Let X denote a set. Then a fuzzy set A in X is a, set of ordered
pairs
A={{x,/j,a(x))}, x e X (2.1)
where the membership function p,A X M C [0,1] is a map from set X
into membership space M. When M contains only 0 and 1, A is nonfuzzy (or
crisp), and its membership function is identical to the characteristic function of
the nonfuzzy set.
A normalized fuzzy set has a maximal membership function value of 1.
A fuzzy number is a convex, normalized fuzzy set on the real number line
whose membership function is at least segmentally continuous and has a unique
3
modal value. For examples of fuzzy numbers, see triagular and symmetric trian
gular numbers, defined below. A fuzzy interval is a fuzzy set on the real number
line with
A fuzzy interval is an uncertain set on the real line with a mean interval
whose elements have p(x) = 1. As in fuzzy numbers, the membership function
must be convex, nomalized, and at least segmentally continuous.
A possibilistic entity, x, is an entity which is known to exist, but whose
exact value is uncertain. The possibility that the entitys value is a member of
a particular interval I is defined by a possibility measure.
Definition 2.2 Let Cl be the set of real numbers. Then a possibility measure,
n, is a setfunction defined on the set of subsets of Cl, which assumes values in
[0,1] and which satisfies the following axioms:
(0 n(0) = 0;
(ii) n(ft) = 1;
(Hi) VI, VAi C Cl,i Â£ I, II(UiejAj) = sup{H(Ai)\i e I).
If an event, E, is certain, then the possibility of an element x with respect to
event E, IIe:(x) = 1 if x D E ^ 0, and II^x) = 0 if x D E = 0. We define a
possibility distribution, Poss(cu), in terms of the membership function of E:
Vcu fi,Poss({u;}) = /ue(u>). (22)
4
When the measure of the underlying space is finite, there is a companion
setfunction to possibility called the necessity, AT, defined by
M(I) = 1U(IC), (2.3)
where Ic denotes the complement of I. A necessity distribution, Nec(w), is
derived from necessity measure M just as the possibility distribution Poss(w) is
derived from possibility measure II in (2.2). Semantically, the necessity of an
event measures the impossibility of the opposite event.
Throughout the paper, a fuzzy set is denoted by a tilde over the letter, and a
possibilistic entity is denoted with a hat. For example, x denotes a fuzzy subset
of X, and x denotes a possibilistic member of X.
The following basic concepts regarding fuzzy numbers will be useful in later
sections of the paper. Because possibility functions are related to membership
functions by (2.2), these concepts can, in general, be extended to possibilistic
entities as well.
Support. The support of a fuzzy number A is a set {x\/jla(x) > 0.}
Equality. Two fuzzy sets, A and B are said to be equal if and only if
^a(x) = Hb{x) for aU x Â£ X. 1
Containment. A fuzzy set A is said to be a subset of fuzzy set B if and
only if pla{x) < Eb{x) for all x G X.
Intersection. The intersection of A and B is defined as the largest fuzzy
set contained in both A and B. The membership function of A A B is
1This definition, given by Bellman and Zadeh [1], is generally accepted. We explore broader
interpretations of A = B in Section 2.2.
5
given by
Vaab(x) = TDin(fiA(x),fMB(x))t x E X. (2.4)
2
Relation. A fuzzy relation R in the product space X x Y is a fuzzy set
characterized by a membership function [jlr which associates with each
ordered pair a grade of membership /xr(x, y) in 5ft.
Extension Principle. For our purposes, the extension principle is used to
extend usual functions to fuzzy setvalued arguments. [32]. Let / :U>V,
and let F be a fuzzy set in U. To obtain the fuzzy set f(F) on V, we first
extend / to set valued arguments via
f(A) = if(u) ueA}. (2.5)
We then apply the extension principle:2 3
f(F)(v) = snp{F(u) : v = f(u)} (2.6)
= Oif{v = /(w)} = 0. (2.7)
Noninteractive. Consider a fuzzy number C where C C SR2 which is a
direct product of two fuzzy numbers A e 5? and B E 3? such that for all
(x,y) E 5ft2,
Vc(x,y) = Va{x) AfjLB(y). (2.8)
2Minimum is only one of a continuum of operators that define intersection, but this dis
cussion is beyond the scope of this paper.
3This is the extension principle for finite domains, with which this paper is concerned. It
may or may not hold in the continuous case, which is beyond the scope of this paper.
6
If condition (2.8) holds for A and B, then A and B, are said to be non
interactive. Semantically, two numbers are noninteractive if they can be
assigned values independently of each other [4]. Every fuzzy and pos
sibilistic model examined in this paper operates on the assumption of
noninteractivity of all uncertain entities.
Symmetric triangular. A symmetric triangular fuzzy number, A = (a, (5),
has a membership function /Ja centered at a value a, with a spread (3 such
that
flA(x) = <
1 for x = a,
0 for x < a /?,
0, for x > a + /?,
1 ba 1 p for x {a P, a + /?)
(2.9)
A symmetric triangular fuzzy number is illustrated in Figure 2.1.
Triangular. A triangular fuzzy number, A = (a,/?,7), has a membership
function 11 a centered at a value a, with a support (/?, 7) such that
1 for x = a,
0 for x < a /3,
o, for x > a + 7,
1 xa 1 (3 3 for x G (a P, a)
1l=x 7 for x e (a, a + 7)
(2.10)
A triangular fuzzy number is illustrated in Figure 2.2.
7
mu(x)
Figure 2.1: Symmetric triangular fuzzy number, x (a, (3).
Trapezoidal. A trapezoidal fuzzy interval, A = (a, ft, 7,5), has a member
ship function Pa with a core (a, ft) and a support (7,5) such that
I* a 0)
1 for x e [a,/?],
0 for x Â£ (7,6),
1 for x E (7,0),
1fEf. torxeW.i).
A trapezoidal fuzzy interval is illustrated in figure 2.3.
(2,11)
LR. [4] A LR fuzzy interval A = (aR, aL, (3R, [3l)lr has a membership
function pa and reference functions L and R. L : [0, inf) > [0,1] and
R : [0, inf) > [0,1] are upper semicontinuous and strictly decreasing in
8
0
alpha beta
alpha
alpha+gamma
Figure 2.2: Triangular fuzzy number, x = (a, (5,7).
the range (0,1], and is defined as follows:
/j,A(x) = <
1 for x G [aL.
0 for x Â£ (aL
U^) for x G [aL
ry( xaR pk for x G [aR
aL)
(2.12)
An LR fuzzy interval is illustrated in figure 2.4.
2.2 Fuzzy Relations
Optimization models usually involve constraints consisting of equalities, in
equalities, or both. For deterministic optimization, the meaning and computa
tion of the constraint set is clear. In optimization under uncertainty, however,
the meaning and computation of equality and inequality must be deter
mined. To this end, Dubois [2] gives a comprehensive analysis of fuzzy relations
9
mu(x)
alpha beta
gamma delta
Figure 2.3: Trapezoidal fuzzy interval x = (a, (3,5,7).
in his 1987 paper. He gives four interpretations of fuzzy equalities, and four
interpretations of fuzzy inequalities. These relations are outlined below, using
crisp intervals M = [m~,m+] and N = [n~,n+] and their fuzzy equivalents, M
and N.
The statement M > N can be interpreted in any of the four following ways:
i. V.x Â£ M, Vy Â£ N, x > y. This is equivalent to m~ > n+. The correspond
ing fuzzy interval relation is Nec(x_ > y+) = inf{max(l ym(u),1
Hn{v))\u < u}.
ii. Vx Â£ M, By Â£ TV, x > y. This is equivalent to m > n~. The cor
responding fuzzy interval relation is Nec(x_ > y~) = inf {sup {max (1
u v
fj,M(u),yN(v))\u > u}}.
iii. Bx Â£ M, Vy Â£ N, x > y. This is equivalent to m+ > n+. The correspond
10
Figure 2.4: An illustration of a LR fuzzy interval A = (aR, aL, f3R, (3l)lr
ing fuzzy interval relation is Poss(x+ > y+) = sup{inf{min(/i^(u), 1
u v
/ijv(u))u < u}}.
iv. 3(x,y) Â£ M x N, x > y. This.is equivalent to m+ > n~. The correspond
ing fuzzy interval relation is Poss(a;+ > y~) = sup{min(/iJv'(ti), hn(v))\u >
w}. The .four inequality relations are illustrated in Figure 2.5. Dubois de
fines four possible interpretations of the statement M = N, here given in
descending order of stringency:
v. Zadehs fuzzy set equality: Hm = Pw (which occurs'when both (vi) and
(vii) are satisfied at the same time.)
vi. Vx Â£ M, 3y Â£ N, x y (which is equivalent to M C N). The corre
sponding fuzzy interval relation is Poss(x = y) II (M, N).
11
Figure 2.5: Ranking indexes for fuzzy intervals. Heavy lines represent hm
and Hn, and light lines represent 1 Hm and 1 /^y. Point a is Nec[x~ > y],
Point b is Poss[a;+ > y~], point c is Poss[(a:+ > y+], and Nec[a;_ > y+] 0 is
not shown.
vii. Vy
sponding fuzzy interval relation is Nec(x y,y) = A/jv(M).
viii. 3(x,y) Â£ M x N, x = y (which is equivalent to iV fl M ^ 0.) The
corresponding fuzzy interval relation is Nec(y = X] x) = J\[m(N).
Inequality relation (i) indicates that x is necessarily greater than y, 4 (iv) in
dicates that x is possibly greater than y, 5 and (ii) and (iii) fall somewhere in
between. Similarly, equality relation (v) indicates that x is necessarily equal to
y, (viii) indicates that x is possibly equal to y, and (ii) and (iii) fall somewhere
4This is the pessimistic view. The decision maker who requires that m~ > n+ in order to
satisfy M > N is taking no chances [14].
5This is the optimistic view. The decision maker who merely requires that m+ > m~ in
order to satisfy M > N has a hopeful outlook.
12
in between.
2.3 Fuzzy Optimization
Optimization aims to allow a decisionmaker to achieve the best possible
outcome, while meeting a particular set of requirements. Traditional (crisp)
optimization treats outcomes and requirements that are definitive and quantita
tive, and stochastic optimization solves problems with parameters characterized
by probability distributions. Fuzzy and possibilistic optimization solve problems
for which the outcomes and/or the requirements are imprecisely or qualitatively
defined.
Fuzzy and possibilistic optimization, while similar, are motivated by very
different semantics [12]. In general, a lack of information about the values of the
model parameters lead to possibilistic optimization. Some researchers [6] prefer
to think of this lack of information as ambiguity. When the decisionmaker
evinces some flexibility regarding the requirements, or the parameters are vague
in nature, fuzzy optimization is in order. It has also been said the that possibilis
tic uncertainty is databased, and fuzzy uncertainty is preferencebased [8]. The
membership grade of a fuzzy goal has been thought of as representing the de
gree of satisfaction, and that of a possibility distribution has been thought of as
representing degree of occurrence [6]. Table 2.1 provides a quick reference for
the differences between fuzzy and possibilistic optimization. The chapters that
follow examine examples of fuzzy programming and possibilistic programming,
as well as programming that combines fuzzy and possibilistic elements.
A key element of the constrained mathematical program is the equality,
inequality, or set containment that describes each constraint. The set of fea
13
Fuzzy Optimization Possibilistic Optimization
fuzzy
membership functions
preferencebased
vague
subjective
degree of satisfaction
imprecise
possibility distributions
informationbased
ambiguous
subjective or objective
degree of occurence
Table 2.1: Characteristics of fuzzy and possibilistic optimization.
sible solutions are those which satisfy these constraints. It is natural, then,
that the concept of feasibility in the fuzzy mathematical program be derived
from the definitions of fuzzy equalities, fuzzy inequalities, and fuzzy set con
tainment. Dubois inequality and equality relations (see Section 2.2), used to
define constraints, lead to the concepts of strong and weak feasibility. The feasi
bility indices introduced below correspond with the relation indices enumerated
in Section 2.2. We assume that the fuzzy programming problem has equality
constraints of the form f(x,a) = g(x,b), where a and b are vectors of LR fuzzy
numbers. We further assume inequality constraints of the form f(x, a) < g(x, b)
where a is a vector of LR fuzzy numbers and b is a fuzzy threshold, as illustrated
in Figure 2.6.
For inequality constraints:
i. Very Strong Feasibility: VSF{x) = J\f(f(x, a)~ > g(x,b)+)
14
Figure 2.6: Illustration of a fuzzy threshold. The membership function of M
is represented by line (a). A fuzzy threshold of type [M, + inf), the fuzzy set
of numbers that is possibly greater than y, is represented by line (b). A fuzzy
threshold of type ]M, + inf), the fuzzy set of numbers that is necessarily greater
than y, is represented by line (c).
ii. Medium Strong Feasibility: MSF(s) = Af(f(x,a)~ > g(x,b)~)
iii. Medium Weak Feasibility: MWF(x) = U(f(x,a) + > g{x,b)+)
iv. Very Weak Feasibility: VWF(x) U(f(x, a))+ > g(x, b)~).
For equality constraints:
v. Weak Feasibility: \fx, WF(x) = II(/(a;,a) = g(x,b))
vi. Strong Feasibility I: \/x, SF(x) = N{g{x,b) = f(x, a); a)
vii. Strong Feasibility II: Vx, SF(x) = Af(f(x, a) = g(x, b)]b).
15
Figure 2.7: Illustration of Weak Feasibility (a), and Strong Feasibility II (b),
and Strong Feasibility I (c). Heavy lines represent fiM and (An, and light lines
represent 1 [jlm and 1 [iN.
As we shall see in the chapters that follow, these feasibility definitions are in
corporated, or can be incorporated, in many proposed fuzzy and possibilistic
programming formulations.
16
3. Review of Selected Possibilistic and Fuzzy Programming
Formulations
There have been numerous formulations set forth for modeling imprecise
objectives and constraints, each with associated input semantics and solution
interpretations. This chapter organizes and reviews a representative selection of
fuzzy and possibilistic formulations, many which will be implemented in Chapter
Four.
3.1 Fuzzy Programming
Vague parameter values leads to fuzzy programming. When the vagueness
represents a willingness on the part of the decisionmaker to bend the constraints,
fuzzy programming can also be called flexible programming. In these cases, the
decisionmaker is willing to lower the feasibility requirements in order to obtain
a more satisfactory objective function value, or, in some cases, simply in order
to reach a feasible solution. Such flexible constraints are commonly referred to
as soft constraints.
In the case of soft constraints, it is the inequality (or equality) itself which is
viewed to be fuzzy (i.e. Ax
hand side has vague value, in which case the right hand side is viewed to be fuzzy
(i.e. Ax
right hand sides, consider the example of a computer dating service. Suppose
woman A specifies that she would like to date a mediumheight man. Woman
A defines medium as a fuzzy set characterized by a triangular membership
17
function, centered at 69, with a spread of 3. This is a hard constraint because
the man is required to be mediumheight, but mediumheight is a fuzzy set. This
is, therefore, an example of a fuzzy righthand side (i.e. Ax = b). Now consider
woman B, who desires to date a man of 69. Unlike woman A, woman B is
willing to compromise a little on the height requirement in order to be matched
with a date who meets some of her other requirements. Her satisfaction level
with a 69 man is 1, with a 68 or 70 man is , and with a 67 or 71 man
is . This is an example of a soft constraint, represented by a fuzzy equality
(i.e. Ax=b). Notice that the membership functions in the two fuzzy cases are
the same (symmetric triangular centered at 69 with a spread of 3), but the
semantics are different. This section contains formulations of fuzzy programs
with both types of fuzziness. We shall see how the semantics of the problem
influence its formulation.
3.1.1 Bellman and Zadeh
In a landmark paper [1], Bellman and Zadeh were the first to propose an
approach to mathematical programming with fuzzy components. Conventional
mathematical programming has three principal components: A set of alterna
tives, X, a set of constraints, C, which limit the choice of alternatives, and an
objective function, z, which measures the desirability of each alternative. Bell
man and Zadeh propose that, in a fuzzy environment, a more natural framework
is one in which goal(s) replace the objective function(s), and a symmetry be
tween goals and constraints erases the differences between them.
Each fuzzy goal, Gj, or constraint, C\, is defined as a fuzzy subset of the
solution alternatives, X, via a membership function (/^(x) or Then a
18
fuzzy decision, D, is defined as a fuzzy set resulting from the intersection of G
and G, and is characterized by membership function nn as follows:
Pd{x) = A AteO) (3.1)
= min \fMc(x),HG(x)].
An optimal decision, in turn, is one which maximizes fiD. Refer to Figure 3.1
for a visual interpretation of /x^, /j,c, and ild.
Figure 3.1: An illustration of a fuzzy decision.
3.1.2 Tanaka, Okuda, and Asai
Tanaka, Okuda, and Asai suggested an implementation of Bellman and
Za,dchs fuzzy decision using acuts [29].1
1In the literature, acut, alevel, and aset are used synonymously.
19
Definition 3.1 IfC is a fuzzy constraint in X, then an acut of C, denoted by
Ca, is the following crisp set in X:
Ca = {x\fac{x) > a} fora e (0,1]
Ca = cls{a;[/i(:r) > 0} fora 0
(3.2)
where cls{A} denotes the closure of set {A}.
Tanaka, et. al make the following remarkable observation:
sup ij,d(x) = supfc* A max/j,G(x)\.
(3.3)
For illustration, consider Figure 3.2. In case (i), a = ai, and Ca is the interval
between the endpoints of the acut, [(7(aa), C+(ai)]. The maximum p,G in
this interval is shown in the example. In this case, a\ < maxGai /iG(x)\, so [a\ A
maxCai Hg(x)\ = ai In case (ii), a2 > maxca2 Hg{x), so [a2 A max^ fJ,G(x)] =
maxca2 1^g{x) In case (iii), a* = max
picture that a* = supxHd (If not, refer to Figure 3.1 for an illustration of po)
In case (iii), a = a* is also supa[a A maxCa fj,G(x)\. For any a < a*, we have
case (i), where [ai A maxoai iaG(x)\ = i < a*; and for any a > a*, we have
case (ii), where [a2 A maxCc(2 (^g(x)} = maxCa2 /uG(x) < a*.2
This result allows the following reformulation of the fuzzy mathematical
problem:
2The formal proof, which follows the reasoning illustrated here pictorially, is omitted for
brevitys sake. However, the interested reader is referred to Tanaka, Okuda, and Asai [29].
Determine (a*,x*) such that
a*Af(x*) = sup[a A max/(x)]
q. xeca
(3.4)
20
(i)
(>i)
(iii)
Figure 3.2: An illustration of supa[a A max^ Hg{x)} = supx fj,D.
The researchers suggest an iterative algorithm which solves the problem. The
algorithm cycles through a series of steps, each of which brings it closer to
a solution to the relation a* = max/(A). When a* is determined to within a
tolerable degree of uncertainty, it is used to find x* such that f(x*) = max/(X).
c'S
This cumbersome solution method is impractical for largescale problems.
It should be noted that Negoita [19] suggested an alternative, but similarly
complex, algorithm for the flexible programming problem based on Tanakas
findings.
3.1.3 Zimmermann
Just two years after Tanaka, et al. suggested the use of acuts to solve the
fuzzy mathematical problem (FMP), Zimmerman published a linear program
ming equivalent to the acut formulation.
21
Beginning with the crisp programming problem,
minimize Z cx (3.5)
subject to Ax < b
x > 0,
the decision maker introduces flexibility in the constraints, and sets a target
value for the objective function, G:
cx < Z (36)
Ax < b
x > 0.
A linear membership function jii is defined for each flexible right hand side,
including the goal Z as
1
&ij %j
^ ^ 1 Q'ij'X'j ^i
} ' di
0
for aijxj A h,
for Y,j aijxi e (,Ji> h +
for aijxi >bi + dt
(3.7)
where di is the decision makers maximum allowable violation of constraint (or
goal) i.
According to the convention set forth by Bellman and Zadeh, the fuzzy
decision D is characterized by
and
Hd = min aijxj)>
3
(3,8)
max mm
x>0 i
tt(Â£ aijX3)
3
(3.9)
22
is the decision with the highest degree of membership.
The problem of finding the solution is therefore
maximize Hd{%) (3.10)
subject to x > 0.
In the membership functions fa from (3,7), Zimmermann substitutes b[ for
bi/di and B[ for JT aij/di He also drops the 1 (which does not change the
solution to the problem) to obtain the simplification fa = b\ B[x. Equation
(3.10) then becomes
maximize min;(f)' B[x) (311)
subject to x > 0.
This is equivalent to the following linear program:
max a (312)
subject to a < b[ B[x V i
x>0.
Equation (3.12) can easily be solved via the simplex or interior point methods.
See Section 4.2.1 for implementation and results.
3.1.4 Verdegay
The solutions examined so far for the fuzzy programming problem have been
crisp solutions. Ralescu [21] first suggested that a fuzzy problem should have
a fuzzy solution, and Verdegay [30] proposes a method for obtaining a fuzzy
23
solution. Verdegay considers a problem with fuzzy constraints,
maximize z = f(x) (3.13)
subject to a; C C,
where the set of constraints have a membership function /ic, with alphacuts
Ca
Ver degay defines xa as the set of solutions that satisfy constraints Ca. Then
a fuzzy solution to the FMP is
maxz=/(rr) (3.14)
X&Ca
Vae[0,i].
Verdegay proposes solving (3.14) parametrically for a G [0,1] to obtain a fuzzy
solution x, with acut Xa, which yields fuzzy objective value z, with acut za.
It should be noted that the (crisp) solution obtained via Zimmermanns
method corresponds to Verdegays solution in the following way: If the objec
tive function is transformed into a goal, with membership function hg, then
Zimmermanns optimal solution, x*, is equal to Verdegays optimal value x(a)
for the value of a which satisfies
= cx*a) = a. (3.15)
In other words, when a particular cncut of the fuzzy solution, (xa) yields an
objective value (za = cxa) whose membership level for goal G, {poiza)) is equal
to the same a, then that solution Xa corresponds to Zimmermann optimal
solution, x*.
24
I propose that, for continuous membership functions, this optimal a corre
sponds to max(a:x(a!) ^ 0).
Observation 3.2 The value of a which satisfies IXg(z^ cx*) = a is max(a:x(c*:) ^
0).
Proof: The a cut, Xm f x(a) is a solution to Ca A Ga. So
max(aix(a!) ^ 0) = max(a(7a A Ga 0). (3.16)
Bellman and Zadeh [1] (see Section 3.1.1) have shown that
max(o;C,a A Ga^ 0) =sup//jD=5A<5(x). (3.17)
X
Tanaka, et al. [29] (see Section 3.1.2) have shown that
supju# = sup[a A max)UG(x)]. (3.18)
X Ot Ca
They have furthermore shown (as illustrated is Section 3.1.2, figure 3.2) that
the sup occurs when
a = max//G(:r). (3.19)
Ca
Because the objective function z = cx is linear, we have
max /zG(x) = max iaG{z = cx). (3.20)
Ca XCa
Since Zimmermans optimal solution x* is, by definition, the x which maxi
mizes hg over Ca,
max/iG(x) = iaG(x*). (3.21)
Ca
We have the result that the value of a which satisfies iic{z*a = cx*) = a is
max(o!x(a:) ^ 0).
25
3.1.5 Surprise
Recently, Neumaier suggested a way to model fuzzy righthand side values
(Ax < b) based on the concept of surprise [20, 13]. Neumaier defines a sur
prise function, s(x\E), which corresponds to the amount of surprise a variable
x produces, given a statement E. For example, if E is the statement on an
airport arrival screen Flight 347 will arrive at 4:00, and x is the actual arrival
time of Flight 347, one would experience quite a bit of surprise if the plane
arrived at 3:25, but not much surprise if the plane landed at 4:05. Correspond
ingly, s(4:05E) would be small and .s(3:251Â£7) would be large. The range of s
is [0, oo), with s = 0 corresponding to an entirely true statement E, and s = oo
corresponding to an entirely false statement E.
The most plausible values are those values of x for which s(x\E) is minimal,
so Neumaier proposes that the best compromise solution for an optimization
problem with fuzzy goals and constraints can be found by minimizing the sum
of the surprise functions of the goals and constraints. Each fuzzy constraint,
(Ax)i
is translated into a fuzzy equality constraint,
(Ax)i = ii, (3.23)
where the membership function im(Q of is the possibility that > Â£. These
membership functions are subsequently translated into surprise functions by
*(0 = MO"1 1)2; (324)
26
and the contribution of all constraints are added to give the total surprise
= (3.25)
i i
Thus, the fuzzy optimization problem is
minimize ]Tb Sj^Ac)*) (3.26)
subject to x >0.
For an example of a surprise function corresponding with a fuzzy goal, see
Figure 3.3. For triangular and trapezoidal numbers, the surprise function is
simple, smooth, and convex, leading to a tractable nonlinear programming
problem.
Figure 3.3: Derivation of surprise function for (Tx)i < b^umor The member
ship function, AtbTurnor((Tx)i) is on the left, and the surprise function, s(x\(Tx)i <
brumor) is on the right.
3.2 Possibilistic Programming
Fuzzy imprecision arises when elements of a set (for instance, a feasible set)
are members of the set to varying degrees, which are defined by the membership
27
function of the set. Possibilistic imprecision arises when elements of a set (say,
again, a feasible set) are known to exist as either full members or nonmembers,
but whether they are members or nonmembers is known with a varying degree
of certainty, which is defined by the possibility distribution of the set. This pos
sibilistic uncertainty arises from a lack of information. As an example, consider
a hikers accent of Mt. Elbert, from the trail head (elevation 9,600 ft) to its
summit (elevation 14,433). It is known for certainty that the hiker crosses the
10,000 foot elevation line at some point; but because of terrain irregularities,
(rocks, snow, water erosion) the hiker will not know with certainty when he
steps over the line.
One result of this distinction is that fuzziness and possibilistic uncertainty
manifest themselves in different regions of the optimization problem. Given the
basic linear program,
maximize cx (3.27)
subject to Ax < b,
x>0,
fuzziness can occur in the righthandside (b), and/or in the inequality (<);3
and possibilistic values can occur in the objective function coefficients c, in the
constraint matrix coefficients A, and/or in the righthandside b.4 In this section,
3To date, no model has been proposed which deals with fuzzy objective function coefficients
of fuzzy constraint matrix coefficients. However, Lodwick [9] provides an example of a fuzzy
constraint. Suppose a constraint, ayx < bi,i G [l,n], is meant to apply only to members of a
particular fuzzy set, Y. Now suppose that the element represented by row i of the constraint
matrix is a member of set Y to a degree defined by hy{v) Then the constraint i should be
multiplied by the membership value of j/i, resulting in fuzzy coefficients. A similar argument
can be applied for objective function coefficients.
4To date, there is neither semantic nor model for a possibilistic inequality.
28
we examine possibilistic programming formulations that account for possibilistic
uncertainty in A, in b, and in c, and well as in combinations of the three.
3.2.1 Buckley
J.J. Buckley has suggested an algorithm for dealing with possibilistic cost
coefficients, constraint coefficients, and right hand sides. Though Buckley pro
vides no motivating example for the development of this possibilistic optimiza
tion, problems with possibilistic A, 6, and c do arise, as we shall see in Chapter
4. Consider the possibilistic linear program
minimize Z = cx (3.28)
subject to Ax > b, x > 0.
where A = [ay] is an m x n matrix of trapezoidal fuzzy numbers ay =
(aya, Oijp, ay7, aijg), b = (6i,...,6m)* is an m x 1 vector of trapezoidal fuzzy
numbers 6* = (bia, bip, b^n bis), and c = (ci,..., cn) is a 1 x n vector of trapezoidal
fuzzy numbers c* = (cia, Cip, ci7, c^). The fuzzy numbers are the possibility dis
tributions associated with the variables, and place a restriction on the possible
values the variables may assume [33]. For example, Poss[5y = a] = /ia(dy) is the
possibility that ay is equal to a. Stated another way, Poss[ay = a]= x means
that, given the current state of knowledge about the value of ay, we believe that
x is the possibility that variable ay could take on value a.
Because this is a possibilistic LP, the objective function will be governed by
a possibilistic distribution, Poss[Z = z\. Let us consider this simple example as
29
we follow Buckleys derivation of Poss[Z = z]:
minimize 2; = dxi + ex2 (3.29)
subject to fx 1 + gx2 < h
ix 1 + jx2 < 0
xux2 > 0
To derive the possibility function, Poss[Z = z], Buckley first specifies the possi
bility that x satisfies the ith constraint. Let
n(hi, h) = mm(fxan(dn),fxain(din), g,bi{bi)), (3.30)
which is the simultaneous distribution of a^, 1 < j < n, and 6*. In our example
(3.29) this corresponds to
n(/,i?,h) =mm(nF(f),fj,G(g),^H(h)). (3.31)
Then the possibility that x > 0 is feasible with respect to the zth constraint
is:
Poss[a; G = sup(II(aj, 6j)ajX > bi). (3.32)
In our example (3.29) this corresponds to
Poss[a: G Fi\ = sup(II(/,g, h)\fxi + gx2 < h) (3.33)
f,9,h
Poss[a: G =sup(II(i, j, k)\ixi + jx2 < k).
i,j,k
Now the possibility that x > 0 is feasible with respect to all constraints is
Poss[a; G J7] = min (Poss[x G Tj\) (3.34)
1
30
In our example (3.29),
Poss[:r G T\ = min(Poss[a: G Ju], Poss[z G TF^) (3.35)
Buckley next constructs PossfZ = z\x], which is the conditional possibility
that the objective function Z obtains a particular value z, given values for x.
The joint distribution of the possibilistic cost coefficients Cj is
11(c) =min(/xCl(c1)>...,/zCn(cll)). (3.36)
In our example (3.29),
11(c) = min fj,E(e)). (3.37)
Therefore
Poss[Z = z\x] = swpc(Jl(c)\cx = z). (3.38)
In our example (3.29), to obtain the possibility that the objective function
will take on a particular value z, given the solution x, (i.e. Poss[Z = zx]), we
consider all the values of d and e for which dxi + ex2 = z. Of these, the pair with
the highest simultaneous possibility values will yield an objective function value
of z with the same possibility level. For example, consider the case illustrated
in Figure 3.4, when d = (2,2) and e = (3,2) are symmetric triangular fuzzy
numbers, (xi,x2) = (1,1) and z = 7. We consider values of d and e such that
d x 1 + e x 1 = 7. We could select d = 2 and e = 5 with /xd(2) = 1 and
A*jE?(5) = 0. The joint possibility of d = 2 and e = 5 is min(^Â£)(2), /Ge(5)) = 0.
The maximum joint possibility of d and e such that dx 1 + e x 1 = 7 occurs
when d = 3 and e = 4. The joint possibility that d = 3 and e = 4, which is
31
Figure 3.4: Example with d = (2, 2) and e = (3,2).
min(/iÂ£)(3), Â£iÂ£i(4)) = .5, we set equal to the possibility that z = 7 given that
(.Xi,X2) = (1,1). So Poss[z = 7(1,1)] = .5]. A combination of (3.34) and (3.38)
yields the possibility distribution of the objective function:
Poss[Z = z] = sup[min(Poss[Z = z\x\, Possja; E .F])]. (3.39)
x>0
In our example (3.29), to obtain the possibility that the objective function will
take on a particular value z (that is Poss[Z = z\, we consider PossfZ = z\x], as
described above, for all positive ads, and select the x which maximizes Poss[Z =
z\x].
This definition of the possibility distribution motivates Buckleys solution
method. Recall that because we are dealing with a possibilistic problem, the
solution will be governed by a possibilistic distribution. Buckleys method de
pends upon a static a, chosen a priori. The decision maker defines an acceptable
level of uncertainty in the objective outcome, 0 < a < 1. For a given a, we
32
define the left and right endpoints of the acut of a fuzzy number x as rE(o;)
and x+(a), respectively. Using these, Buckley defines a new objective function:
minimizeZ(o:) = c~{a)x (3.40)
subject to A+(a)x > b~{a).
Since both x and at are variables, this is a nonlinear programming problem.
When a is fixed in advance, however, it becomes linear. We can use either the
simplex method or an interior point method to solve for a given value of a. If a
maximal value for a is desired, the LP method must be applied iteratively. For
implementation and results, see Section 4.3.1.
It should be noted that this LP is constrained upon the bestcase scenario.
That is, for a given alevel, each variable is multiplied by the largest possible
coefficient ad.(a), and is required to be greater than the the smallest possible
right hand side b^(a). We should interpret z(a) accordingly. If the solution to
the LP is implemented, the possibility that the objective function will attain the
level z(a) is given by a. Stated another way, the bestcase scenario is that the
objective function attains a value of z(a), and the possibility of the best case
scenario occurring is a.
3.2.2 Tanaka, et al.
In the mid 1980s, Tanaka and Asai [26] and Tanaka, Ichahashi, and Asai
[28] proposed a technique for dealing with ambiguous coefficients and right hand
sides based upon a possibilistic definition of greater than zero. The reader
will note that this approach bears many similarities to the flexible programming
proposed by Tanaka, Okuda, and Asai a decade earlier, which was discussed
33
in Section 3.1.2. Indeed, the 1984 publications refer to fuzzy variables. This
approach has subsequently been classified [6],[8] as possibilistic programming
because the imprecision it represents stems from a lack of information about the
values of the coefficients.
Consider a programming problem with noninteractive possibilistic A, b, and
c, whose possible values are defined by fuzzy matrix A, fuzzy vector b, and fuzzy
vector c, respectively:
minimize cx (341)
subject to Ax
x > 0.
Tanaka, et al. transform the problem in several steps: First, the objective
function is viewed as a goal. As in flexible programming, the goal becomes a
constraint, with the aspiration level for the objective function on the righthand
side of the inequality. Next, a new decision variable Xq is added. Finally, the
bs are moved to the left hand side, so that the possibilistic linear programming
problem is
a&>0 (3.42)
i = 1,2, ...,m,
a/>0
where x' = (1,3:*)* = (1, Â£1,3:2, and a* = (bi,an, ...ain).
Note that all the parameters, A, b, and c are now grouped together in the
new constraint matrix A. Because the objective function(s) have become goals,
34
the cost coefficients, c, are part of the constraint coefficient matrix A. Fur
thermore, because the righthand side values have been moved to the lefthand
side, the b's are also part of the constraint coefficient matrix. The righthand
sides corresponding the former objective functions are the aspiration levels of
the goals.
Each constraint becomes
Yi = dix > 0,
(3.43)
where
^in)
Yi is almost positive, denoted by % > 0, is defined by
Yi > 0 (3.44)
fJYii0) < 1 h,
xlai > 0.
The measure of the nonnegativity of Yt is h: The greater the value of h,
the stronger the meaning of almost positive (see Figure 4.12). Actually, h is
1 a, where a is the level of membership used by Bellman and Zadeh.
Tanaka and Asai [26] developed this theory for triangular fuzzy numbers,
and Tanaka et al. extended it to trapezoidal fuzzy numbers in [28]. Inuiguchi
et al. [6] generalized the approach for LR fuzzy numbers. For the sake of
simplicity, our discussion here will deal with trapezoidal fuzzy numbers, denoted
x = (7, a, (3, S), where (7,5) is the support of x, and (a, /?) is the core of x.
35
0
Figure 3.5: An illustration of 0.
Using (3.44), we can rewrite each constraint from (3.43) as
Mri( 0) = 1
cfix
(i nYx
< 1 h
a\x > 0,
where x > 0. Then (3.45) reduces to
(a* h(a# 7i))Jz > 0.
Since we wish to find the largest h that satisfies these conditions,
maximize h = h*
subject to (ccj h(aci 7j))hc > 0, Vi
h Â£ [0,1].
(3.45)
(3.46)
the LP becomes
(3.47)
36
Since both x and h are variables, this is a nonlinear programming problem.
When h is fixed, it becomes linear. We can use the simplex method or an interior
point algorithm to solve for a given value of h. If the decision maker wishes to
maximize h, the LP method must be applied iteratively. For application and
implementation, see section 4.3.2.
3.2.3 Fuzzy Max
Dubois and Prade [3] suggested that the concept of fuzzy max could be
applied to constraints with fuzzy parameters. The fuzzy max concept was used
to solve possibilistic linear programs with triangular possibilistic coefficients by
Tanaka et al. [27]. Ramik and Rimanek [22] applied the same technique to LR
fuzzy numbers. For consistency, we will discuss the fuzzy max technique with
respect to trapezoidal numbers.
The fuzzy max, illustrated in Figure 3.6, C = max[A,B\, is the extended
maximum operator between real numbers, and defined by the extension principle
(2.6) as
Hc{c)= max min^(o),/x5(6)].
{a,o:c=mai(a,o))
Using fuzzy max, we can define an inequality relation as
A > B
max(^4, B) = A.
(3.48)
(3.49)
37
Figure 3.6: Illustration of fuzzy max, c = max(d, b). Fuzzy max might be
thought of as the rightmost of the two membership functions.
Applying (3.49) to a fuzzy inequality constraint:
f(x,a)>g(x,b) (3.50)
max(/(x, a),g(x, b)) = f(x,a).
Observe that the inequality relation defined by (3.49) yields only a partial or
dering. That is, it is sometimes the case that neither A> B nor B > A holds.
To improve this, Tanaka et al., introduce a level h, corresponding to the decision
38
makers degree of optimism. They define an /ilevel set of A > B as
A >h B
max [A] a > max[jB]a
and
min [A] a > min [Â£]
for all a
This definition of A > B requires that two of Dubois inequalities from Section
2.2, (ii) and (iii), hold at the same time, and is illustrated in Figure 3.7.
Figure 3.7: Illustration of A >h B.
Tanaka, et al. suggest a similar treatment for a fuzzy objective function. A
problem with the single objective function, Maximize z(x, c), becomes a multi
39
objective problem with objective functions:
inf(;z(x,c)a) a G [0,1]
(3.52)
maximize
sup(z(a;, c)a) a G [0,1].
Clearly, since a can assume an infinite number of values, (3.52) has an
infinite number of parameters. Since (3.52) is not tractable, Inuiguchi et al. [6]
suggest using the following approximation using a finite set of G [0,1]:
3.2.4 Jamison and Lodwick
Jamison and Lodwick ([7, 11]) develop a method for dealing with possibilistic
right hand sides that is a possibilistic generalization of the recourse models
in stochastic optimization. Violations of constraints are allowable, at a cost
determined a priori by the decision maker.
Jamison and Lodwick choose the utility (that is, valuation) of a given in
terval of possible values to be its expected average (a concept defined by Yager
Definition 3.3 The expected average (or EA) of a possibilistic distribution of
a is defined to be
It should be noted that the expected average of a crisp value is the value
(3.53)
[31].)
(3.54)
itself. Since a (a) = a^a) = a, EA(a) = \ f^(a + a)dx = a dx = a.
40
Jamison and Lodwick start from the following possibilistic LP:
maximize cTx (3.55)
subject to Ax <6,
x > 0.
By subtracting a penalty term from the objective function, they transform 3.55
into the following possibilistic NLP:
maximize cTx +pT max(0, Ax b) (3.56)
subject to x > 0.
The max is taken componentwise, and each pi < 0 is the cost per unit violation
of the right hand side of constraint i. The utility, (i.e. the expected average), of
the objective function is chosen to be minimized. The possibilistic programming
problem becomes
minimize z = cTx + piL4(max(0, Ax b)) (3.57)
subject to x [0,17].
A closedform objective function (for the purpose of differentiating when
solving) is achieved in ([10]) by replacing
max(0, Ax b) (3.58)
with
Ax b) + e2 + Ax b
2
(3.59)
Jamison and Lodwicks method can be extended, ([7]), to account for possi
bilistic values for A, b, c, and even the penalty coefficients p with the following
formulation:
maximize EAf(x) =
~ Jg1 (c~(a)Tx + c+(o')Tx (p+(a)) max(0, A~(a)x fc+(a:)) (p_(i
max(0, _A+(a:)a; b~(a)))da. (3j
This approach differs significantly from the others weve examined in several
regards. First, many of the approaches weve seen have incorporated the objec
tive function(s) as goals into the constraints in the Bellman and Zadeh tradition.
Jamison and Lodwick, on the other hand, incorporate the constraints into the
objective function. Bellman and Zadeh create a symmetry between constraints
and objective, while Jamison and Lodwick temper the objective with the con
straints. A second distinction of the expected average approach is the nature
of the solution. The other formulations we have examined to this point have
produced either (1) a crisp solution for a particular value of a, (namely, the
maximal value of a), or, (2) a fuzzy/possibilistic solution which encompasses all
possible a values. The Jamison and Lodwick approach provides a crisp solution
via the expected average utility which encompasses all alpha values. This may
be a desirable quality to the decision maker who wants to account for all pos
sibility levels and still reach a crisp solution. Fuzzy and crisp solutions and the
implications for implementation will be further discussed in Chapter 5.
3.2.5 Luhandjula
42
Luhanjulas [15] formulation of the possibilistic mathematical program de
pends upon his concept of more possible values. He first defines a possibility
distribution with respect to constraint F as
II x = (3.61)
where (j,f(u) is the degree to which the constraint F is satisfied when u is the
value assigned to the solution X.
Then the set of more possible values for X, denoted by VP(X), is given by
VJX) = n^maxHx^)). (3.62)
U
In other words, VP(X) contains elements of U which are most compatible with
the restrictions defined by Ux It follows from intuition, and from Luhanjulas
formal proof [15] that when is convex, VP(X) is a realvalued interval, and
when Ylx is strongly convex, VP(X) is a single real number.
Luhandjula considers the mathematical program
maximize cx (3.63)
subject to Ai
x>0.
By replacing the possibilistic numbers c,A;, and b{ with their more possible val
ues, Vp(c), Vp(Ai), and Vp{bi), respectively, Luhandjula arrives at a deterministic
43
equivalent to equation (3.63):
maximize kx
subject to hi G Vp(ci)
Xi < Si
t{ G Vp(dij)
Si G bjj(6j)
x > 0.
This formulation varies significantly from the other approaches considered
thus far: The possibility of each possibilistic component is maximized individu
ally. Other formulations have required that each possibilistic component Cj, Aij,
and bi achieve the same possibility level defined by a.
This formulation also has a distinct disadvantage over the others weve con
sidered. To date there is no proposed computational method for determining
the more possible values, Vj,, so there is no way to solve the deterministic MP.
3.3 Programming With Fuzzy and Possibilistic Components
Sometimes the values of an optimization problems components are am
biguous and the decisionmakers are vague (or flexible) regarding feasibility
requirements. This sections explores a couple of approaches for dealing with
such fuzzy/possibilistic problems. .
3.3.1 Robust Programming
The goal of robust optimization, which has its roots in stochastic optimiza
tion, is to produce a solution whose quality will withstand a wide variety of
parameter realizations. Robust optimization seeks to mitigate the effects of
(3.64)
44
uncertainty rather that merely anticipating it. Hence, robustness reflects a ten
dency to hedge against uncertainty, sacrificing some performance in order to
avoid excessive volatility [17]. Robust fuzzy optimization seeks a solution which
guarantees compliance with constraints at every possibility level.
A fuzzy robust programming problem [18, 3] is a mathematical program
with possibilistic constraint coefficients d^ whose possible values are defined by
fuzzy numbers of the form d^:
max cx (3.65)
subject to d'x' C bi (3.66)
x' = (1, xl)t > 0.
Zadeh [32] defines the setinclusion relation M C N as /xAq(r) < /i^(r)Vr e R.
(Recall that Dubois [2] interprets the setinclusive constraint a^x' C % is inter
preted as a fuzzy extension of the crisp equality constraint. Robust program
ming, however, interprets the setinclusive constraint to mean that the region
in which a^x' can possibly occur is restricted to bi, a region which is tolerable
to the decision maker. Therefore, the left side of (3.66) is possibilistic, and the
right side is fuzzy.
Negoita [18] defines the fuzzy right hand side as follows:
k = {rr > bi}.
(3.67)
45
As a result, we can interpret a(x' C 6* as an extension of an inequality constraint.
The setinclusive constraint (3.66) is reduced to
at(a)x
bi(a)
for all a e (0,1],
If we abide by Negoitas definition of b (3.67), bf = oo for all values of
a, so we can drop the first constraint in (3.68). Nonetheless, we still have an
infinitely (continuum) constrained program, with two constraints for each value
of a 6 (0,1]. Inuiguchi observes [6] that if the lefthand side of the membership
functions for a^, (Hi, din, h are identical for all i, and the righthand side of the
membership functions for aiQ,an, ...ain,bi are identical for all i, the constraints
are reduced to the finite set,
bi,a (3.69)
bifi
di,& bitS
As per our usual notation, (7,6) is the support of the fuzzy number, and
(a, P) is its core. In application, constraint formulation (3.69) has limited utility
because of the narrowly defined sets of memberships functions it admits. For
example, if a,iQ, an,...ain, bt are defined by trapezoidal fuzzy numbers, they must
all have the same spread, and therefore the same slope, on the righthand side;
46
and they must all have the same spread, and therefore the same slope, on the
lefthand side if (3.69) is to be implemented. Recall that in fuzzy robust pro
gramming, the ays are possibilistic, reflecting a lack of information about their
values, and the 6;s are fuzzy, reflecting the decision makers degree of satisfac
tion with their possible values. It is possible that n possibilistic components
and 1 fuzzy component will share identicallyshaped distribution, but it is not
something likely to happen with great frequency.
3.3.2 Delgado, et al.
Delgado, Verdegay, and Villa propose [16] the following formulation for deal
ing with ambiguity in the constraint coefficients and righthand sides, as well as
vagueness in the inequality relationship:
maximize cx , (3.70)
subject to Ax
x>0.
In addition to (3.70), membership functions are defined for the possible val
ues of each possibilistic element of A, membership functions are defined for
the possible values of each possibilistic element of b, and membership function Hi
gives the degree to which the fuzzy constraint i is satisfied. Stated another way,
Hi is the membership function of the fuzzy inequality. The uncertainty in the
aijS and the fr;S is due to ambiguity concerning the actual value of the param
eter, while the uncertainty in the
regarding the necessity of satisfying the constraints in full.
Delgado, et al. do not define the fuzzy inequality, but leave that specification
to decision maker. Any ranking index that preserves ranking of fuzzy numbers
47
when multiplied by a positive scalar is allowed. For instance, one could select
any of Dubois four inequalities from Section 2.2. Once the ranking index is
selected, the problem is solved parametrically, as in Verdegays earlier work [30]
(see Section 3.1.4.)
To illustrate this approach, let us choose Dubois pessimistic inequality (i),
which interprets A > b to mean VT G A, Vy G B,x > y. This is equivalent to
a+ > b~. Then 3.70 becomes
maximize cx (371)
subject to afj(a)x < b^(a)
x>0.
3.3.3 Comments on Fuzzy/Possibilistic Programming
What does the solution to the mixed fuzzy/possibilistic problem mean? In
both the approaches covered in this section, the same a. cuts define the level
of ambiguity in the coefficients, and the level at which the decisionmakers
requirements are satisfied. We interpret the solution to mean that for any value
a G (0,1], we have a possibility a of obtaining a solution that satisfies the
decision maker to degree a.
Using the same a value for both the possibilistic and fuzzy components of
the problem is convenient, but does not necessarily provide a meaningful model
of reality. For example, suppose a = .8. If the solution, x(.8) is implemented,
there is a .8 possibility that it will provide the decision maker with at least .8
degree of satisfaction. On the other hand, if z(.3) is implemented, there is a
.3 possibility that it will provide the decision maker with at least .3 degree of
48
satisfaction. It is apparent that the quality of the solution diverges rapidly. A
more useful approach might be to require a higher level of satisfaction for a
solution with a smaller possibility of attaining feasibility. For example, if I were
to accept a solution that had only .3 possibility of occurring, I might require
that it provide me with at least .9 degree of satisfaction. This is akin to the idea
of being willing to take greater risks for the possibility of greater rewards.
49
4. Application of Possibilistic and Fuzzy Programming Techniques
to the Radiation Therapy Problem
4.1 The Radiation Therapy Planning Problem
The use of particle beams to treat tumors is called the radiation therapy
planning (RTP) problem [5]. Beams of particles, usually photons or electrons,
are oriented at a variety of angles and with varying intensities to deposit radia
tion dose (measured as energy/unit mass) to the tumor. The goal is to deposit
a tumorcidal dose to the tumor while minimizing damage to surrounding non
tumor tissue.
The process begins with the patients computed tomography (CT) scan.
Each CT image is examined to identify and contour the tumor and normal
structures. The image is subsequently vectorized into pixels. Likewise, can
didate beams are discretized into beamlets. This study is concerned with a
twodimensional image of a head tumor at a resolution of 64 x 64 pixels, and
ten beam angles, each divided into 10 beamlets. This is potentially a 642 x 102
problem. However, since not all pixels are in the potential paths of radiation,
and not all beamlets have the potential to deliver dose to the tumor, we set
these to zero a priori and remove them from the analysis. This corresponds to
blocking the beam, which is always done in practice.
The RTP in practice proceeds as follows. After the oncologist delineates the
tumor and critical structures, a candidate set of beam intensities are obtained by
some optimization technique or purely by human choice. These beam intensities
50
are used as inputs to a Federal Drug Administration (FDA) approved dose
calculator to produce the graphical depiction of the dose deposited at each pixel.
A histogram of the percentage of pixels receiving a particular radiation, dose,
called the dose volume histogram (DVH), is one such graphical depiction. In
general, the radiation oncologist never consults the objective function value; just
the final dose distribution is viewed.
A treatment plan is the identification of a set of beam angles and weights
that provides a lethal dose to the tumor cells while sparing healthy tissue, with
a resulting dose distribution acceptable to the radiation oncologist. A dose
transfer matrix A, specific to the patients geometry, represents how a unit of
radiation in beamlet j is deposited in body pixel i. The components of A are
determined by the fraction of pixel i which intersects with beamlet j, attenuated
by the distance of the pixel from the place where the beam enters the body. A
can be divided into a matrix T which contains dose transfer information to tumor
pixels only, matrices C\ through Ck which contain does transfer information to
pixel in critical organs 1 though K, and body matrix B which contains dose
transfer information for all nontumor and noncriticalorgan pixels in the body.
The variable vector x represents the beamlet intensities, and the right hand side
vector b represents the dosage requirements. The vector b can be divided into
btumor j bCl, , bCk and &b0
In the RTP literature, there is no agreement on what the objective function
function of the RTP problem should be [10]. In this paper, the chosen objective
function is minimizing total weighted radiation. Since each element aij of the
attenuation matrix A represents radiation delivered to body pixel i by one unit of
51
radiation from beamlet j, the sum of the elements of the column a* will give the
entire amount of radiation delivered to the body by one unit of radiation from
beamlet j. We therefore select Cj = JT as the objective function coefficient
of Xj.
The crisp formulation of the RTP is
minimize cx (4.1)
subject to B < bb0dy
C, < bCl
Ck < bcK
T <
b
tumor
T <
b
tumor
(4.2)
Two inequalities represent the T = btumor equality constraint. This split
representation has a natural interpretation. T > btumor (which is equivalent
to T < bturnor) represents the requirement that the dose delivered to tumor
pixels is high enough to kill them. T < bturnor represents the requirement that
no body tissue, tumorous or otherwise, be burned.
Uncertainty might occur in four areas of the mathematical programming
problem: the righthand side, the inequality, the A matrix coefficients, and the
cost coefficients. Let us consider what each of these means in the context of the
radiation therapy problem:
Uncertainty in the right hand side implies uncertainty about what levels of
radiation are tumorcidal and what levels are safe for the body and critical
52
organs.
Uncertainty in the inequality implies the oncologists flexibility regarding
the satisfaction of the dosage requirements.
Uncertainty in coefficient ay implies uncertainty about the degree of radi
ation beamlet j imparts to pixel i. One (possibilistic) source of this type
of uncertainty is patient positioning and breathing. 1 In addition, we have
possibilistic uncertainty in the A matrix because the attenuation matrix
is a linearizion of the nonlinear attenuation process. Finally, some uncer
tainty in the a^s results from the fuzzy nature of the pixels themselves.
Biologically, cells are not strictly tumor or nontumor cells. Some cells
take on partial tumor characteristics, which results in a fuzzy membership
function.
Uncertainty in cost coefficient Cj implies ambiguity about the impact of
the intensity of beamlet j on the overall radiation to the body.
4.2 Fuzzy RTP Problem
4.2.1 Flexible Programming
Recall that in flexible programming, the decision maker relaxes the param
eters of the constraints and the goals in order to achieve a more desirable result
(or sometimes simply to obtain a feasible result.) Flexible programming can
1 Ambiguity caused by patient positioning and breathing will not be noninteractive from,
one body pixel to another, since the body is positioned and breathes as a whole. Since all
the fuzzy/possibilistic formulations in this study make the assumption of noninteractivity, we
have a less than ideal modeling situation.
53
be interpreted as fuzziness in the inequality relationship. Algorithms for flexi
ble programming were proposed by Tanaka, et al. [29], by Negoita and Sularia
[19], and by Zimmermann [34]. We implement Zimmermanns solution here be
cause he was the first to propose a linear programming equivalent to the flexible
programming problem. Given the following fuzzy programming problem,
cx < Z (43)
Ax < b
x>0,
Zimmermann calls for the solution of the crisp linear program:
maximize a (44)
subject to a < b[ B[x V i
x>0,
where di is the maximum tolerable violation of constraint or goal i, 5' = 6,/d,,
and B[ = ]>T a^/dj.
Flexible programming in the radiation therapy application represents will
ingness on the part of the radiation oncologist to exceed the radiation levels
deemed safe for the body, critical organs, or tumor, or to fall short of a level be
lieved to be lethal for the tumor. With this in mind, each right hand side of the
RTP problem is redefined as a trapezoidal fuzzy number: biJ0(jy = (0,30,0,35),
bcritical = (0,20,0,25), and btumor = (56,60,60,64). Zimmermanns formulation
of the RTP problem defines no goals, only constraints.
This formulation was programmed in MATLAB and solved with the opti
mization toolbox. MATLAB calculated a solution to Zimmermanns problem
54
in microseconds, with an optimal a* = .22, and optimal x* represented in the
following dose volume histogram (figure 4.1). This DVH leaves something to be
Figure 4.1: Dose volume histogram for RTP problem solved via Zimmer
manns fuzzy technique. Optimal a* = .22.
desired, as 85% of the tumor region recieves significantly less that the ideal 60
Gy radiation level.
4.2.2 Verdegay
Recall that Verdegay proposes a parametric solution to the fuzzy problem
mayiz f{x) (4.5)
XSCqc
VaG [0,1],
55
where Ca is the a cut of the constraint set C. For the radiation therapy
problem, Verdegays method looks like this:2 *
minimize E* <46)
i
subject to B < bb0dy(a)
Ci < bCl(a)
Cn < bcn(a)
T E btumor{pi)
T E btumor((x)
V a Â£ [0,1]
Rather than solve parametrically to obtain a solution into which specific
values of a may be inserted, I wrote a MATLAB program which accepts an a
value as input and returns a crisp solution. The fuzzy solution is obtained by
running the program with multiple a values. This program is extremely fast
executing in under a second. The maximum a value of the fuzzy solution is 0.22.
As expected, this corresponds with the results for the Zimmermann method,
which maximizes a. Indeed, the Verdegay solution at a = 0.22 is exactly the
same as the Zimmermann solution. The reader will recall from section 3.1.4 that
this is the avalue for which Verdegays solution and Zimmermanns solutions
were shown to respond, so this is the expected result.
2The objective function in 4.6 is to minimize overall radiation, i.e. the cvector is all ones,
so z = cx = x.
56
The dose volume histograms below show snapshots of the fuzzy solution at
cl = 0 (Figure 4.2) and a = 0.22 (Figure 4.3). Notice that the the DVH line
become steeper as a increases this is because the constraints are more strictly
satisfied.
Vcrdegay Fuzzy Optimization 64x64 10 beams, alpha = 0
Figure 4.2: Dose volume histogram for RTP problem solved via Verdegays
fuzzy technique. Snapshot taken at a = 0.
4.2.3 Surprise
The surprise approach handles fuzziness in the righthand sides of the goals
and constraints. Semantically, this means that radiation oncologist is vague
about the maximum acceptable dosage for the body and critical organs, and
about the dose required to kill the tumor. Lodwick et al. [13] and Lodwick
and Bachman [10] previously implemented the surprise approach for the RTP
problem.
To define appropriate surprise functions for the righthand sides, let us re
view the DVH. Figure 4.4 shows the perfect tumor DVH, a typical DVH, and
57
Verdegay Fuzzy Optimization 64x64 10 beams, alpha = 0.22
Figure 4.3: Dose volume histogram for RTP problem solved via Verdegays
fuzzy technique. Snapshot at a = .22.
the differences between them that we wish to penalize.
The surprise function for the tumor right hand side, represented in Figure
3.3, is
s(x\(Tx)i = bTumor) = (7t l)2 (4.7)
Oi ICC)
wher ea = ftw((Ti)i).
The surprise functions for the other right hand sides are similarly defined.
The fuzzy optimization using surprise functions is:
total pixels N
minimize z = J2 SiiY^OijXj) (4.8)
i=1 j=1
subject to x > 0.
58
Figure 4.4: Dose Volume Histogram for tumor dosage. The heavy line repre
sents ideal dosage distribution, dotted line represents typical dosage distribution,
and hatched areas represent areas the surprise function should penalize.
The MATLAB implementation required the addition of a hard constraint to 4.8
to require that minimum dosage is delivered to the tumor pixels. That is,
total pixels N
minimize z = (49)
i= 1 j1
Subject tO YsjLl Tijxj > train
X >0.
The DVH obtained from the implementation of (4.9) appears below (Figure
4.5). The tumor dosage in this DVH is particularly nice it is almost square
at the 60 Gy dosage level. This quality, however, comes at the expense of the
other constraints. Critical organs receive a much higher overall dosage than they
did in solutions calculated by the flexible programming techniques. The body
dosage exceeds 50 Gy in some places 20% higher than the maximum allowable
59
dosage. The surprise approach is more flexible in its handling of soft constraints
Fuzzy Optimization with surprise 64x64 10 beams
Figure 4.5: Dose volume histogram obtained with fuzzy surprise methodology.
than traditional flexible programming since it does not require that every soft
constraint be satisfied at the same (maximal) alevel, but minimizes the sum of
constraint violations, allowing a tradeoff between constraint violations. As we
have seen, however, the method can sometimes afford too much flexibility: in
our RTP implementation, we had to add hard constraints in order to obtain a
reasonable DVH. A fuzzy methodology that requires the introduction of crisp
constraints in order to achieve a viable solution is clearly disadvantageous.
4.3 Possibilistic RTP Problem
4.3.1 Buckley
60
Recall Buckleys LP equivalent to the possibilistic minimization problem:
minimize Z = c~x, (410)
subject toA+x > b~
x > 0.
In the previous chapter, we dealt with possibilistic righthand sides and
fuzzy inequalities. In order to validate the model, I first programmed Buckleys
method with crisp c and A (model BuckelyB3), and possibilistic b, to simulate the
Zimmermann program used in the previous section. I then added uncertainty by
using possibilistic A and b and crisp c (model BuckleyAB.) Lastly, I programmed
Buckleys method with possibilistic A, b, and c (model BuckleyABC.)
As expected,the BuckleyB model with crisp A and c and possibilistic b
agrees with the Zimmermann model. Recall that Zimmermann maximizes a,
while BuckleyB receives a as an input parameter. The maximum value attained
by a in the Zimmermann model, to two digits, is 0.22. When 0.22 is the a
value input for BuckleyB, the results are identical (See Figures 4.6 and 4.1.)
Furthermore, the BuckleyB model is infeasible for a > 0.22. This agrees with
the Zimmerman model, which maximizes a. As a decreases, the quality of our
solution increases (for an example, see Figure 4.7, where a = 0.) Consider, for
example, the C2 line, which represents dosage to critical organ 2. When a = 0,
no part of the C2 organ is receiving more than 10 units of radiation. When
a = .22, however, 40% of the C2 organ receives more than 10 units of radiation.
Recall that Buckley provides the most optimistic solution, so as the level of
3A11 MATLAB code is available from the author upon request.
61
Buckleyfuzzy B 64x64 10 beams, alpha .22
Figure 4.6: Dose volume histogram for RTP problem with possibilistic right
hand side, Buckleys method, a = 0.22.
uncertainty, represented by a, increases, the quality of the bestcase scenario
increases. Put another way, the greater the risk, the larger the potential reward.
The DVH for a = 0.22 (Figure 4.6) represents an outcome z, whose possibility
of occurring, once the solution to BuckelyB is implemented, is only 0.22. Now
consider what happens to the solution when we add uncertainty in A. When
the elements of A are defined by a possibility distribution rather than a crisp
number, we change the ith constraint of the problem from
a,iiXi + ... + ainxn < bf(a) (411)
to
a7i(a)xi + +a,n{a)xn < bf(a). (4.12)
Because there is now variability on both sides of the inequality, it is easier to
satisfy the constraint. Since the possibility that x is feasible increases, Poss[Z =
62
Buckleyfuzzy B 64x64 10 beams, alpha = 0
Figure 4.7: Dose volume histogram for RTP problem with possibilistic right
hand side, Buckleys method a = 0.
z] increases (see (3.34) and (3.39)). Accordingly, we expect a higher maximum
a level. Indeed, the maximum feasible a in BuckleyAB is 0.70, whereas the
maximum feasible a. in BuckleyB was 0.22. Additionally, the resultant DVH
from the maximum a in BuckleyAB indicates a better clinical result than the
DVH which results from the maximum a in BuckleyB To see this, compare the
BuckelyB result in Figure 4.8 with the BuckleyAB results in Figure 4.6. The
tumor dosage is similar in the two DVHs but the critical organ and body doses
prescribed by the BuckleyAB model are lower.
Some qualities of the potential bestcase scenario solution increase as a
decreases. (See Figures 4.9 4.11.) The total dosage delivered to the body
and critical organs decrease as a decreases. But consider what happens to the
tumor dosage. When a = .7, 55 Gy is the smallest dosage any part of the
tumor receives. But when a = 0, the maximum tumor dosage is less than 55
63
BuckleyAB 64x64 10 beams, alpha .70
Figure 4.8: Dose volume histogram for RTP problem with possibilistic right
hand side and A matrix, a = 0.70.
Gy, with most of the tumor receiving less than 60 Gy. Recall that the crisp
constraint requires a tumor dosage of 60 Gy, and the fuzzy triangular number
which represents the tumor dosage is (60,56,64). In other words, at a = 0,
we require tumor dosage to exceed 56 Gy. This begs the question, how have
we arrived at a DVH which displays a tumor dosage of 50 Gy? Recall that we
are dealing with possibilistic right hand sides, and possibilistic coefficients on he
lefthand side. Buckleys formulation selects the lower end of the left hand side
acut, and requires it to be larger than the upper end of the right hand side
acut. That is, Aj(a)x < b'Ka). This generous formulation describes a feasible
area larger than its crisp counterparts. To evaluate the optimal solution selected
from this generous feasible region, we multiply by the crisp attenuation matrix.
This crisp evaluation of a possibilistic solution produces puzzling results.
64
BuckleyAB 64x64 10 beams, alpha .50
Figure 4.9: Dose volume histogram for RTP problem with possibilistic right
hand side and A matrix, a = 0.50.
Now consider adding uncertainty to c. Recall that Cj is defined to be JT a^.
When the a^s are possibilistic, as in the BuckleyAB model, the Cj s are sums of
possibilistic numbers, and have an inherent possibility distribution built in. It
is not surprising then, that the results of model BuckleyABC, with possibilistic
parameters A, b, and c are identical to the results of BuckelyAB, with only
possibilistic A and possibilistic b. The dose volume histograms are omitted here
because they are redundant.
The Buckley technique requires that all parameters carry the same degree
of uncertainty. Since uncertainty in the A matrix represents something entirely
different from uncertainty in the b vector (variation in the degree of radiation
delivered to a pixel versus uncertainty about what level of radiations constitutes
a lethal dose), this is a significant limitation of the models. It is conceivable, even
likely, that the decision maker would choose to account for a different level of
65
BuckleyAB 64x64 10 beams, alpha = .20
Figure 4.10: Dose volume histogram for RTP problem with possibilistic right
hand side and A matrix, a = 0.20.
vagueness in the dosage limitations (represented by the bvector) than the level
of ambiguity in the attenuation matrix (represented by the Avector). Another
potential disadvantage of Buckleys technique is that it accepts a as a parameter.
If the goal is to find the solution with the lowest level of uncertainty, the Buckley
method is at a disadvantage finding a maximum a requires tedious trialand
error. On the other hand, the Buckley method offers a significant advantage to
the decision maker who wishes his solution to reflect the flexibility allowed by a
particular level of uncertainty.
4.3.2 hLevel Technique
The hlevel technique, proposed by Tanaka and Asai [26] and Tanaka, et al.
[28] solves the optimization problem,
> 0, i = l,2, ...,m, x' > 0, (413)
66
BuckleyAB 64x64 10 beams, alpha = .00
Figure 4.11: Dose volume histogram for RTP problem with possibilistic right
hand side and A matrix, a = 0.00.
where x' = (i,x*), t = (1, x\, x%, ..., xny, and a* = (Â£>*, da, ...din). In this for
mulation, the 6s have been moved to the right hand side, and any objectives
have been reframed as goals and appear with the objectives. So the possibilistic
matrix A captures ambiguity in A, 6, and ct.
The resultant nonlinear program for fuzzy trapezoidal coefficients =
(7,a,P,5) is
maximize h (414)
subject to
(ctj h(ai 7i))Ta: > 0, i = (0,..., m)
he [0,1] (4.15)
We solve this as iterative linear programs, each with a different fixed h. The
result for the 64 pixel, 10 beam head problem is h* = 0.21, which means that
67
the optimal solution, when implemented, has a 0.21 possibility of meeting the
feasibility requirements set forth by the radiation oncologist. The DVH shows
the resultant beamlet intensities multiplied by a crisp attenuation matrix. This
is a particularly attractive DVH because the tumor dose distribution is almost
square.
hlevel 64x64 10 beams, h = .21
Figure 4.12: Dose volume histogram for RTP problem solved via Tanakas
hlevel technique, h* = .21.
It should be noted that the right hand sides in this formulation have a mean
ing different from that of the right hand sides when the problem was formulated
for Zimmermanns flexible programming technique in Section 3.1.3. In Zimmer
mans fuzzy problem, the membership function of an element of the righthand
side, represents the radiation oncologists subjective degree of satis
faction with radiation at level x. In Tanaka, et aVs possibilistic problem, the
possibility distribution of bi (Poss(Bi = x)) represents the objective possibility
that the maximal tolerable radiation is actually x.
68
4.3.3 Fuzzy Max Technique
Recall that the fuzzy max definition of A >h B is
A >h B (4.16)
maxfAjo, > maxfRjo;
and
min[A]a > min[B]&
for all a [h, 1].
In practice, this means that a fuzzy mathematical .program with n con
straints is solved with 2n crisp inequalities, one for the upper bound of each
acut, and one for the lower bound.
We have heretofore dealt with crisp equivalents to the fuzzy or possibilistic
RTP problem in the form
B bbody
Ci bCl
Ck < bcK
T btumor
T btumor
Remember that the tumor constraints appear twice because the dosage must
be both high enough to kill the tumor and low enough to prevent burning the
patient.
69
Using Tanaka et aVs formulation with the level /iwhich represents the
decision makers certainty requirementsprovided a priori [27], the fuzzy max
version of the fuzzy inequality yields the following crisp form:
1 to 1 bBh
eft
n A
Th <
Bn bBh
Cin bclh
Cfch bckh
Tn bTh
,~bTh.
Note that the four constraints on T now require equality at both the alpha
cuts, that, is and = 6^. This is an unlikely scenario, particularly
given the arbitrary nature of membership functions and possibility distributions.
As expected, the implementation of this formulation was infeasible for all values
of h because it is so heavily constrained.
To relax this constraint a bit, we replace Bellman and Zadehs equality
definition (Dubois equality definition (v)),
= ^g(b,x) > (417)
70
with Dubois strong equality definition (vi),
f(a,x)Cg(b,x). (4.18)
In practice, this meant reducing the four tumor constraints to just two,
(4.19)
which resulted in a feasible problem.
When the possibility distributions for the elements of the constraint ma
trix A are triangular, with a spread of 10% of the value of a^, (e.g. dij =
(0.9ay,Oy, 1.1 ay),) the FLP is only feasible for h < .01. Decreasing the spread
to 5% results in feasible solutions of h < .05, as represented by Figure 4.13.
(Compare this to Buckley: a = .71 for a spread of 10%. This fuzzy max format
is very unforgiving.) The DVH is very attractive. The average tumor dose is
closer to 60 Gy than either the hlevel or the Zimmermann results. Also, the
body and critical organ dosage are lower overall than in the Zimmermann of
/ilevel solutions. One byproduct requirement of the hlevel technique is that if
the constraint is met at bounds (that is, for nonbasic variables), then the spread
of A must be at least as big as the spread of b, as illustrated is Figure 4.14.
This might be useful to a decisionmaker who wishes to ensure similar possibility
distributions on both sides of the inequality. This implementation is slower than
some of the other programs, taking up to 330 seconds for MATLAB to execute.
This is due, in part, to the fact that the fuzzy max formulation requires twice
as many constraints as any of the other formulations implemented to this point.
Running time might also be affected by MATLAB memory allocation problems.
1 I < bTh
~Tk 1 1 cr l
71
1
Fuzzy Max 64x64 10 beams, h .05
Figure 4.13: Dose volume histogram for RTP problem solved via the fuzzy
max technique, h* = .05.
4.3.4 JamisonLodwick Penalization Method
The Jamison and Lodwick formulation for the RTP problem with possibilis
tic right hand side is as follows:
min z = cTx + pbEA ^max ^0, Bx bbody^ ^ +
Ef=lPCiEA (max (0, Qx q)) \pTEA Tx tj (4.20)
subject to 0 < x < U,
Let the body and critical organ right handside values be defined a trape
zoidal fuzzy numbers bbody = (0,35,0,40) and = (0,25,0,30) respec
tively. Let the tumor tolerances be defined by the triangular fuzzy number
bx = (60,56,64). Semantically, this means that the best evidence from research
or empirical experience indicates that as long as the body pixels receive 35 units
72
Figure 4.14: In the figure on the left, A
or less, the radiation oncologist is 100% satisfied with dosage at the body pixels.
However, body pixels that receive more than 40 units axe completely unaccept
able to the oncologist. The same interpretation holds for the critical organs. For
the tumor pixels, the radiation is to be held within the range of 56 and 64 with
60 being the ideal dosage.
Thus, the left and right membership functions are, respectively, for the body:
Kody (a) = 0
Hedy M = 40 5q:
The left and right membership functions for the critical organs are
bCi (a) = 0
&J,. (ck) = 30 5a,
and the left and right membership functions for the tumor are
73
bT (a) = 56 + 4a
bf (a) = 64 4a.
The penalty functions for the Jamison and Lodwick method associated with
the RTP are as follows, where p > 0 are vector weights.
i. Body:
PbodyEA (max (o, Bx fy^)) =
\pbody (/o max (o, b^ody (a) Bx) da + ff max (o, Bx bjody (a)) da) =
\Vbody Jo max (0, Bx (40 5a)) da,
(4.21)
since b^ody (a) Bx = Bx < 0.
ii. Critical organs
PCiEA (max (o, Qx 5^)) =
\pCi (fo max (0,65. (a) Cix) da + ff max (0, CiX bÂ£. (a)) da) =
\pci fo max (0, Qx (30 5a)) da,
(4.22)
since bq. (a) CiX = CiX < 0.
iii. Tumor
PtEA (Tx bx )
\px (fo max (0, b^ (a) Tx) da + max (o,Tx bf (a)) da) =
Px (fo max (0,56 + 4a Tx) da + fJ,1 max (0, Tx (64 4a)) da) .
(4.23)
74
Using the following substitution for the maximum in the above integrals gives a
closed form of the integral:
max(0, x) = i{Vx2+ E + x}
where E > 0 is very small.
(4.24)
Thus, the penalty, b(x), for the body pixels is,
1 f1
b(x) Pbody / max (o, Bx 40 + 5a) da
2 Jo
= ^Pbody J  V(Bx 40 + 5a)2+ E + Bx 40 + 5cc j da
= ~Pbody ^Bx 37.5 + J \/25a2 + 10(Bx 40)a + (Bx 40)2+ Eda
= ^Pbody  Bx 37.5 + 5 J <\Ja2 + y (Bx 40)a + ^ [(Bx 40)2+ E]da
= Pbody \Bx 37.5 + 5
f ^
Jo
a2 + &ia + Cirfa j
^Pbody ^Bx 37.5 + 5 y/l + b\ + Ci + &i^ ^biy/ci+
2Cl^ 0n In ^2 + + 2\/l + 6x + Cx
where,
(4.25)
since
&i = (.Ba; 40)
5
Ci=^(Bx40)2+^'
/ \J a2 + bia + ci da
Jo
\+ i61) l6iv/5T
+ fg&i ~Cll In ^2 + f>i + 2\/l + 6i + c\
(4.26)
(4.27)
75
The critical organ and tumor constriants are derived similarly. The possibilistic
RTP becomes:
K
mm2: = cTx + b(x) + '^2/Ci{x) + t\(x) + t2(x) (4.28)
i=1
subject to 0 < x < U.
This convex nonlinear program takes significantly longer to execute than
many of the linear programming equivalents we have examined so far. For
instance, using MATLAB to solve the head problem with 64 x 64 resolution
takes more than 15 minutes when the Jamison and Lodwick method is used. In
contrast, the same problem with Zimmermanns flexible programming algorithm
in MATLAB takes less than a second.
A dose volume histogram of a solution obtained via Jamison and Lodwicks
penalization method is shown in Figure 4.15. The attractive solution shows a
very consistent tumor dosage centered around 60 Gy, as evidenced by the near
vertical tumor dose line in the DVH. Body and critical organ dosages are higher
than in the DVHs of other techniques weve used. A particular problem is the
body dosage, which exceeds 50 Gy.
A notable feature of the Jamison and Lodwick method is that it imposes no
hard boundaries for constraints (like when a = 0 in some of the other methods
weve seen.) Constraint violations are penalized, but could be huge if the penalty
is small enough. Further experimentation should explore the impact of other
penalty values on the solution quality.
4.3.5 Robust Programming
Recall that robust programming seeks a solution that satisfies the con
straints at all crlevels by solving the following problem with possibilistic A
76
Jamison and Lodwick Possibilistic Optimization 64x64 10 beams
Figure 4.15: Dose volume histogram for Jamison Lodwick method.
matrix and fuzzy inequalities represented by a membership function for the b
vector:
maximize cx (4.29)
subject to af{a)x < b*(a)
a(a)x>b~(a)
for all a G (0,1].
Since a can take on a infinite number of values, the program is infinitely con
strained. Recall that Inuiguchi suggests approximating the infinitely constrained
program by solving for a discrete set of {cq} G (0,1]. For this implementation,
the interval (0,1] is divided into a variable number of equalsized intervals, and
the constraints in (4.29) applied for a at the right endpoint of each subinterval.
For example, four values for a result in the following formulation for the RTP:
77
max cx
(4.30)
subject to al(.2h)x < &d(.25)
ar(.25)a:>6r(.25)
a+(.50)i
ar(.50)a;>6r(.50)
af(.75)x
ar(.7b)x>b(.75)
af(l)x
aT(l)x>bi(l)
(4.31)
)
It should be evident that if the crisp problem has n constraints, the robust
formulation has 2rn constraints, where r, or resolution, is the number of dis
crete a values considered. Large r values, which more closely approximate the
infinitely constrained robust program, are computationally costly.
78
For the radiation therapy problem, the robust implementation is
B+(a)x < b+ody(a)
B~(ot)x > bfody(a)
Cf(a)x
C(a)x>bc.(a)
T+(a)x < b+umor(a)
T (ot)x ^ b^UTnor{(x)
tumor
(4.32)
(4.33)
Each element of the A matrix is given a triangular possibility distribution with a
spread of ten percent of its value. That is % = (%, .9%, 1.1%). Semantically,
this means the the decision maker is confident that the level of radiation delivered
to pixel % by beamlet j lies somewhere in the interval (.9%, 1.1%), and that.
Poss[% = %] = 1. In previous formulations of the RTP problem involving fuzzy
goals (see Sections 4.2.1, 4.2.2, and 4.2.3) we have represented the right hand
sides as the following fuzzy numbers:
Because robust programming requires feasibility at all possibility levels, includ
ing a = 1, it is very restrictive. For instance, using btumor = (60,60,56,64),
results in the following requirements, (by substituting a = 1 into the last two
bbody (0,35,0,40)
bcritical ~ (0, 25, 0, 30)
btumor = (60,60,56,64).
tumor
79
lines of 4.32):
%x < 60
T{X > 60
Vi.
As expected, the MP proved to be infeasible with these parameters. For the
purpose of investigating the robust formulation, the parameters were relaxed as
follows to create a feasible problem:
hody = (0, 40, 0, 50)
bcritical = (0, 30, 0, 40)
humor (55,65,45, 75).
MATLAB solved the problem for the resolution r = 1 in 11 seconds, for
r = 2 in 23 seconds, for r = 3 in 34 seconds, and for r = 4 in 103 seconds.
It could not solve the problem for r = 5. Since the increase is linear for the
first three resolution values, it is likely that MATLAB encountered a memory
allocation problem for larger problems.
Regardless of the resolution, the robust program produced the same solution
with the given parameters. This is not surprising, since a 1 should be the
most restrictive a level, and a = 1 occurs no matter the resolution. It is not
necessarily the case, however, that no r over 1 need be calculated, if Ax has a
bigger spread than b, it will not be contained in b for small a values.
Figure 4.16 is a representative DVH for the robust solution. It is very similar
the DVH produced by the fuzzy max technique, which in turn had superior
80
Robust Optimization 64x64 10 beams, r = 3
Figure 4.16: Dose volume histogram for robust formulation, with a resolution
of 3.
results to Zimmerman and the hlevel technique.
4.4 The Crisp RTP Problem
As reference point for fuzzy and possibilistic optimization methodologies,
an implementation of the crisp RTP problem is included. Recall the crisp RTP:
minimize cx (4.34)
subject to Bx < bsody
C\x < bCl
CiX < bCi
Tx 5: bxumor
Tx ^ b^umor
X >0.
81
As in other implementations addressed by this paper, A is the attenuation
matrix, and c is the column sum of A so that the objective function is to minimize
total (weighted) radiation. Previously, b was defined as the following trapezoidal
fuzzy numbers: for the body, (0,35,0,40), for the critical organs, (0,25,0,30),
for the tumor, (60,60,56,64). The crisp formulation is infeasible if the core
values (corresponding to a = 1) of b are assigned to b. When the support
values (corresponding to a = 0) of b are assigned to b, the crisp formulation
produces the solution represented in Figure 4.17. The crisp solution does not
vary significantly from the robust and fuzzy max solutions, which were among
the best.
Crisp Optimization 64x64 lO beams
Figure 4.17: Dose volume histogram for the crisp formulation.
5. Comparison of Optimization of Uncertainty and Techniques
5.1 Formulation Selection
As we have seen, fuzzy and possibilistic optimization formulations are abun
dant. The choice of formulation depends firstly upon the semantics of the prob
lem, secondly upon the desired format of the outcome, and the lastly on the
performance of a particular algorithm. What follows are several dichotomies
which should aid a decision maker in the selection of an appropriate algorithm.
5.1.1 Hard or Soft Constraints
If the constraints in question are required to be met exactly as stated, the
decision maker should select a formulation that uses hard constraints, such as
Tanakas hlevel technique, the fuzzy max technique, Neumaiers surprise tech
nique, or Buckleys technique. On the other hand, if the constraints indicate
preference, are flexible, or can be violated at a cost in order to obtain a more
desirable outcome, a formulation with soft constraints is appropriate. Methods
we have examined with soft constraints are Zimmermanns, Verdegays, Jami
son and Lodwicks, Luhandjulas, Delgados, and the fuzzy robust approach. It
should be noted that in the case of hard constraints, there may not be feasible
solutions in practice.
5.1.2 Fuzzy or Possibilistic Parameters
Another determinant of the appropriate formulation is the nature of the
problem parameters. If the problem contains parameters whose values are known
83
with precision, and the uncertainty lies in the inequality, a problem with crisp
parameters, such as Zimmermanns or Verdegays should be chosen. If the prob
lem contains, parameters which are imprecise because there is not sufficient in
formation to fully determine their values, a model that incorporates possibilistic
parameters should be chosen, possibly from the following: Luhandjulas tech
nique, the fuzzy max technique, Tanakas hlevel technique, Buckleys technique,
fuzzy robust optimization Jamison and Lodwicks expected average technique,
or Delgados technique. If, however, the value of some parameter is not fully
defined because the parameter assumes multiple values to various degrees, the
decision make should choose a formulation with fuzzy parameters. The only
model fully examined in this paper which incorporates fuzzy parameters is Neu
maiers surprise approach.
5.1.3 Feasibility or Optimality
Model selection will also depend on the decision makers goal. The problem
might be to find a feasible solution (often in these cases, a crisp program is in
feasible). In this case, the decision maker should select a model which seeks the
most feasible solution, such as Zimmermanns technique, fuzzy robust optimiza
tion, surprise, Tanakas hlevel technique, or Delgados technique. If there is an
objective function to be optimized, the decision maker might choose to reformu
late the objective function as an additional constraint, as in goal programming.
In this case, the problem again becomes one of finding the most feasible solu
tion, and'a formulation chosen from those listed above is advised. On the other
hand, the decision maker might desire to optimize the objective function for a
given feasibility level, choosing a model which seeks an optimal solution such
84
as Buckleys, Luhandjulas, Verdegays, or the fuzzy max model. The Jamison
and Lodwick approach alone among those implemented in this paper affords the
decision maker the ability to trade off between optimality and feasibility.
5.1.4 Fuzzy or Crisp Solution
A formulation might also be selected because a crisp solution or a fuzzy/possibilistic
solution is desired. In general, algorithms that require a input give fuzzy so
lutions, and those which do not take an a level a priori give crisp solutions.
Jamison and Lodwick differ greatly from other formulations, by providing a
crisp solution that is not given for a particular a value, but is averaged over all
a values.
The question to consider is when to defuzzify. At some point, all solutions
must become crisp, because one cannot implement a fuzzy or possibilistic solu
tion. For instance, in a transportation problem, the decisionmaker might send
the order to ship 3 = (3,2,4) widgets, but in the end, the dockworker will load
a crisp 2, a crisp 3, or a crisp 4 on the truck. The tradeoff is generally between
flexibility and control. A decisionmaker who selects a formulation which yields
a crisp solution knows what solution will be implemented, for there is only one.
A decisionmaker who selects a formulation with a fuzzy or possibilistic solution
must at some point choose (or allow another, like the dockworker in the example
above, to choose) a crisp instantiation of the solution to implement. Further
more, this choice must be made without the aid of an optimization program.
Leaving the defuzzifying to a person rather than a program grants flexibility,
but relinquishes control over the quality of the solution.
85
When choosing a fuzzy or crisp solution, it is important to keep the semantics
of the problem in mind. As an example, consider Verdegays fuzzy optimization
technique and Buckleys possibilistic approach, in the case that b is possibilistic
and A and c are crisp. The two programs give noncrisp solutions that are
calculated in exactly the same way, but semantically means very different things.
It makes sense to give an a priori a level for Buckley, because it corresponds
with the decision makers level of comfort with risk. In Verdegay, the a level
corresponds to the decision makers level of constraint satisfaction. Would the
decision maker not always wish to maximize the level of constraint satisfaction,
which corresponds to using Zimmermanns method?
Bearing this in mind, the following formulations produce noncrisp solu
tions: Verdegays, Buckleys, fuzzy max, and Delgados. For a crisp solution,
good choices are Zimmermanns, surprise, hlevel, Luhandjulas, Jamison and
Lodwicks, and fuzzy robust optimization.
5.1.5 Complexity and Speed
If multiple formulations fit the semantics of the problem and meet the needs
of the decision maker, efficiency may be considered a deciding factor. The list is
a ranking, from fastest to slowest, of the techniques tested in this paper: crisp,
Verdegays, Zimmermanns, Buckleys, robust, A level, fuzzy max, Jamison and
Lodwicks.
5.2 Additional Thoughts
The focus of this chapter has been the selection of a fuzzy or possibilistic
optimization algorithm. The broader question is whether a fuzzy/possibilistic
algorithm should be selected at all.
86
For instance, in the radiation therapy problem, the solution to the problem is
multiplied by a crisp attenuation matrix to calculate delivered dose and produce
the dose volume histogram. A possibilistic solution, calculated with possibilistic
A, might be feasiblei.e. produce an acceptable dose when multiplied by the
possibilistic A, but be infeasible when multiplied by a crisp A. For an example
of this, refer to the Buckley implementation is Section 4.3.1. Since the DVH
is the tool used by radiation oncologists to judge the quality of the solution,
this is a problem. To accurately display the possibilistic solution, the DVH
should be calculated using a possibilistic A. Even if the radiation oncologist
employed a possibilistic DVH, the Federal Drug Administration (FDA) regulates
radiation levels with crisp guidelines. A possibilistic solution that produces a
feasible DVH based on a possibilistic attenuation matrix might still violate FDA
regulations, which are evaluated crisply. I put forth that optimization problems
whose solutions will be held to crisp standards should be solved with crisp
optimization.
An additional question that must be considered is whether or not the quality
of the solution to the fuzzy problem does not significantly surpass the quality of
its crisp counterpart. If it does, the marginal quality must be enough to outweigh
the computational inefficiency of optimization under uncertainty. In the case of
the RTP implementations in this paper, the crisp formulation produced a DVH
of quality comparable to the most successful fuzzy and possibilistic techniques
in much less time.
87
6. Avenues for Further Research
This thesis is by no means a conclusive work, having presented at least as
many issues as it resolved. Several areas of research that I would like to pursue
in the future are outlined below:
i. Fuzzy robust programming depends on setinclusive constraints rather
than inequality constraints. Setinclusive constraints with respect to crisp
sets (intervals) were considered by Soyster, [23], [24], [25]. He treated the
case of interactive parameters, and developed a dual of the linear program
with setinclusive constraints. He proposed solving the set inclusive LP
via the dual. I would like to study the possible extension of this solution
method to fuzzy linear programs with setinclusive constraints.
ii. Using the same a value for both the possibilistic and fuzzy components of
the problem is convenient, but does not necessarily provide a meaningful
model of reality. A more useful approach might be to require a higher
level of satisfaction for a solution with a smaller possibility of attaining
feasibility. This is akin to the idea of greater risktaking for the possibility
of greater rewards.
iii. Every fuzzy and possibilistic technique covered in this paper assumes non
interactivity of the parameters. As in the case of the RTP problem, in
teractive parameters arise in practice. How does using a model designed
for noninteractive parameters impact the solution of a problem whose
88
parameters are, in actuality, interactive? Noninteractivity in fuzzy and
possibilistic optimization is closely related to independence in stochastic
optimization, so similar research regarding the assumption of independent
parameters in stochastic models might illuminate the issue.
iv. At the heart of any optimization under uncertainty is the vehicle by which
the uncertainty is defined be it the endpoints of the interval, the fuzzy
membership function, the possibility distribution, or the probability distri
bution. Membership functions and possibility distributions are subjective
in nature, which raises many questions. How should they be elicited and
defined? How does the definition affect the solution of the problem? Is it
theoretically and practically sound to adjust the membership function or
possibility distribution in order to obtain a more desirable solution? If so,
how should it be done?
89
