- Email: [email protected]
Low-dose-rate extrapolation of data from animal carcinogenicity experiments— analysis of a new statistical technique
Low-Dose-Rate from
Animal
Analysis
Extrapolation Carcinogenicity
of a New
Statistical
of Data ExperimentsTechnique
HARRY A. GUESS National Institute of Environmental Health Sciences, National Institutes of Health, Research Triangle Park, North Carolina AND
KENNY
S. CRUMP
Louisiana
Technical
Universip,
Ruston, Louisiana
ABSTRACT Use of data from animal experiments to estimate the human long-term exposure to very low doses of chemicals in the environment
cancer risk due to poses a number of
biological and statistical problems. One of the statistical problems is to extrapolate the animal dose-response relation from the high dose range where test data are “available to low doses, which humans might encounter. We develop the mathematical properties of a new technique for performing this extrapolation. (Strictly speaking, it is an interpolation, since some animals are always tested at the background dose as a control.) The technique is based on maximum-likelihood estimation of parameters in dose-response relations whose form is derived
from a general
multistage
of stages is itself an unknown, the number although only finitely many will be nonzero
carcinogenesis
model.
Since the number
of parameters to be estimated is infinite, for any given set of data. Existence and
uniqueness properties of the estimates are proved using results of Karlin and Studden and of Krein and Rehtman from the theory of interpolation of polynomials with nonnegative coefficients. The algorithm for calculating the estimates is based on reducing an infinite set of Kuhn-Tucker conditions to an equivalent finite set of conditions, by exploiting properties of lexicographic ordering. This paper deals with the mathematical properties of the estimation technique. Statistical properties will be discussed in a subsequent paper.
1.
BACKGROUND
The procedure currently used to assess the human cancer risk associated with various chemicals may be thought of as consisting of four parts, involving biological, statistical, pharmacological, and to some extent epidemiological considerations (Hoe1 et al. [l]). First, data are obtained from MATHEMATICAL
0 American
Elsevier
BIOSCIENCES Publishing
32, 15-36 (1976)
Company,
Inc., 1976
15
16
HARRY
A. GUESS AND KENNY
S. CRUMP
laboratory experiments in which groups of animals are exposed to the compounds under investigation. This step includes evaluating the biological appropriateness of the experimental procedure and assessing the quality of the data. Second, the data are analyzed statistically to obtain an animal dose-response relationship applicable to the low-dose range. Since it is typically impractical to test more than a few hundred animals at a given dose level, the doses at which most of the tests are conducted must be high enough to produce positive responses (e.g., the appearance of a specific type of tumor) in about 1 percent or more of the animals tested. It is necessary to use these data to estimate the maximum dose levels at which the predicted responses should be below some specified upper limit, typically in the range 10m5 to lo-‘. This problem, which is essentially one of interpolating animal dose-response relationships in the range between the lowest positive dose tested and the background dose, is commonly referred to in the cancer literature as the low-dose extrapolation problem. After a low-dose extrapolation has been obtained for the animal dose-response relation, interspecies conversion factors based on biological and pharmacological considerations are applied to obtain an equivalent human dose-response relationship for the low-dose range. Finally, these results are combined with other information (e.g., the existence of susceptible human population subgroups) to obtain an overall assessment of human risk. Further discussion of all of these steps is given by Hoe1 et al. [l] and in a news article in Science (Gillette [2]). The focus of this paper is on the mathematical properties of a new technique for performing the low-dose extrapolation that constitutes the second part of the process discussed above. The statistical properties of the technique, such as strong consistency and asymptotic normality, will be discussed in a subsequent paper. Although at present there is no universally accepted method for performing this extrapolation, a common technique which has been used in developing environmental standards is the Mantel-Bryan procedure (Mantel and Bryan [3]; Mantel et al. [4]). This technique uses a probit model of the form
(for the case where the response at background dose may be taken to be negligible) to set “virtually safe” dose levels based on evaluation of dichotomous response data. Here @ is the normal cumulative distribution function, d is the dose rate, pd is the response probability at dose rate d, /I is a parameter which is conservatively assumed to be 1 in most applications, and CYis a parameter to be conservatively estimated from the data. Abbott’s correction [5] is used to adjust for background.
EXTRAPOLATION
2.
OF ANIMAL
INTRODUCTION
CARCINOGENICITY
DATA
17
AND SUMMARY
The approach we will take to the low-dose extrapolation problem is to base estimates on curves whose general parametric form is derived from a class of models which reflect what a number of investigators have judged to be plausible mechanisms of carcinogenesis. Specifically, the technique involves calculating maximum-likelihood estimates of coefficients in a general multistage carcinogenesis model developed by Armitage and Doll [6], extended by Peto (unpublished work) and Crump et al. [7], and investigated by Brown [8]. In this model it is assumed that K, > 1 different events must occur in a single cell before cancer is initiated in the cell and that the times-to-occurrence of these events are independently distributed according to exponential distributions with rates LY~ + &dt per unit time, for i = 1.2,. . . , K,. Here d is the dose rate; LY~ is the rate at which an initiating event for the ith stage occurs due to spontaneous (i.e., not dose-related) carcinogenesis, and & is the rate per /ith power of dose at which an initiating event for the ith stage occurs due to dose-related carcinogenesis. From these assumptions it is shown by Crump et al. [7] that if the time from cancer initiation in a single cell until an observable cancer develops in a tissue is assumed to be stochastically independent of the time to cancer initiation in a single cell and functionally independent of the dose rate, then the probability P(d, t) of an observable response by time t in a tissue exposed to a constant dose rate d is given by P(d,t)=l-exp[-Q(d)F(1)1,
where F is some function
Q(d)=
of time and Q is a polynomial
01,z 0.
;i (a,+/$+ i-l
of the form
p, > 0.
Observe that the degree of Q may be less than K, if some of the b, are zero. The degree of Q will be called the number of dose-related stages associated with the carcinogenic agent under test and will be denoted by K Q K,. In dichotomous data, which we consider, the time is taken to be constant and on a par with the animal life span. Thus, in the subsequent analysis we will assume that the response probability function P(d) has the form P(d)=l-exp
i I=1 1 -
In the sequel we will not be concerned
2
(q+&d4)
with the parameters
(1)
(Y,and pi as
18
HARRY
such. Instead,
A. GUESS
AND
KENNY
S. CRUMP
we rewrite (1) in the form P,(d)ZP(d)=l-e-‘(d),
(2)
where the function f has the form
5 _I?‘,
f(d)=
f; > 0.
i=o
Although this functional form of the response probability function P(d) was derived from a model involving several assumptions about carcinogenesis, it is actually somewhat more general than would seem to be indicated by the detailed discussion above. In fact, this functional form simply amounts to assuming that the response probability function has the form (2) where f is a polynomial with nonnegative coefficients. We regard the coefficients f; > 0 as parameters to be estimated from the data. The degree off, K, which is also unknown, may be interpreted as the number of dose-related stages in the carcinogenic process associated with the chemical under test. We will use the method of maximum likelihood for the estimation. Observe that since K is unknown, there are actually infinitely many parameters to be estimated in this model. In the experimental setting which we consider, a population of laboratory animals is randomly subdivided into n + 1 treatment groups, each consisting of n, animals. The animals in the ith treatment group are exposed to a dose rate di of a chemical compound whose carcinogenicity is to be evaluated. We denote by X, the number of animals in the ith treatment group that are observed to have a positive response. Then the likelihood function of the experimental outcome is given by
L=fi ; j-0 t ’
1[
P(cq]'[
1- P(q)y,
assuming the animals respond independently. When P, f,and L are given by (2). (3) and (4) respectively, we will write L(f)= L. Taking the natural logarithm of L(f) gives lnL(f)=C+
i j=O
xjln[l-
.-‘@)I-
i:
(nj-xj)f(4),
(5)
J=o
where C is a constant that does not depend on f. As long as f(d,) > 0 for x, > 0, In L(f) is well defined. If x, =O, it is consistent with (1) to take x,In[l-e-‘(d,)]=O even if f(c$)=O. It is easy to see that InL(f) is well defined for any f that maximizes L(f), so maximizing L(f) is equivalent to maximizing In L(f).
EXTRAPOLATION
OF ANIMAL
CARCINOGENICITY
19
DATA
It is easily seen that lnL(f) is a concave function of f,,,_f,, . . ,fK for fixed values of nj, x1, and 4, so for each fixed value of K the problem of obtaining maximum-likelihood estimates of fO,fi , . . . ,fK is the standard nonlinear-programming problem of maximizing a concave function subject to the linear constraints A. > 0, 0 < i Q K. In our case, however, K is unknown, so the problem is an infinite-dimensional nonlinear-programming problem. In the next section we investigate the problem of maximizing the likelihood function globally over all K= 1,2,. . . . First we introduce some additional definitions and notation. 3.
MODEL
FORMULATION
AND CHARACTERIZATION
A function f is said to be absolutely monotonic (a.m.) on the interval [0, b), where 0 < b < cc, whenever f(x) = Z~=,,fkxk has radius of convergence greater than or equal to b and fk > 0 for all k > 0 (cf. Widder [9], pp. 144 ff.). Let AM denote the set of all a.m. functions on [0, cc) and let AMP denote the set of all polynomials that are a.m. on [0, cc). In this terminology our initial problem is to find a polynomial Q, E AMP such that L( Q,) = sup[ L( Q)] Q E AMP]. We will need the following lemmas. LEMMA
I
If O
then there exists
QEAMP
Proof: If d,,> 0 this is contained in Karlin and Studden [ 10, p. 230, Theorem 10.11, so we assume d,=O. The function g(d)=f (d)-f (0) is a.m., and g(O)=O. Without loss of generality we may assume n to be even, so, by the theorem referred to above, there exists P E AMP such that P (di) = g(dJ for 1
lemma is an immediate
consequence
of Lemma
1.
2
sup[L(Q)lQ EAMP]=sup[l(f)lf attained, so is the other.
EAM],
and if one of the suprema is
We next give an example which shows that this common supremum can fail to be attained, illustrates geometrically what goes wrong, and suggests what to do about it. Suppose the dose-response data have the form shown in Fig. 1, where
y,=-ln
i
l-5
1
n, .
20
HARRY
A. GUESS AND KENNY
S. CRUMP
FIG. I.
Since an a.m. function is either strictly convex or else linear, no a.m. function can pass through all four points (&v,). However, these points can be approximated arbitrarily closely by a.m. polynomials. Furthermore it is easy to see that sup[L(_/)lf EAM] is equal to the value of L(f) which would be calculated by setting f (c$)=yj; any other values of f (4) will yield a strictly lower value for L(f). Thus, in this example sup[L(f)]f~AM] cannot be attained. The geometric reason why the supremum is not attained in this case has to do with what is essentially an end effect at the highest dose point. If we take the maximum likelihood over all a.m. polynomials of degree K or less and then let K increase, we see that it becomes possible to fit the point (d4,y4) more and more closely by a term of increasingly high order which approaches zero for d-c d4 as K approaches infinity. Since there are no data at doses than d4, there is nothing to prevent the curve from approaching a discontinuous positive jump at d4 as polynomials of increasingly high degree are admitted to the class over which the maximum is taken. This phenomenon should not be regarded as being unique to this rather special example; the point is that if the response at highest dose lies above the general trend shown by responses at lower doses, then as a.m. polynomials of increasingly high degree are admitted, it becomes possible to fit the response at highest dose simply by having a curve which hooks up sharply at d,. In other words, if we admit a.m. polynomials of arbitrarily high degree, then we are implicitly permitting the maximum amount of upturn in the dose-response curve at highest dose to be unrestrained by the requirements imposed by responses at lower doses.
EXTRAPOLATION
OF ANIMAL
CARCINOGENICITY
21
DATA
Fortunately, it is possible to avoid these end effects by enlarging the class of functions over which the maximum is taken to include functions which are a.m. on the interval [O,d,) and are permitted to take a positive jump at d,.Let
sd.(d)=O,
= 1,
d< 4,
(6)
d > d,,
and let (7)
AMJ=[flf(d)=Q(d)+c&(d),QEAMP,c>O].
The dose range of greatest interest for purposes of assessing the risk of cancer in humans is typically far below the highest dose tested in animal experiments. Tests on laboratory animals are usually conducted with very large doses in order to be reasonably sure of getting a response in the small groups of animals tested. The typical problem in constructing a dose-response curve is to extrapolate the data from the high dose range, where most of the tests are conducted, to the very low dose range, which is of greatest interest for purposes of setting human exposure limits. In view of this it makes sense biologically, as well as mathematically, to permit the response at the highest dose to be accommodated, if necessary, by means of a jump at d,.This is effectively what happens anyway if we allow K to be treated as the unknown that it is. We will see that with this formulation the maximization problem typically has a global solution and that the solutions can be computed. Figure 2 shows three dose-response curves of the form P,(d)=1 - e-f@‘) obtained by maximum likelihood estimation using data on the long-term exposure of female CF-1 mice to DDT in the diet.’ The algorithm discussed in Sec. 4 was used to compute the estimates. The lowest curve, labeled K= 00,comes from the function f, (d)= 0.04478497423 + 0.002052495669d + 0.05407064 (d), which maximizes the likelihood function over all curve, labeled K=20, comes from the polynomial
f EAMJ.
f2a (d)=0.04478497424+0.002052495667d+ (4.0139x
The middle
10-40)d’3
+(1.7200x
10-35)d’4+(8.2571
x 10-38)d15+(3.4130x
+(1.3741x
10-42)d’7+(5.5031
x
+ (8.8076 x
10-50)d20,
10-40)d’6
10-45)d'8+(2.2018x 10-47)d’9
‘The data are for liver hepatomas in females; data from the parental and Fl generations are pooled. These data are taken from Tables 3 and 4 of Tomatis et al. [ 1 I].
22
HARRY
-
observed of
0.W
AND
KENNY
S. CRUMP
frequency
liver
hepatomos
in female 05
A. GUESS
mice
I
0.0-L
1
50
0
100
150
200
I 250
d Dose (ppm)
FIG. 2.
which
maximizes
the likelihood
20 or less. The top curve,
function
labeled
over all a.m. polynomials
K = 4, comes
f4 (d) = 0.04483 + O.O02038d+ which
maximizes
the likelihood
function
curves
are nearly
coincident
of degree
the polynomial
(1.39 x 10-9)d4,
over all a.m. polynomials
4 or less. Here d is the dose rate, ppm of DDT of the mice. The three
from
of degree
in the diet over the life span
for doses
below
about
100 ppm.
The curve labeled K = 20 hooks up sharply at the highest dose, the upturn being caused by very high-order terms whose value is negligible except in the vicinity of the highest dose. Note that the low-order coefficients, which are the ones of importance for the low-dose extrapolation, are the same as the lower-order coefficients of the K= cc curve, up to 9 significant figures. In the low dose range [i.e., where the increased risk over background, P,(d)- PJO), is about 10m3 or less], the estimates of increased risk over background based on the curves K= 20 and K= cc agree to 9 significant figures. This illustrates with actual data the point discussed theoretically above, namely, that maximizing the likelihood function over all f E AMJ is equivalent to considering the limit of the maximizing polynomials as higher-
EXTRAPOLATION
and
OF ANIMAL
higher-degree
maximum
polynomials
is taken.
We will need
approach
and
the following
of Karlin
LE;IMA 3 (Karlin
are admitted
The former
cally equivalent to the latter implement computationally. theorem
CARCINOGENICITY
is much
lemma,
23
DATA
to the class
is both easier
which
over
theoretically to treat
is a minor
which
the
and practi-
theoretically modification
and of a
and Studden. and Studden)
Let O< d, < . . < d,, and suppose that c ,, . . ,c, are nonnegatice constants. Then a necessary and sufficient condition for there to exist an absolute& monotonic function g on [0, CC) mch that g(0) =O, g(d,) = c, for 1 Q i < n - 1, and g(d,J < c, is that
2 aid:+’ [ ,=I for any set of constants
>O
forallk>O]
*
[;,a,c,>O]
(8)
a,, . . ,a,.
Proof This follows from the proof of Theorem 10.3 on pp. 23 I-232 of Karlin and Studden [lo] because the functions zi,(k) = d,k+’ for i = 1.2,. . , n and k = 0, 1,2,. . . form a Tchebycheff system on the nonnegative integers, for the same reasons that the functions u,(k)= d,k do. Equation (8) above expresses the duality between the moment space and polynomial space for the system (6,) in the same way that Eq. (10.10) on p. 232 of Karlin and Studden expresses the duality for the system {u,}. THEOREM
I
Zf x, < n,, then there exists f EAMJ =sup[L(h)Jh
EAM].
such that L(f)=max[L(g)lgEAMJ]
Zf x, = n,, then for each ‘EAMJ,
L(f)
g)j gE
AMJI, Proof
We will give proof
for the case d, = 0; the proof
for the case d,, > 0
and let { f(“)}E=, be a sequence of is similar but easier. Let x,O. Since L(f(“)) < Cexp[ - (n,, - x,)f(“)(d,J] for some constant C, we have limsup, f(“‘)($.) < limsupnf(m)(d,)< cc for each 0 < j < n. Thus. by replacing the original sequence with a subsequence if necessary, we may assume that the limits cj = lim m--rmf (“)(d,) exist for each 0 < j < n. Let g(“)(d) = f c”)(d) -f(“)(O); g(“)(d,) = c, - c0 for each 0 < j < n. then gem) is a.m., g(“)(O) = 0, and lim,_, Hence by Lemma 3, if a,,. . ..a, are such that IX;= ,a,dF+ ’ > 0 for all then x7= ,a,g(“)(di) > 0. Taking limits, we see that k=O,l,...,
2 a,dk+’ > 0 i=l
for all k > 0
1 i5 +
a,(c,--a)>0
,=I
I
HARRY
24
A. GUESS AND KENNY
S. CRUMP
So by Lemma 3 there exists a function g absolutely monotonic on [0, co) such that g(0) = 0, g(d,) = c, - ce for 1 < i < n - 1, and g(d,) $ c, - cg. Using Lemma 1, let Q EAMP be such that Q (~I~)=g(d~)+ ce for 0 < i < n, let c= c,- Q(d,) > 0, and letf(d)= Q(d)+c&,(d). ThenfEAMJ, andf(d,)= c, for 0 < j < n. By construction, L(f) =sup[L(h)lh EAM]. But sup[L(h)lh EAM]= sup[L( g)l g EAMJ] because the functions g,(d)= Q(d)+ ~(d/d,)~ are in AM and g,(d,)+_f(d,) for 0 $ j < n, as k-co. This proves the part of the theorem pertaining to x, < n,. When x,,= n, it is easy to see that sup[L( g)l g E AMJ] cannot be attained, because, for g(d) = Q(d) + ~6~ (d), L(g) is a strictly increasing function of c.
Our next theorem and corollary characterize tions to this maximization problem are unique. THEOREM
the extent to which solu-
2
Zff,gEAMJ
and L(f)=L(g)=sup[L(h)lhEAMJ],
for each
g(d,)=f(d,)
2
C
njg(d,)=
then
~ESP, (9)
"jf(d,>y
JESZ
jESZ
where
SP= [jlO< j< SZ=[jJO<
n and x,>O],
j
(10)
Proof As discussed earlier, maximizing L is equivalent to maximizing lnL, so we have lnL(f)=lnL( g)=sup[lnL(h)lhEAMJ]. Let G(y)=ln(le-~)andletO
lnL(ef+(l-/3)g)-sup[lnL(h)lhEAMJ] =lnL(Bf+(l-8)g)-[@lnL(f)+(l-B)lnL(g)] =
x
x,T,,
jESP
where
Since f and g maximize L, it is easy to see that f(4) > 0 and g(4) > 0 for each jESP. Since G is strictly convex on (0, co), we have Tj > 0 for each jESP, and if ,f(d,J#g(+o) for some j0 E SP, then TiO> 0. This would imply In L( of+ ( 1- (3) g) > sup[ln L( h)l h E AMJ], which is impossible. Hence f(d,) = g(4) for all j E SP, and since In L(f) = In L( g), Eq. (9) follows.
EXTRAPOLATION
OF ANIMAL
CARCINOGENICITY
DATA
25
Using Theorem 2 and some results of Krein and Rehtman [ 121, we may formulate conditions for the likelihood function to have a unique maximum in AMJ. To state the conditions we will need to introduce the notion, taken from approximation theory, of the index of a function in AMJ. The definition below is adapted from Krein and Rehtman [12, p. 1271. DEFINITION
Let ~EAMJ be given by f(d)= Q(d)+ c&(d), where Q(d)=q,+q,d’ + . . . + qKdK. The index of ,f is defined as follows. Let k, < k, < . . . < k,,, be the subscripts of coefficients q, such that q, > 0. Divide the numbers k, beginning with k,,,, into groups in the following manner. If k, - k,_, > 1, then the first group consists only of k,; otherwise the first group consists of the two numbers k, and k,,,_,. Remove the first group from the sequence 1.) The last group consists of k, and k, in case (a) and of k, alone in case (b). Let s denote the number of groups. We say that the index of the polynomial Q is 2s unless the last group consists of the number 0 alone, in which case the index of q is defined to be 2s - 1. If c = 0, we define the index of f to be equal to the index of Q. If c > 0, the index off is defined to be one plus the index of Q. Essentially, the index of a function in AMJ is a weighted count of the number of its positive coefficients. Our uniqueness condition states that if a maximum-likelihood estimate f turns out to have an index less than the number of positive doses 4 for which the response X, is also positive, thenf is unique, i.e., L(f)> L(g) for all gEAMJ such that g#f. It is worth remarking that this uniqueness condition usually seems to be satisfied in practice as long as there are about three or more positive doses with positive responses. The condition may be stated precisely as follows: COROLLARY
I
Let NP denote the number of subscripts j > 0 for which both 4 > 0 and xj > 0. Let f E AMJ maximize L otter AMJ, and let I(f) denote the index off as defined above. If Z(f) Q NP- 1. then f is the unique function in AMJ maximizing L ouer AMJ. If d,= 0 and either x0 > 0 or NP= n, then uniqueness holds wheneoer the index of the function q, + q,d’ + ’ * . + q,dk-’ + cC3,Jd) is less than or equal to NP - 1. Proof: We will follow the argument end we first introduce some notation.
in Krein and Rehtman [12]. To this Let j, < . .
26
HARRY
A. GUESS
such that both d, > 0 and x, > 0. Let jEAMJ AMJ], and let Z(j) Q NP - 1. We may write
where Q,E AMP and c, > 0. For dNp_ , = d, let Let
z= {O,l)...)
AND
KENNY
satisfy L(f)=
qf,
=
cJ;
S. CRUMP
sup[l(h)]h
otherwise let
qf,
E
=
0.
co},
s={lE.qqf>o}, and
where 6;,, = 1 if i = j and S,,, = 0 otherwise. In this notation we may write
c,=
c qjwo
O
Its
where c;=f($,) for 0 < i Q NP- 1. We now proceed to show that the argument in Krein and Rehtman [12] applies to our case. It is easy to see that our definition of the index of the function f coincides with the definition p. 127 in Krein and Rehtman of the index of the sequence {q, : 1 E S }. This in turn is the same as their definition on p. 128 of the index of the representation {c,}~~;‘. Thus what we call Z(f) is what they call the index of the representation {ci}yL;‘. Hence by Theorems 4.3 and 3.1 and Eq. 3.4 of Krein and Rehtman, since I(f) < NP - 1, there is a unique polynomial Q EAMP satisfying
Q(d,)=cc
O
Q(d,p,p_,)~CNP-l.
Since Q, also satisfies these conditions, QJ- Q. Now let g be any other function in AMJ such that L(g)= L(f)=sup[L(h)]h EAMJ]. We may write
EXTRAPOLATION
OF ANIMAL CARCINOGENICITY
DATA
27
where Q, EAMP and cp > 0. It follows from Theorem 2 that Q,(d,) = Q,(<,) =c, for 0~ i 0, then by Theorem 2 g(d,,) =f(d,) and, since Q, = Q,, we have cg = c,. If X, =O, then it is not hard to show directly that cg = c,= 0. Hence in all cases cn = c, and Q, z Q,, so g-f. This completes the proof of the first assertion. Now let d,= 0 and either x,>O or NP= n. Then by Theorem 2,f(O) and hence 40’are uniquely determined, and so are
forO
Let f&d)=q(+q:d’+‘-
+q’dK-’ K
’
To prove the second assertion it suffices to use the above argument with 5, playing the role of c, and using Q, in place of Q,. This completes the proof. When the maximum-likelihood estimatef fails to be unique, there will be a (possibly infinite) set of low-dose extrapolations each of which maximizes the likelihood function. In such cases we seek a canonical way to choose one of them. One way which is biologically meaningful and also easy to handle mathematically is to seek the maximum probability of additional response over background, P,(d,) - P,(O). at a given environmental dose rate d, that can be calculated using anyf EAMJ such that f.(f)=sup[L( g)l gE AMJ]. Whenever da=0 and a positive response is observed at this dose, it follows from Theorem 2 that P,(O) is unique. In this frequently occurring case, maximizing P,(d,)- P/(O) is equivalent to maximizing P/(d,). Since Pf(d,) is an increasing function of f(~$,), this is equivalent to solving the problem maximize
f( d, )
(11)
subject to ~EAMJ
and
L(f)=sup]l(g)lgEAMJ],
where d, is a given dose rate. Once we have obtained one function sup[l(h)lh EAMJ], it follows from Theorem (11) and (12) is equivalent to the problem maximize
(12)
g EAMJ such that L(g)= 2 that the problem defined by
5 l;d,’ r=O
(13)
HARRY
28
A. GUESS
AND
KENNY
S. CRUMP
subject to
;$oi.d:+cs, (4)=s(d,)t
jESP,
I51; jSSZ c n,d,i = c n&d,), 1jESZ
i=o
(14)
(W2
[
J;>o,c>o,
(16)
and KE[O,1,2,
. ..I.
(17)
where d, and g(d,) (0~ j < n) are given, d,< d,,, and L(g)=max[l(h)]h~ AMJ]. For each fixed value of K this is a linear-programming problem and therefore is readily solvable by standard algorithms. We will show in Theorem 3 that a global maximum (over all K) exists. It is possible to modify a standard linear-programming algorithm to yield an algorithm which would converge to the global maximum. For most practical purposes, however, such a modification appears unnecessary. In typical applications d,
EAMJI}.
(18)
if and only if X, < n,. In this case we have
3
If x, < n,
and d, < d,,, then there exists f E M such that f(d,)
= sup[ g(d,)l
gEMI. ProoJ: By Theorem
1, M is nonempty,
M = { flf
E AMJ andf
‘If n E SZ, then c = 0; we omit the proof.
so let g E M; then by Theorem satisfies (9)).
2
(19)
EXTRAPOLATION
OF ANIMAL CARCINOGENICITY
Since de < d,, we have f(d,) < f(d,J,
DATA
29
and it follows from (9) that
)ISE M 1
sup[f(d,
It follows from (20) that there exists a sequence f’“‘(d)= and such that the following
cw {f(m)}z3,
in M such that
Qcm)(d)+ccm)& n (d)
(21)
limits exist:
lim f’“‘(d,)=sup[h(d,)Ih m--t’x
EM l=c, (22) O
,Il= f(“‘( 4. ) = c,,
In the remainder of the proof we will assume that do=0 and d,# d, for each j. When the contrary is true the proof is similar but easier. Let g(“)(d)
= Q(“)(d)
- Qcm)(0); then gCrn)~ AM, gCm’(0) = 0,
and
g’“‘(d n ) < f’“‘(d n )-f’“‘(O). Hence by Lemma 3 and Eq. (22) i ajO, j=l
all
k>O I
j=l n
2 a,(c,-c,)+a,(c,-co)>0
=+
[
Using Lemma 3 in the other direction, such that
we conclude
h (de>= ce- co,
h(4) = cj h(d,)
< c,-
co,
co.
O
1. that there exists h EAM
30
HARRY
By Lemma 1 there exists Q EAMP h(d,)+ ca for 0 < j < n. Let
A. GUESS
AND
KENNY
S. CRUMP
such that Q (d,)= h(d,)+ c,, and Q(4)=
c=c,-Q(d,,)=c,-c,-h(d,,)>O, and let (23) Then f E AMJ,
and for O
&-py~).
f(d,)=c,=
Since each f(“) satisfies (9) so does f. Hence f~ M, and f(d,)= hE M]. This completes the proof. 4.
GLOBAL
MAXIMIZATION
sup[h(d,)]
ALGORITHM
In this section we develop an algorithm for computing a function ~EAMJ such that L(f)=sup[l(h)]h EAMJ]. We outline the general structure of the algorithm and prove convergence. Each f E AMJ has the form
f(4=Q(d)+4in(d)a
(24)
where Q is a polynomial of arbitrary degree with nonnegative coefficients, c is a nonnegative number, and Sdn is the jump at d,, defined by (6). We may write Q in the form Q(d)=
q, > 0,
5 q;d’, r=O
(25)
where all but finitely many of the coefficients q, are zero. Since maximizing L(fl is equivalent to maximizing lnL(f), and since lnL(f) is a concave function of the vector (c,q)=(c,q,,q,, .), we will work with lnL(f) rather than with L(fl. Our problem may be stated as follows maximize
(26)
lnL(f)
subject to the constraints c>o
q, >O
for O< i< K,
and
q,=O
for i>K
(27)
EXTRAPOLATION
OF ANIMAL
CARCINOGENICITY
31
DATA
and K=O,l,2
,...,
(28)
where lnL(f) is given by (5) f 1sgiven by (24) and Q is given by (25). By Theorem 1 this problem has a solution if and only if X, < n,. For each fixed value of K, (26) and (27) define a concave nonlinear-programming problem with linear constraints, and hence the Kuhn-Tucker conditions are necessary and sufficient for a vector (c, qO,. . . , qK) to yield a maximum. For future reference we state this formally as a lemma. LEMMA
4 (Kuhn-Tucker)
For each fixed K the following conditions are necessary and sufficient for a function f to soke the nonlinear-programming problem defined by (26) and (27), where lnL(f) is g&en by (5), f by (24), Q by (25), and P by (2): aInL(f)
-i
84,
[-n,+&]$
=0
if
q,>O, (29)
J
J=o
GO
if
q,=O
for all 0 < i < K, and
alnL(f) -[
X,
ac
-
-%+
P(d,)
=o
if
c >O,
GO
if
c=O.
I
(30)
ProojI Equations (29) and (30) are the Kuhn-Tucker conditions for the nonlinear-programming problem defined by (26) and (27). Since the constraints are linear and the function to be maximized is concave, the Kuhn-Tucker conditions are necessary and sufficient. (See Zangwill [13], Theorems 2.14 and 2.15 and Exercise 2.13.) We will also have occasion to refer to the nonlinear Kth-dimensional problem where the jump is not allowed at the highest dose. This problem may be stated maximize
lnL(f)
(31)
0 d i Q K,
(32)
subject to the constraints Aa0
for
where lnL(f) is given by (5) and f is given by (3). This is also a concave nonlinear-programming problem, and the corresponding necessary and
32
HARRY
sufficient
Kuhn-Tucker alWf> 56
conditions = i
A. GUESS
AND
S. CRUMP
for a function f to be a solution
[ -n,+&]$=O
are
if J;>O, J
j=O
GO
for all 0~ i < K, where lnL(f) LEMMA
KENNY
is given by (5) andf
if
J=O
is given by (3).
5
In order that f EAMJ sati.& In L(fl=sup[lnL(h)/h sary and sufficient that (29) and (30) hold for all i > 0.
EAMJ]
it is neces-
Proof: Necessity is an obvious consequence of Lemma 4. If the conditions hold for all i > 0, then from (30) x, < n,, and so, by Theorem 1, there exists g EAMJ such that lnL(g)=sup[lnL(h)Jh EAMJ]. We may write g(d)= 0 (d)+ E&(d) where Q is an a.m. polynomial. It follows from Lemma 4 that lnL(f) > lnL(g). Hence lnL(f)=sup[lnL(h)]h EAMJ]. This proves sufficiency.
The obvious analog of this lemma holds for the case where no jump at d, is permitted. In this case a global maximum need not exist, but (29) still is necessary and sufficient. The exact result is as follows: COROLLARY In order that QEAMP satisb lnL(Q)=sup[lnL(Q)lQEAMP] necessary and sufficient that (29) hold for all i > 0. (Here
c =0
it is in all
formulas.). Proof Necessity is obvious. To prove sufficiency (29) (and that c = 0). Then
lnL(Q)=sup[lnL(Q)IQEAMP,(degreeof
suppose that Q satisfies
Q)
(33)
holds for each K, by the Kuhn-Tucker theorem. Since (29) holds for all - n,, + x,/ P(d,,) < 0 and hence x, < n,. By Theorem 1 there exists gEAMJ such that lnL(g)=sup[lnL(h)]h EAMJ]. Let g(d)= &(d)+ c&,(d), where Q EAMP. For each k let Qk be defined by Qk(d)= o(d)+ c(d/dJk. Then Q,EAMP and limk_,Qk(<)=g($ for O< j < n, and hence lim ,,,lnL(Q,)=lnL( g). It follows from (33) that i > 0, we have
lnL(Q,)g for each k, and so lnL(Q)=
lnL(Q)
< lnL(g)
lnL( g). This proves sufficiency.
We now have a set of infinitely many conditions which together are necessary and sufficient for a function to globally maximize the likelihood
EXTRAPOLATION
OF ANIMAL CARCINOGENICITY
DATA
33
function. The key step in the development of our algorithm is to replace this set by a set which has only finitely many conditions and yet is necessary and sufficient. To do this we will introduce the following lexicographic order relation. DEFINITION
A vector (u,, . . . , u,) will be said to be iexico-negative if the first nonzero term in the sequence u,,, u,_ ,, . . . , u. (read starting with u,) is negative. [Thus the vectors (O,O,. . . ,O) and (2,2,. . , ,2, - 1) are lexico-negative.] LEMMA
6
(a) The vector (u,, . . . un) defined by uj =
I
-T+
xj P(4)
(34) 1
is lexico-negative if and only if there exists some integer i, such that ClIn L(f)/ aq,< 0 for all i > i,. (b) If (u,, . . , u,,) defined by (34) is Mexico-negative, then 8 In L( f)/ aq, < 0 for all i > i,, where i, is given by i,=O
if all
U, GO,
OG i < n,
(35)
and otherwise
In i,=l+
j@axo~j~jo-I[-9+X,/p(4)l I-ni0+x/0/p(40))) I
40 In &-I i-1
(36)
where j, is the index of the first nonzero term in the sequence u,, u,_ I,...,t+l (starting with u,). Proof Both (a) and (b) follow from straightforward will be omitted. THEOREM
computations,
which
4
The following conditions are necessary and sufficient for In L(f) = sup[lnL(g)(gEAMJ], where InL(f) is given by (5), f by (24), Q by (25), and P by (2): q,=O
foraIli>
K,
for some K,
(uo,..., u,) is lexico-negative, Eqs. (29) and (30) hold for all i < K,
(37) (38) (39)
HARRY
34
A. GUESS
AND
KENNY
S. CRUMP
and
alnuf)
ho
forall
K+l
a4i
Here the uj are defined by (34). and i, is defined i, < K + 1, then (40) is CacuousIy satisfied.)
by (35) and
(36).
(If
Proof. We first prove necessity. Let ~EAMJ, and let lnL(f)= sup[lnl(h)lh EAMJ]. Then (39) holds by Lemma 4 and, since f EAMJ, we have (37). It follows from Lemmas 5 and 6 that (38) holds. By (37) and Lemma 5 we get (40). This proves necessity. We now prove sufficiency. By (38), (40) and Lemma 6 we have a lnL( f)/ aq, G 0 for all i > K + 1. This, together with (37) and (39) implies that the conditions of Lemma 5 are satisfied. Hence we have sufficiency. COROLLARY
Equations (37) through (40) omitting the requirement of (30) in Eq. (39) and taking c = 0 throughout, are necessary and sufficient conditions for there to exist
Q EAMP
such that lnL(Q)=sup[lnL((i)lQ
The proof is similar to that of Theorem
Proof
Theorem 4 forms the basis for our algorithm. algorithm is as follows:
EAMP]. 4 and is omitted. The general scheme of the
(a) Starting with a given K, solve the nonlinear-programming problem defined by (26) and (27). (b) Once a maximum for fixed K has been obtained in step (a), Eqs. (38) and (40) are checked. If both are satisfied, then (by Theorem 4) the current value of (c, qO,. . . , qk) is globally optimal and the algorithm stops. If either (38) or (40) fails to hold, we replace K by K+ 1 and return to step (a), using our current vector (c, q,,, . . , qk, 0) as a starting point for performing the new maximization. Since a global maximim exists when X, < n, (by Theorem 1) and since the conditions checked by the algorithm each time K is increased are necessary and sufficient for a vector (c,q,, . . ,qk) to yield a global maximum, the algorithm will eventually converge to a global maximum. This proves convergence of the algorithm to a global maximum. We now state briefly the method we have used to solve the finite-dimensional nonlinear-programming problem defined by (31) and (32), i.e., step (a) above. We will need the following lemma. LEMMA
7
Suppose for a fixed K that a solution f’ to the nonlinear-programming problem defined by (31) and (32) satisfies the condition x,/n, > P(f’(d,)). If
EXTRAPOLATION
OF ANIMAL
CARCINOGENICITY
f* is any solution to the same nonlinear-programming
replaced by n - 1, then x,/n,
35
DATA
problem except that n is
> P(,F(d,,)).
Proof: Let
and
We then have e,(f’)98,_,(_~‘)<8,_,(f*) and e,,(f’)=e,_,(f’)+a,(f’). Now suppose the conclusion of the theorem does not hold, that is, x,/n, < P(f*(d,)), and define g,=cuf*+(l -LY)~), 0~ a g 1. It is easy to see that g,EAMP (O h(a), since Q,-,(f*) is a maximum. Thus h(a) is increasing on 0 < cx < 1, and consequently Q,_,(f’)=h(O)< h(a,)=8,_,(ga,)< h(l)=Q,_,(f*). On the other hand, we have the strict inequality a,( g,,)> afi(f’) because x,/n, is the unique value of P(f(d,)) that yields the global maximum for a,. Therefore, B,~,(s,,)+a,(g,l)>8,-,(f’)+a,if’), tradicts our assumption that Q,(j)
i.e.,
Qn(g,,)>Q,,(f’). and This completes
IS a maximum.
Now suppose we have solved the nonlinear-programming fined by (3 1) and (32) and obtained the solution f’ EAMP,
this conthe proof.
problem dewhich must
the Kuhn-Tucker conditions (29). If, in addition, x,/n, < seen that f(d) =J’(d) + c6dn(d) EAMJ, where c = 0, must satisfy the Kuhn-Tucker conditions (29) and (30) for the problem defined by (26) and (27). If, instead, x,/n, > P(f’(d,,)); we proceed to obtain f* EAMP that solves the problem defined by (31) and (32) with n necessarily
P(f’(d,
replaced
satisfy
)), it is easily
by
n - 1.
We
are
x,/n, > P(fc(d,)). It follows c&(d) EAMJ, where c= x,/n,
guaranteed easily
from
by this
the
above
inequality
- P(fc(d,, )) satisfies
lemma
that
that f(d)=f*+
the conditions
(29) and
(30) and thus is a solution to the problem defined by (26) and (27). With this procedure the two problems defined by (26) and (27) and by (31) and (32) respectively can both be solved using a single algorithm. The algorithm we have used for this purpose is basically the method of Davidon, Fletcher, and Powell; see Fletcher and Powell [14]. It was necessary to modify this method, which applies only to unconstrained optimization problems, to take into account our nonnegativity constraints.
36 5.
HARRY
FUTURE
A. GUESS
AND
KENNY
S. CRUMP
WORK
In a subsequent paper we will investigate the statistical aspects of this technique. We will show that in general the maximum-likelihood estimates of the response probabilities P(d) are strongly consistent and unbiased for all d< d,,. Under some additional hypotheses the estimates will be shown to be asymptotically normal. Using asymptotic normality, confidence intervals will be computed and compared with simulation results. Finally, risk estimates obtained using this technique will be compared with risk estimates using several widely known techniques. We thank Kevin Deal for his efficient programming tions regarding computational methods.
and for useful sugges-
REFERENCES D. G. Hoel, D. W. Gaylor, R. L. Kirschstein, U. Saffiotti, and M. A. Schneiderman, Estimation of risks of irreversible, delayed toxicity, J. Toxicol. and Enoiron. Health 1, 133-151(1975). R. Gillette, Cancer and the environemnt (II): Groping for new remedies, Science 186, 242-245 (1974). N. Mantel and W. R. Byran, Inst. 27, 455470 (1961).
Safety
testing
of carcinogenic
agents,
J. Nail. Cancer
N. Mantel, N. Bohidar, C. Brown, J. Ciminera, and J. Tukey, An improved “MantelBryan” procedure for “safety” testing of carcinogens, Cancer Res. 35, 865-872 (1975). W. S. Abbott, A method of computing the effectiveness of an insecticide, J. Econ. Entomol. 18, 265-267 (1925). P. Armitage and R. Doll, Stochastic models for carcinogenesis, Proc. Fourfh Berkeley Symp. Math. Stat. Probab. Univ. of Calif. Press, Berkeley, 1961, Vol. 4, pp. 19-38. K. S. Crump, D. G. Hoel, C. H. Langley, and R. Peto, Fundamental carcinogenic processes and their implications for low-dose risk assessment, Cancer Res. (1976), to 8 9 10 11 12
13 14
appear. C. C. Brown, Statistical aspects of extrapolation of dichotomous dose-response data, submitted for publication. D. V. Widder, The Laplace Transform, Princeton U. P., Princeton, 1941. S. Karlin and W. J. Studden, Tchebycheff Systems: With Applications in Analysis and Statistics, Wiley, New York, 1966. L. Tomatis, V. Turusov, N. Day, and R. T. Charles, The effect of long-term exposure to DDT on CF-1 mice, Inr. J. Cancer 10, 489-506 (1972). M. G. Krein and P. G. Rehtman, Development in a new direction of the CebysevMarkov theory of limiting values of integrals, Usp. Mat. Nauk (N.S.) 10, No. 1 (63), 67-78 (1955); Am. Math. Sot. Transl. Ser. 2, 12, 123-135 (1959). W. I. Zangwill, Nonlinear Programming, A Unified Approach, Prentice-Hall, Englewood Cliffs, N. J., 1969. R. Fletcher and M. J. D. Powell, A rapidly convergent descent method for minimization, Comput. J. 6, 163-168
(1963).
Animal
Analysis
Extrapolation Carcinogenicity
of a New
Statistical
of Data ExperimentsTechnique
HARRY A. GUESS National Institute of Environmental Health Sciences, National Institutes of Health, Research Triangle Park, North Carolina AND
KENNY
S. CRUMP
Louisiana
Technical
Universip,
Ruston, Louisiana
ABSTRACT Use of data from animal experiments to estimate the human long-term exposure to very low doses of chemicals in the environment
cancer risk due to poses a number of
biological and statistical problems. One of the statistical problems is to extrapolate the animal dose-response relation from the high dose range where test data are “available to low doses, which humans might encounter. We develop the mathematical properties of a new technique for performing this extrapolation. (Strictly speaking, it is an interpolation, since some animals are always tested at the background dose as a control.) The technique is based on maximum-likelihood estimation of parameters in dose-response relations whose form is derived
from a general
multistage
of stages is itself an unknown, the number although only finitely many will be nonzero
carcinogenesis
model.
Since the number
of parameters to be estimated is infinite, for any given set of data. Existence and
uniqueness properties of the estimates are proved using results of Karlin and Studden and of Krein and Rehtman from the theory of interpolation of polynomials with nonnegative coefficients. The algorithm for calculating the estimates is based on reducing an infinite set of Kuhn-Tucker conditions to an equivalent finite set of conditions, by exploiting properties of lexicographic ordering. This paper deals with the mathematical properties of the estimation technique. Statistical properties will be discussed in a subsequent paper.
1.
BACKGROUND
The procedure currently used to assess the human cancer risk associated with various chemicals may be thought of as consisting of four parts, involving biological, statistical, pharmacological, and to some extent epidemiological considerations (Hoe1 et al. [l]). First, data are obtained from MATHEMATICAL
0 American
Elsevier
BIOSCIENCES Publishing
32, 15-36 (1976)
Company,
Inc., 1976
15
16
HARRY
A. GUESS AND KENNY
S. CRUMP
laboratory experiments in which groups of animals are exposed to the compounds under investigation. This step includes evaluating the biological appropriateness of the experimental procedure and assessing the quality of the data. Second, the data are analyzed statistically to obtain an animal dose-response relationship applicable to the low-dose range. Since it is typically impractical to test more than a few hundred animals at a given dose level, the doses at which most of the tests are conducted must be high enough to produce positive responses (e.g., the appearance of a specific type of tumor) in about 1 percent or more of the animals tested. It is necessary to use these data to estimate the maximum dose levels at which the predicted responses should be below some specified upper limit, typically in the range 10m5 to lo-‘. This problem, which is essentially one of interpolating animal dose-response relationships in the range between the lowest positive dose tested and the background dose, is commonly referred to in the cancer literature as the low-dose extrapolation problem. After a low-dose extrapolation has been obtained for the animal dose-response relation, interspecies conversion factors based on biological and pharmacological considerations are applied to obtain an equivalent human dose-response relationship for the low-dose range. Finally, these results are combined with other information (e.g., the existence of susceptible human population subgroups) to obtain an overall assessment of human risk. Further discussion of all of these steps is given by Hoe1 et al. [l] and in a news article in Science (Gillette [2]). The focus of this paper is on the mathematical properties of a new technique for performing the low-dose extrapolation that constitutes the second part of the process discussed above. The statistical properties of the technique, such as strong consistency and asymptotic normality, will be discussed in a subsequent paper. Although at present there is no universally accepted method for performing this extrapolation, a common technique which has been used in developing environmental standards is the Mantel-Bryan procedure (Mantel and Bryan [3]; Mantel et al. [4]). This technique uses a probit model of the form
(for the case where the response at background dose may be taken to be negligible) to set “virtually safe” dose levels based on evaluation of dichotomous response data. Here @ is the normal cumulative distribution function, d is the dose rate, pd is the response probability at dose rate d, /I is a parameter which is conservatively assumed to be 1 in most applications, and CYis a parameter to be conservatively estimated from the data. Abbott’s correction [5] is used to adjust for background.
EXTRAPOLATION
2.
OF ANIMAL
INTRODUCTION
CARCINOGENICITY
DATA
17
AND SUMMARY
The approach we will take to the low-dose extrapolation problem is to base estimates on curves whose general parametric form is derived from a class of models which reflect what a number of investigators have judged to be plausible mechanisms of carcinogenesis. Specifically, the technique involves calculating maximum-likelihood estimates of coefficients in a general multistage carcinogenesis model developed by Armitage and Doll [6], extended by Peto (unpublished work) and Crump et al. [7], and investigated by Brown [8]. In this model it is assumed that K, > 1 different events must occur in a single cell before cancer is initiated in the cell and that the times-to-occurrence of these events are independently distributed according to exponential distributions with rates LY~ + &dt per unit time, for i = 1.2,. . . , K,. Here d is the dose rate; LY~ is the rate at which an initiating event for the ith stage occurs due to spontaneous (i.e., not dose-related) carcinogenesis, and & is the rate per /ith power of dose at which an initiating event for the ith stage occurs due to dose-related carcinogenesis. From these assumptions it is shown by Crump et al. [7] that if the time from cancer initiation in a single cell until an observable cancer develops in a tissue is assumed to be stochastically independent of the time to cancer initiation in a single cell and functionally independent of the dose rate, then the probability P(d, t) of an observable response by time t in a tissue exposed to a constant dose rate d is given by P(d,t)=l-exp[-Q(d)F(1)1,
where F is some function
Q(d)=
of time and Q is a polynomial
01,z 0.
;i (a,+/$+ i-l
of the form
p, > 0.
Observe that the degree of Q may be less than K, if some of the b, are zero. The degree of Q will be called the number of dose-related stages associated with the carcinogenic agent under test and will be denoted by K Q K,. In dichotomous data, which we consider, the time is taken to be constant and on a par with the animal life span. Thus, in the subsequent analysis we will assume that the response probability function P(d) has the form P(d)=l-exp
i I=1 1 -
In the sequel we will not be concerned
2
(q+&d4)
with the parameters
(1)
(Y,and pi as
18
HARRY
such. Instead,
A. GUESS
AND
KENNY
S. CRUMP
we rewrite (1) in the form P,(d)ZP(d)=l-e-‘(d),
(2)
where the function f has the form
5 _I?‘,
f(d)=
f; > 0.
i=o
Although this functional form of the response probability function P(d) was derived from a model involving several assumptions about carcinogenesis, it is actually somewhat more general than would seem to be indicated by the detailed discussion above. In fact, this functional form simply amounts to assuming that the response probability function has the form (2) where f is a polynomial with nonnegative coefficients. We regard the coefficients f; > 0 as parameters to be estimated from the data. The degree off, K, which is also unknown, may be interpreted as the number of dose-related stages in the carcinogenic process associated with the chemical under test. We will use the method of maximum likelihood for the estimation. Observe that since K is unknown, there are actually infinitely many parameters to be estimated in this model. In the experimental setting which we consider, a population of laboratory animals is randomly subdivided into n + 1 treatment groups, each consisting of n, animals. The animals in the ith treatment group are exposed to a dose rate di of a chemical compound whose carcinogenicity is to be evaluated. We denote by X, the number of animals in the ith treatment group that are observed to have a positive response. Then the likelihood function of the experimental outcome is given by
L=fi ; j-0 t ’
1[
P(cq]'[
1- P(q)y,
assuming the animals respond independently. When P, f,and L are given by (2). (3) and (4) respectively, we will write L(f)= L. Taking the natural logarithm of L(f) gives lnL(f)=C+
i j=O
xjln[l-
.-‘@)I-
i:
(nj-xj)f(4),
(5)
J=o
where C is a constant that does not depend on f. As long as f(d,) > 0 for x, > 0, In L(f) is well defined. If x, =O, it is consistent with (1) to take x,In[l-e-‘(d,)]=O even if f(c$)=O. It is easy to see that InL(f) is well defined for any f that maximizes L(f), so maximizing L(f) is equivalent to maximizing In L(f).
EXTRAPOLATION
OF ANIMAL
CARCINOGENICITY
19
DATA
It is easily seen that lnL(f) is a concave function of f,,,_f,, . . ,fK for fixed values of nj, x1, and 4, so for each fixed value of K the problem of obtaining maximum-likelihood estimates of fO,fi , . . . ,fK is the standard nonlinear-programming problem of maximizing a concave function subject to the linear constraints A. > 0, 0 < i Q K. In our case, however, K is unknown, so the problem is an infinite-dimensional nonlinear-programming problem. In the next section we investigate the problem of maximizing the likelihood function globally over all K= 1,2,. . . . First we introduce some additional definitions and notation. 3.
MODEL
FORMULATION
AND CHARACTERIZATION
A function f is said to be absolutely monotonic (a.m.) on the interval [0, b), where 0 < b < cc, whenever f(x) = Z~=,,fkxk has radius of convergence greater than or equal to b and fk > 0 for all k > 0 (cf. Widder [9], pp. 144 ff.). Let AM denote the set of all a.m. functions on [0, cc) and let AMP denote the set of all polynomials that are a.m. on [0, cc). In this terminology our initial problem is to find a polynomial Q, E AMP such that L( Q,) = sup[ L( Q)] Q E AMP]. We will need the following lemmas. LEMMA
I
If O
then there exists
QEAMP
Proof: If d,,> 0 this is contained in Karlin and Studden [ 10, p. 230, Theorem 10.11, so we assume d,=O. The function g(d)=f (d)-f (0) is a.m., and g(O)=O. Without loss of generality we may assume n to be even, so, by the theorem referred to above, there exists P E AMP such that P (di) = g(dJ for 1
lemma is an immediate
consequence
of Lemma
1.
2
sup[L(Q)lQ EAMP]=sup[l(f)lf attained, so is the other.
EAM],
and if one of the suprema is
We next give an example which shows that this common supremum can fail to be attained, illustrates geometrically what goes wrong, and suggests what to do about it. Suppose the dose-response data have the form shown in Fig. 1, where
y,=-ln
i
l-5
1
n, .
20
HARRY
A. GUESS AND KENNY
S. CRUMP
FIG. I.
Since an a.m. function is either strictly convex or else linear, no a.m. function can pass through all four points (&v,). However, these points can be approximated arbitrarily closely by a.m. polynomials. Furthermore it is easy to see that sup[L(_/)lf EAM] is equal to the value of L(f) which would be calculated by setting f (c$)=yj; any other values of f (4) will yield a strictly lower value for L(f). Thus, in this example sup[L(f)]f~AM] cannot be attained. The geometric reason why the supremum is not attained in this case has to do with what is essentially an end effect at the highest dose point. If we take the maximum likelihood over all a.m. polynomials of degree K or less and then let K increase, we see that it becomes possible to fit the point (d4,y4) more and more closely by a term of increasingly high order which approaches zero for d-c d4 as K approaches infinity. Since there are no data at doses than d4, there is nothing to prevent the curve from approaching a discontinuous positive jump at d4 as polynomials of increasingly high degree are admitted to the class over which the maximum is taken. This phenomenon should not be regarded as being unique to this rather special example; the point is that if the response at highest dose lies above the general trend shown by responses at lower doses, then as a.m. polynomials of increasingly high degree are admitted, it becomes possible to fit the response at highest dose simply by having a curve which hooks up sharply at d,. In other words, if we admit a.m. polynomials of arbitrarily high degree, then we are implicitly permitting the maximum amount of upturn in the dose-response curve at highest dose to be unrestrained by the requirements imposed by responses at lower doses.
EXTRAPOLATION
OF ANIMAL
CARCINOGENICITY
21
DATA
Fortunately, it is possible to avoid these end effects by enlarging the class of functions over which the maximum is taken to include functions which are a.m. on the interval [O,d,) and are permitted to take a positive jump at d,.Let
sd.(d)=O,
= 1,
d< 4,
(6)
d > d,,
and let (7)
AMJ=[flf(d)=Q(d)+c&(d),QEAMP,c>O].
The dose range of greatest interest for purposes of assessing the risk of cancer in humans is typically far below the highest dose tested in animal experiments. Tests on laboratory animals are usually conducted with very large doses in order to be reasonably sure of getting a response in the small groups of animals tested. The typical problem in constructing a dose-response curve is to extrapolate the data from the high dose range, where most of the tests are conducted, to the very low dose range, which is of greatest interest for purposes of setting human exposure limits. In view of this it makes sense biologically, as well as mathematically, to permit the response at the highest dose to be accommodated, if necessary, by means of a jump at d,.This is effectively what happens anyway if we allow K to be treated as the unknown that it is. We will see that with this formulation the maximization problem typically has a global solution and that the solutions can be computed. Figure 2 shows three dose-response curves of the form P,(d)=1 - e-f@‘) obtained by maximum likelihood estimation using data on the long-term exposure of female CF-1 mice to DDT in the diet.’ The algorithm discussed in Sec. 4 was used to compute the estimates. The lowest curve, labeled K= 00,comes from the function f, (d)= 0.04478497423 + 0.002052495669d + 0.05407064 (d), which maximizes the likelihood function over all curve, labeled K=20, comes from the polynomial
f EAMJ.
f2a (d)=0.04478497424+0.002052495667d+ (4.0139x
The middle
10-40)d’3
+(1.7200x
10-35)d’4+(8.2571
x 10-38)d15+(3.4130x
+(1.3741x
10-42)d’7+(5.5031
x
+ (8.8076 x
10-50)d20,
10-40)d’6
10-45)d'8+(2.2018x 10-47)d’9
‘The data are for liver hepatomas in females; data from the parental and Fl generations are pooled. These data are taken from Tables 3 and 4 of Tomatis et al. [ 1 I].
22
HARRY
-
observed of
0.W
AND
KENNY
S. CRUMP
frequency
liver
hepatomos
in female 05
A. GUESS
mice
I
0.0-L
1
50
0
100
150
200
I 250
d Dose (ppm)
FIG. 2.
which
maximizes
the likelihood
20 or less. The top curve,
function
labeled
over all a.m. polynomials
K = 4, comes
f4 (d) = 0.04483 + O.O02038d+ which
maximizes
the likelihood
function
curves
are nearly
coincident
of degree
the polynomial
(1.39 x 10-9)d4,
over all a.m. polynomials
4 or less. Here d is the dose rate, ppm of DDT of the mice. The three
from
of degree
in the diet over the life span
for doses
below
about
100 ppm.
The curve labeled K = 20 hooks up sharply at the highest dose, the upturn being caused by very high-order terms whose value is negligible except in the vicinity of the highest dose. Note that the low-order coefficients, which are the ones of importance for the low-dose extrapolation, are the same as the lower-order coefficients of the K= cc curve, up to 9 significant figures. In the low dose range [i.e., where the increased risk over background, P,(d)- PJO), is about 10m3 or less], the estimates of increased risk over background based on the curves K= 20 and K= cc agree to 9 significant figures. This illustrates with actual data the point discussed theoretically above, namely, that maximizing the likelihood function over all f E AMJ is equivalent to considering the limit of the maximizing polynomials as higher-
EXTRAPOLATION
and
OF ANIMAL
higher-degree
maximum
polynomials
is taken.
We will need
approach
and
the following
of Karlin
LE;IMA 3 (Karlin
are admitted
The former
cally equivalent to the latter implement computationally. theorem
CARCINOGENICITY
is much
lemma,
23
DATA
to the class
is both easier
which
over
theoretically to treat
is a minor
which
the
and practi-
theoretically modification
and of a
and Studden. and Studden)
Let O< d, < . . < d,, and suppose that c ,, . . ,c, are nonnegatice constants. Then a necessary and sufficient condition for there to exist an absolute& monotonic function g on [0, CC) mch that g(0) =O, g(d,) = c, for 1 Q i < n - 1, and g(d,J < c, is that
2 aid:+’ [ ,=I for any set of constants
>O
forallk>O]
*
[;,a,c,>O]
(8)
a,, . . ,a,.
Proof This follows from the proof of Theorem 10.3 on pp. 23 I-232 of Karlin and Studden [lo] because the functions zi,(k) = d,k+’ for i = 1.2,. . , n and k = 0, 1,2,. . . form a Tchebycheff system on the nonnegative integers, for the same reasons that the functions u,(k)= d,k do. Equation (8) above expresses the duality between the moment space and polynomial space for the system (6,) in the same way that Eq. (10.10) on p. 232 of Karlin and Studden expresses the duality for the system {u,}. THEOREM
I
Zf x, < n,, then there exists f EAMJ =sup[L(h)Jh
EAM].
such that L(f)=max[L(g)lgEAMJ]
Zf x, = n,, then for each ‘EAMJ,
L(f)
g)j gE
AMJI, Proof
We will give proof
for the case d, = 0; the proof
for the case d,, > 0
and let { f(“)}E=, be a sequence of is similar but easier. Let x,
2 a,dk+’ > 0 i=l
for all k > 0
1 i5 +
a,(c,--a)>0
,=I
I
HARRY
24
A. GUESS AND KENNY
S. CRUMP
So by Lemma 3 there exists a function g absolutely monotonic on [0, co) such that g(0) = 0, g(d,) = c, - ce for 1 < i < n - 1, and g(d,) $ c, - cg. Using Lemma 1, let Q EAMP be such that Q (~I~)=g(d~)+ ce for 0 < i < n, let c= c,- Q(d,) > 0, and letf(d)= Q(d)+c&,(d). ThenfEAMJ, andf(d,)= c, for 0 < j < n. By construction, L(f) =sup[L(h)lh EAM]. But sup[L(h)lh EAM]= sup[L( g)l g EAMJ] because the functions g,(d)= Q(d)+ ~(d/d,)~ are in AM and g,(d,)+_f(d,) for 0 $ j < n, as k-co. This proves the part of the theorem pertaining to x, < n,. When x,,= n, it is easy to see that sup[L( g)l g E AMJ] cannot be attained, because, for g(d) = Q(d) + ~6~ (d), L(g) is a strictly increasing function of c.
Our next theorem and corollary characterize tions to this maximization problem are unique. THEOREM
the extent to which solu-
2
Zff,gEAMJ
and L(f)=L(g)=sup[L(h)lhEAMJ],
for each
g(d,)=f(d,)
2
C
njg(d,)=
then
~ESP, (9)
"jf(d,>y
JESZ
jESZ
where
SP= [jlO< j< SZ=[jJO<
n and x,>O],
j
(10)
Proof As discussed earlier, maximizing L is equivalent to maximizing lnL, so we have lnL(f)=lnL( g)=sup[lnL(h)lhEAMJ]. Let G(y)=ln(le-~)andletO
lnL(ef+(l-/3)g)-sup[lnL(h)lhEAMJ] =lnL(Bf+(l-8)g)-[@lnL(f)+(l-B)lnL(g)] =
x
x,T,,
jESP
where
Since f and g maximize L, it is easy to see that f(4) > 0 and g(4) > 0 for each jESP. Since G is strictly convex on (0, co), we have Tj > 0 for each jESP, and if ,f(d,J#g(+o) for some j0 E SP, then TiO> 0. This would imply In L( of+ ( 1- (3) g) > sup[ln L( h)l h E AMJ], which is impossible. Hence f(d,) = g(4) for all j E SP, and since In L(f) = In L( g), Eq. (9) follows.
EXTRAPOLATION
OF ANIMAL
CARCINOGENICITY
DATA
25
Using Theorem 2 and some results of Krein and Rehtman [ 121, we may formulate conditions for the likelihood function to have a unique maximum in AMJ. To state the conditions we will need to introduce the notion, taken from approximation theory, of the index of a function in AMJ. The definition below is adapted from Krein and Rehtman [12, p. 1271. DEFINITION
Let ~EAMJ be given by f(d)= Q(d)+ c&(d), where Q(d)=q,+q,d’ + . . . + qKdK. The index of ,f is defined as follows. Let k, < k, < . . . < k,,, be the subscripts of coefficients q, such that q, > 0. Divide the numbers k, beginning with k,,,, into groups in the following manner. If k, - k,_, > 1, then the first group consists only of k,; otherwise the first group consists of the two numbers k, and k,,,_,. Remove the first group from the sequence
I
Let NP denote the number of subscripts j > 0 for which both 4 > 0 and xj > 0. Let f E AMJ maximize L otter AMJ, and let I(f) denote the index off as defined above. If Z(f) Q NP- 1. then f is the unique function in AMJ maximizing L ouer AMJ. If d,= 0 and either x0 > 0 or NP= n, then uniqueness holds wheneoer the index of the function q, + q,d’ + ’ * . + q,dk-’ + cC3,Jd) is less than or equal to NP - 1. Proof: We will follow the argument end we first introduce some notation.
in Krein and Rehtman [12]. To this Let j, < . .
26
HARRY
A. GUESS
such that both d, > 0 and x, > 0. Let jEAMJ AMJ], and let Z(j) Q NP - 1. We may write
where Q,E AMP and c, > 0. For dNp_ , = d, let Let
z= {O,l)...)
AND
KENNY
satisfy L(f)=
qf,
=
cJ;
S. CRUMP
sup[l(h)]h
otherwise let
qf,
E
=
0.
co},
s={lE.qqf>o}, and
where 6;,, = 1 if i = j and S,,, = 0 otherwise. In this notation we may write
c,=
c qjwo
O
Its
where c;=f($,) for 0 < i Q NP- 1. We now proceed to show that the argument in Krein and Rehtman [12] applies to our case. It is easy to see that our definition of the index of the function f coincides with the definition p. 127 in Krein and Rehtman of the index of the sequence {q, : 1 E S }. This in turn is the same as their definition on p. 128 of the index of the representation {c,}~~;‘. Thus what we call Z(f) is what they call the index of the representation {ci}yL;‘. Hence by Theorems 4.3 and 3.1 and Eq. 3.4 of Krein and Rehtman, since I(f) < NP - 1, there is a unique polynomial Q EAMP satisfying
Q(d,)=cc
O
Q(d,p,p_,)~CNP-l.
Since Q, also satisfies these conditions, QJ- Q. Now let g be any other function in AMJ such that L(g)= L(f)=sup[L(h)]h EAMJ]. We may write
EXTRAPOLATION
OF ANIMAL CARCINOGENICITY
DATA
27
where Q, EAMP and cp > 0. It follows from Theorem 2 that Q,(d,) = Q,(<,) =c, for 0~ i
forO
Let f&d)=q(+q:d’+‘-
+q’dK-’ K
’
To prove the second assertion it suffices to use the above argument with 5, playing the role of c, and using Q, in place of Q,. This completes the proof. When the maximum-likelihood estimatef fails to be unique, there will be a (possibly infinite) set of low-dose extrapolations each of which maximizes the likelihood function. In such cases we seek a canonical way to choose one of them. One way which is biologically meaningful and also easy to handle mathematically is to seek the maximum probability of additional response over background, P,(d,) - P,(O). at a given environmental dose rate d, that can be calculated using anyf EAMJ such that f.(f)=sup[L( g)l gE AMJ]. Whenever da=0 and a positive response is observed at this dose, it follows from Theorem 2 that P,(O) is unique. In this frequently occurring case, maximizing P,(d,)- P/(O) is equivalent to maximizing P/(d,). Since Pf(d,) is an increasing function of f(~$,), this is equivalent to solving the problem maximize
f( d, )
(11)
subject to ~EAMJ
and
L(f)=sup]l(g)lgEAMJ],
where d, is a given dose rate. Once we have obtained one function sup[l(h)lh EAMJ], it follows from Theorem (11) and (12) is equivalent to the problem maximize
(12)
g EAMJ such that L(g)= 2 that the problem defined by
5 l;d,’ r=O
(13)
HARRY
28
A. GUESS
AND
KENNY
S. CRUMP
subject to
;$oi.d:+cs, (4)=s(d,)t
jESP,
I51; jSSZ c n,d,i = c n&d,), 1jESZ
i=o
(14)
(W2
[
J;>o,c>o,
(16)
and KE[O,1,2,
. ..I.
(17)
where d, and g(d,) (0~ j < n) are given, d,< d,,, and L(g)=max[l(h)]h~ AMJ]. For each fixed value of K this is a linear-programming problem and therefore is readily solvable by standard algorithms. We will show in Theorem 3 that a global maximum (over all K) exists. It is possible to modify a standard linear-programming algorithm to yield an algorithm which would converge to the global maximum. For most practical purposes, however, such a modification appears unnecessary. In typical applications d,
EAMJI}.
(18)
if and only if X, < n,. In this case we have
3
If x, < n,
and d, < d,,, then there exists f E M such that f(d,)
= sup[ g(d,)l
gEMI. ProoJ: By Theorem
1, M is nonempty,
M = { flf
E AMJ andf
‘If n E SZ, then c = 0; we omit the proof.
so let g E M; then by Theorem satisfies (9)).
2
(19)
EXTRAPOLATION
OF ANIMAL CARCINOGENICITY
Since de < d,, we have f(d,) < f(d,J,
DATA
29
and it follows from (9) that
)ISE M 1
sup[f(d,
It follows from (20) that there exists a sequence f’“‘(d)= and such that the following
cw {f(m)}z3,
in M such that
Qcm)(d)+ccm)& n (d)
(21)
limits exist:
lim f’“‘(d,)=sup[h(d,)Ih m--t’x
EM l=c, (22) O
,Il= f(“‘( 4. ) = c,,
In the remainder of the proof we will assume that do=0 and d,# d, for each j. When the contrary is true the proof is similar but easier. Let g(“)(d)
= Q(“)(d)
- Qcm)(0); then gCrn)~ AM, gCm’(0) = 0,
and
g’“‘(d n ) < f’“‘(d n )-f’“‘(O). Hence by Lemma 3 and Eq. (22) i aj
all
k>O I
j=l n
2 a,(c,-c,)+a,(c,-co)>0
=+
[
Using Lemma 3 in the other direction, such that
we conclude
h (de>= ce- co,
h(4) = cj h(d,)
< c,-
co,
co.
O
1. that there exists h EAM
30
HARRY
By Lemma 1 there exists Q EAMP h(d,)+ ca for 0 < j < n. Let
A. GUESS
AND
KENNY
S. CRUMP
such that Q (d,)= h(d,)+ c,, and Q(4)=
c=c,-Q(d,,)=c,-c,-h(d,,)>O, and let (23) Then f E AMJ,
and for O
&-py~).
f(d,)=c,=
Since each f(“) satisfies (9) so does f. Hence f~ M, and f(d,)= hE M]. This completes the proof. 4.
GLOBAL
MAXIMIZATION
sup[h(d,)]
ALGORITHM
In this section we develop an algorithm for computing a function ~EAMJ such that L(f)=sup[l(h)]h EAMJ]. We outline the general structure of the algorithm and prove convergence. Each f E AMJ has the form
f(4=Q(d)+4in(d)a
(24)
where Q is a polynomial of arbitrary degree with nonnegative coefficients, c is a nonnegative number, and Sdn is the jump at d,, defined by (6). We may write Q in the form Q(d)=
q, > 0,
5 q;d’, r=O
(25)
where all but finitely many of the coefficients q, are zero. Since maximizing L(fl is equivalent to maximizing lnL(f), and since lnL(f) is a concave function of the vector (c,q)=(c,q,,q,, .), we will work with lnL(f) rather than with L(fl. Our problem may be stated as follows maximize
(26)
lnL(f)
subject to the constraints c>o
q, >O
for O< i< K,
and
q,=O
for i>K
(27)
EXTRAPOLATION
OF ANIMAL
CARCINOGENICITY
31
DATA
and K=O,l,2
,...,
(28)
where lnL(f) is given by (5) f 1sgiven by (24) and Q is given by (25). By Theorem 1 this problem has a solution if and only if X, < n,. For each fixed value of K, (26) and (27) define a concave nonlinear-programming problem with linear constraints, and hence the Kuhn-Tucker conditions are necessary and sufficient for a vector (c, qO,. . . , qK) to yield a maximum. For future reference we state this formally as a lemma. LEMMA
4 (Kuhn-Tucker)
For each fixed K the following conditions are necessary and sufficient for a function f to soke the nonlinear-programming problem defined by (26) and (27), where lnL(f) is g&en by (5), f by (24), Q by (25), and P by (2): aInL(f)
-i
84,
[-n,+&]$
=0
if
q,>O, (29)
J
J=o
GO
if
q,=O
for all 0 < i < K, and
alnL(f) -[
X,
ac
-
-%+
P(d,)
=o
if
c >O,
GO
if
c=O.
I
(30)
ProojI Equations (29) and (30) are the Kuhn-Tucker conditions for the nonlinear-programming problem defined by (26) and (27). Since the constraints are linear and the function to be maximized is concave, the Kuhn-Tucker conditions are necessary and sufficient. (See Zangwill [13], Theorems 2.14 and 2.15 and Exercise 2.13.) We will also have occasion to refer to the nonlinear Kth-dimensional problem where the jump is not allowed at the highest dose. This problem may be stated maximize
lnL(f)
(31)
0 d i Q K,
(32)
subject to the constraints Aa0
for
where lnL(f) is given by (5) and f is given by (3). This is also a concave nonlinear-programming problem, and the corresponding necessary and
32
HARRY
sufficient
Kuhn-Tucker alWf> 56
conditions = i
A. GUESS
AND
S. CRUMP
for a function f to be a solution
[ -n,+&]$=O
are
if J;>O, J
j=O
GO
for all 0~ i < K, where lnL(f) LEMMA
KENNY
is given by (5) andf
if
J=O
is given by (3).
5
In order that f EAMJ sati.& In L(fl=sup[lnL(h)/h sary and sufficient that (29) and (30) hold for all i > 0.
EAMJ]
it is neces-
Proof: Necessity is an obvious consequence of Lemma 4. If the conditions hold for all i > 0, then from (30) x, < n,, and so, by Theorem 1, there exists g EAMJ such that lnL(g)=sup[lnL(h)Jh EAMJ]. We may write g(d)= 0 (d)+ E&(d) where Q is an a.m. polynomial. It follows from Lemma 4 that lnL(f) > lnL(g). Hence lnL(f)=sup[lnL(h)]h EAMJ]. This proves sufficiency.
The obvious analog of this lemma holds for the case where no jump at d, is permitted. In this case a global maximum need not exist, but (29) still is necessary and sufficient. The exact result is as follows: COROLLARY In order that QEAMP satisb lnL(Q)=sup[lnL(Q)lQEAMP] necessary and sufficient that (29) hold for all i > 0. (Here
c =0
it is in all
formulas.). Proof Necessity is obvious. To prove sufficiency (29) (and that c = 0). Then
lnL(Q)=sup[lnL(Q)IQEAMP,(degreeof
suppose that Q satisfies
Q)
(33)
holds for each K, by the Kuhn-Tucker theorem. Since (29) holds for all - n,, + x,/ P(d,,) < 0 and hence x, < n,. By Theorem 1 there exists gEAMJ such that lnL(g)=sup[lnL(h)]h EAMJ]. Let g(d)= &(d)+ c&,(d), where Q EAMP. For each k let Qk be defined by Qk(d)= o(d)+ c(d/dJk. Then Q,EAMP and limk_,Qk(<)=g($ for O< j < n, and hence lim ,,,lnL(Q,)=lnL( g). It follows from (33) that i > 0, we have
lnL(Q,)g for each k, and so lnL(Q)=
lnL(Q)
< lnL(g)
lnL( g). This proves sufficiency.
We now have a set of infinitely many conditions which together are necessary and sufficient for a function to globally maximize the likelihood
EXTRAPOLATION
OF ANIMAL CARCINOGENICITY
DATA
33
function. The key step in the development of our algorithm is to replace this set by a set which has only finitely many conditions and yet is necessary and sufficient. To do this we will introduce the following lexicographic order relation. DEFINITION
A vector (u,, . . . , u,) will be said to be iexico-negative if the first nonzero term in the sequence u,,, u,_ ,, . . . , u. (read starting with u,) is negative. [Thus the vectors (O,O,. . . ,O) and (2,2,. . , ,2, - 1) are lexico-negative.] LEMMA
6
(a) The vector (u,, . . . un) defined by uj =
I
-T+
xj P(4)
(34) 1
is lexico-negative if and only if there exists some integer i, such that ClIn L(f)/ aq,< 0 for all i > i,. (b) If (u,, . . , u,,) defined by (34) is Mexico-negative, then 8 In L( f)/ aq, < 0 for all i > i,, where i, is given by i,=O
if all
U, GO,
OG i < n,
(35)
and otherwise
In i,=l+
j@axo~j~jo-I[-9+X,/p(4)l I-ni0+x/0/p(40))) I
40 In &-I i-1
(36)
where j, is the index of the first nonzero term in the sequence u,, u,_ I,...,t+l (starting with u,). Proof Both (a) and (b) follow from straightforward will be omitted. THEOREM
computations,
which
4
The following conditions are necessary and sufficient for In L(f) = sup[lnL(g)(gEAMJ], where InL(f) is given by (5), f by (24), Q by (25), and P by (2): q,=O
foraIli>
K,
for some K,
(uo,..., u,) is lexico-negative, Eqs. (29) and (30) hold for all i < K,
(37) (38) (39)
HARRY
34
A. GUESS
AND
KENNY
S. CRUMP
and
alnuf)
ho
forall
K+l
a4i
Here the uj are defined by (34). and i, is defined i, < K + 1, then (40) is CacuousIy satisfied.)
by (35) and
(36).
(If
Proof. We first prove necessity. Let ~EAMJ, and let lnL(f)= sup[lnl(h)lh EAMJ]. Then (39) holds by Lemma 4 and, since f EAMJ, we have (37). It follows from Lemmas 5 and 6 that (38) holds. By (37) and Lemma 5 we get (40). This proves necessity. We now prove sufficiency. By (38), (40) and Lemma 6 we have a lnL( f)/ aq, G 0 for all i > K + 1. This, together with (37) and (39) implies that the conditions of Lemma 5 are satisfied. Hence we have sufficiency. COROLLARY
Equations (37) through (40) omitting the requirement of (30) in Eq. (39) and taking c = 0 throughout, are necessary and sufficient conditions for there to exist
Q EAMP
such that lnL(Q)=sup[lnL((i)lQ
The proof is similar to that of Theorem
Proof
Theorem 4 forms the basis for our algorithm. algorithm is as follows:
EAMP]. 4 and is omitted. The general scheme of the
(a) Starting with a given K, solve the nonlinear-programming problem defined by (26) and (27). (b) Once a maximum for fixed K has been obtained in step (a), Eqs. (38) and (40) are checked. If both are satisfied, then (by Theorem 4) the current value of (c, qO,. . . , qk) is globally optimal and the algorithm stops. If either (38) or (40) fails to hold, we replace K by K+ 1 and return to step (a), using our current vector (c, q,,, . . , qk, 0) as a starting point for performing the new maximization. Since a global maximim exists when X, < n, (by Theorem 1) and since the conditions checked by the algorithm each time K is increased are necessary and sufficient for a vector (c,q,, . . ,qk) to yield a global maximum, the algorithm will eventually converge to a global maximum. This proves convergence of the algorithm to a global maximum. We now state briefly the method we have used to solve the finite-dimensional nonlinear-programming problem defined by (31) and (32), i.e., step (a) above. We will need the following lemma. LEMMA
7
Suppose for a fixed K that a solution f’ to the nonlinear-programming problem defined by (31) and (32) satisfies the condition x,/n, > P(f’(d,)). If
EXTRAPOLATION
OF ANIMAL
CARCINOGENICITY
f* is any solution to the same nonlinear-programming
replaced by n - 1, then x,/n,
35
DATA
problem except that n is
> P(,F(d,,)).
Proof: Let
and
We then have e,(f’)98,_,(_~‘)<8,_,(f*) and e,,(f’)=e,_,(f’)+a,(f’). Now suppose the conclusion of the theorem does not hold, that is, x,/n, < P(f*(d,)), and define g,=cuf*+(l -LY)~), 0~ a g 1. It is easy to see that g,EAMP (O
i.e.,
Qn(g,,)>Q,,(f’). and This completes
IS a maximum.
Now suppose we have solved the nonlinear-programming fined by (3 1) and (32) and obtained the solution f’ EAMP,
this conthe proof.
problem dewhich must
the Kuhn-Tucker conditions (29). If, in addition, x,/n, < seen that f(d) =J’(d) + c6dn(d) EAMJ, where c = 0, must satisfy the Kuhn-Tucker conditions (29) and (30) for the problem defined by (26) and (27). If, instead, x,/n, > P(f’(d,,)); we proceed to obtain f* EAMP that solves the problem defined by (31) and (32) with n necessarily
P(f’(d,
replaced
satisfy
)), it is easily
by
n - 1.
We
are
x,/n, > P(fc(d,)). It follows c&(d) EAMJ, where c= x,/n,
guaranteed easily
from
by this
the
above
inequality
- P(fc(d,, )) satisfies
lemma
that
that f(d)=f*+
the conditions
(29) and
(30) and thus is a solution to the problem defined by (26) and (27). With this procedure the two problems defined by (26) and (27) and by (31) and (32) respectively can both be solved using a single algorithm. The algorithm we have used for this purpose is basically the method of Davidon, Fletcher, and Powell; see Fletcher and Powell [14]. It was necessary to modify this method, which applies only to unconstrained optimization problems, to take into account our nonnegativity constraints.
36 5.
HARRY
FUTURE
A. GUESS
AND
KENNY
S. CRUMP
WORK
In a subsequent paper we will investigate the statistical aspects of this technique. We will show that in general the maximum-likelihood estimates of the response probabilities P(d) are strongly consistent and unbiased for all d< d,,. Under some additional hypotheses the estimates will be shown to be asymptotically normal. Using asymptotic normality, confidence intervals will be computed and compared with simulation results. Finally, risk estimates obtained using this technique will be compared with risk estimates using several widely known techniques. We thank Kevin Deal for his efficient programming tions regarding computational methods.
and for useful sugges-
REFERENCES D. G. Hoel, D. W. Gaylor, R. L. Kirschstein, U. Saffiotti, and M. A. Schneiderman, Estimation of risks of irreversible, delayed toxicity, J. Toxicol. and Enoiron. Health 1, 133-151(1975). R. Gillette, Cancer and the environemnt (II): Groping for new remedies, Science 186, 242-245 (1974). N. Mantel and W. R. Byran, Inst. 27, 455470 (1961).
Safety
testing
of carcinogenic
agents,
J. Nail. Cancer
N. Mantel, N. Bohidar, C. Brown, J. Ciminera, and J. Tukey, An improved “MantelBryan” procedure for “safety” testing of carcinogens, Cancer Res. 35, 865-872 (1975). W. S. Abbott, A method of computing the effectiveness of an insecticide, J. Econ. Entomol. 18, 265-267 (1925). P. Armitage and R. Doll, Stochastic models for carcinogenesis, Proc. Fourfh Berkeley Symp. Math. Stat. Probab. Univ. of Calif. Press, Berkeley, 1961, Vol. 4, pp. 19-38. K. S. Crump, D. G. Hoel, C. H. Langley, and R. Peto, Fundamental carcinogenic processes and their implications for low-dose risk assessment, Cancer Res. (1976), to 8 9 10 11 12
13 14
appear. C. C. Brown, Statistical aspects of extrapolation of dichotomous dose-response data, submitted for publication. D. V. Widder, The Laplace Transform, Princeton U. P., Princeton, 1941. S. Karlin and W. J. Studden, Tchebycheff Systems: With Applications in Analysis and Statistics, Wiley, New York, 1966. L. Tomatis, V. Turusov, N. Day, and R. T. Charles, The effect of long-term exposure to DDT on CF-1 mice, Inr. J. Cancer 10, 489-506 (1972). M. G. Krein and P. G. Rehtman, Development in a new direction of the CebysevMarkov theory of limiting values of integrals, Usp. Mat. Nauk (N.S.) 10, No. 1 (63), 67-78 (1955); Am. Math. Sot. Transl. Ser. 2, 12, 123-135 (1959). W. I. Zangwill, Nonlinear Programming, A Unified Approach, Prentice-Hall, Englewood Cliffs, N. J., 1969. R. Fletcher and M. J. D. Powell, A rapidly convergent descent method for minimization, Comput. J. 6, 163-168
(1963).