- Email: [email protected]
A parametric independence test for clustered binary data
STATISTICS& ELSEVIER
Statistics & Probability Letters 41 (1999) 1 7
A parametric independence test for clustered binary data Dale Bowman * Department of Mathematics, University of Mississippi, University, MS 38677, USA Received February 1997; received in revised form January 1998
Abstract
This paper proposes an independence test for a set of clustered binary observations, such as might be encountered in developmental toxicity studies. An exchangeable binary model is employed, under an assumption of exchangeability among cluster elements, to model probability of positive response. With a Weibull form assumed for response, the independence test is equivalent to testing whether a parameter value is unity. This Weibull form also allows for parametric tests for a covariate-response relationship and for covariate effects on correlation. The procedure is illustrated using data obtained from a developmental toxicity study. (~) 1998 Elsevier Science B.V. All rights reserved Keywords:
Exchangeability; Developmental study; Correlated data
1. Introduction
In many experimental situations binary responses may be elicited in clusters. It is often noted that responses among individuals within the same cluster are dependent while responses across clusters are independent. Such an assumption is usually based on a priori belief and on observable overdispersion relative to a model based upon an assumption o f independence. George and Kodell (1996) develop a likelihood ratio test for independence based on an assumption of exchangeability, which is realistic in many situations. Using the definition of exchangeability given by de Finetti (1974), a set of random variables, X1,X2 . . . . . is exchangeable if for any n, and any vector ( x l . . . . . x n ) , P(X~(1) : X l , . . .
,Xrc(n) =
Xn ) = P ( X , = x l , . . . , X n = xn ),
(1)
for any permutation n(1) . . . . . rffn) of 1. . . . . n. Bowman and George (1995) develop maximum likelihood estimation procedures for response rates for clustered binary data under an assumption of exchangeability. Specifically, suppose X = ( X 1 . . . . . X n ) is a
* Tel.: +1 601 232 7409; e-mail: [email protected]. 0167-7152/99/$ - see front matter (~) 1998 Elsevier Science B.V. All rights reserved PII: S01 6 7 - 7 1 5 2 ( 9 8 ) 0 0 0 8 6 - 8
2
D. Bowman I Statistics & Probability Letters 41 (1999) 1 - 7
cluster of n exchangeable binary random variables. Define 2k =P(X1 = 1..... Xk = 1).
(2)
Note that since the random variables are exchangeable, 2k will be the same for any subset of k of the random variables. Note also that 21 is the marginal response rate. Define R to be the sum of the X/'s. Using an inclusion, exclusion argument, it can be shown that
P(R=r)=
(7)
Z(_l)k
n
k=0
r 2r+k.
(3)
Using this model, Bowman and George (1995) develop maximum likelihood estimates of the parameters 2k and George and Kodell (1996) construct a likelihood ratio test of independence among cluster elements based on the fact that if the observations are independent the probability distribution given above reduces to the ordinary binomial distribution. An alternative approach to testing for independence is derived here by assuming a Weibull-type model for the marginal probabilities 2k. This allows for a parametric test for independence. The estimation of the parameters of the Weibull form is based on the procedure developed in George and Bowman (1995) where a folded logistic form is assumed for 2k. The Weibull-type model is developed and hypothesis tests are derived in Section 2. Section 3 discusses the estimation procedure and Section 4 illustrates the model and testing procedures on data from a developmental study.
2. Parametric model and a test for independence
Let Xl ..... Xn be a cluster of exchangeable binary random variables and let 2k be as before, where, due to the assumption of exchangeability, 2k has the same value for any subset of the random variables of size k. Then it can be shown that the probability that ~ i n l X, = r is given by (-1)k k=O
() n-r k
2r+k.
(4)
George and Bowman (1995) term this the exchangeable binary model. In this model the marginal response rate of individual cluster elements is 21. Correlations between cluster elements may be obtained using this model. In many instances primary interest is in correlations of order two, but correlations of all order may be obtained. Defining Pk as (var(X/))k/2pk = E (0(1 - #)(X2 - p) .... (Ark -- g)),
(5)
where E(X,)= 21 = #, an expression for Pk is given by E j = o ( - 1 ) k-' (~)2~-sJ.y Pk = [21(1 - - /~1)] k/2 k
•
"
(6)
Specifically, the correlation of order 2, often called intra-cluster correlation, P2 is given by (22 - 22) P2 - - ~ 1 ( 1 _ 21 )
(7)
D. Bowman I Statistics & Probability Letters 41 (1999) 1-7
3
In order to estimate a response curve using the exchangeable binary model, George and Bowman (1995) assume a folded logistic model for 2k. That is ,~ =
2 1 + efl~g(k)+fl2zg(k)'
(8)
where fll and t2 are parameters, z is a covariate, and 9(k) is a function of k specified prior to estimation. George and Bowman selected g(k) to be log(k + 1). The resulting expression for 2k is a completely monotonic function when the exponent fllO(k)+ fl2zg(k) is non-negative. The restriction of complete monotonicity on the form assumed for 2k is necessary to ensure that all estimated response probabilities are non-negative. The parameters of the folded logistic model control the effect of covariates on both mean response, 21, and higher-order correlation structures. In some cases this may not be flexible enough to capture all effects of covariates. A more flexible model in this respect is the Weibull-type model where 2k = e -°k',
(9)
where 0 and 7 are functions of covariates. Using the Weibull-type model, mean response is given by 21 = e -°, demonstrating that only 0 influences mean response, while both 0 and a influence the correlation structure. Based on the definition of 2k it is clear that the exchangeable random variables are independent if and only if 2k = (21)k. Under the assumption of independence, the exchangeable model reduces to the binomial distribution. Using the Weibull form, 2k = (21)k if and only if a = 1. Hence a test for independence of the exchangeable observations reduces to a hypothesis test of H0 : a = 1 versus H1 : a ~ 1. It can also be shown that positive correlation among cluster elements corresponds to the case in which a < 1, allowing for a test for positive correlation to be constructed as well.
3. Maximum likelihood estimation Let nl ..... nm be the cluster sizes of m clusters of exchangeable binary observations. Further let x~..... Xm be the number of positive responses within each exchangeable cluster. Further let zi and vi be vectors of covariates effecting mean response and correlation, respectively. Further let 0 and a be vectors of parameters. Then the likelihood function of 0 and a may be written
Z(O,a)= i= 1
ni xi
Z(_I)j
ni~xi
e_O,z(xi+j),,'~
.
(10)
j=0
Maximum likelihood estimates of 0 and a are obtained using an algorithm based on the Nelder Mead simplex method (Nelder and Mead, 1965). The simplex method obtains the minimum of a function of n parameters (here the function is the negative log-likelihood) by comparing the function values at n + 1 vertices of a general simplex. The vertex associated with the highest function value is replaced by a new point in each iteration by contraction, extension or reflection, ultimately achieving the function minimum. The advantage of the simplex method over Newton-Rhapson-type procedures is that no evaluation of gradients or Hessians is required, making the estimation process less sensitive to ill-conditioning. Estimates of the covariance matrix for the parameter estimates are obtained using a fit quadratic surface as suggested by Nelder and Mead (1965). Using a transformation of the form Cov = D r V D where D is a matrix of partial derivatives and V is a covariance matrix, estimates of standard errors of 21 and p may be obtained. Estimates of 0 and a must be such that the Weibull-type model is completely monotonic. At the minimum, this requires that both O'z and a'v be non-negative. More precise requirements for complete monotonicity depend on the covariates z and v. Constrained optimization may be required to force the parameter estimates to remain in the feasible region. With the simplex method the estimates can be constrained by assigning
4
D. Bowman I Statistics & Probability Letters 41 (1999) 1-7
arbitrarily large function values to points outside the feasible region. The simplex method then replaces the non-feasible point within that iteration.
4. Example The procedure is illustrated using a data set from a developmental toxicity study designed to measure the teratogenic effect of diethylene glycol dimethyl ether (DYME) in CD-1 mice (Price et al., 1987). In this experiment five dose groups of pregnant females were administered dose levels of 0,0.0625,0.125,0.25, and 0.5 g/kg of body weight. A dose scaling of log(dose + 1) was used to facilitate the estimation procedure. Each litter represents a cluster and it is recognized (e.g. Williams, 1975; Ryan 1992) that littermates do not respond independently. However, George and Bowman (1995) assume exchangeability of littermate response, in which case the exchangeable binary model is appropriate. The binary endpoint under consideration in this analysis is fetal death. Table 1 gives the frequency distribution of litter size and number of deaths within each litter. The model considered for 2 is
~i)= exp(-(Ol ÷
(11)
02di) * (k)a'+az*di),
where di is the dose level in the ith dose group, 01 + 02di>~O, and al +a2di>>-Ofor i = 1..... 5. This model includes a dose term in both the mean response covariate vector and the correlation covariate vector. This allows for a test of dose effect on both effects as follows. If 02 = 0 then there is no dose effect on mean response rate. If 02 = 0 and a2 is non zero then there is an effect of dose on correlation only. It is also possible with this model to test for independence within each dose group by testing H0:al + a2di = 1 for i = 1. . . . . 5. The simplex estimation procedure described in Section 4 resulted in the following estimates of the parameters in Eq. (11 ), &~= (3.1947, -6.0114, 1.0038, -0.5247) with covariance matrix
A
Cov(~) =
0.1755
-0.0874
0.1167
-0.0061
-0.0874
0.5243
0.1439
-0.1408 /
0.1167
0.1439
0.1749
-0.0607 /
-0.0061
-0.1408
-0.0607
0.0461 /
for o~t=(Ol,O2,al,a2). Corresponding z-scores for the parameter estimates are 7.62,-8.30,2.40 and -2.44, respectively. Note that all estimates are significant indicating a dose effect on both mean and correlation. In Table 2, the corresponding estimates of mean response, 2~ and correlation, P2 are given along with observed response rates. These estimates are compared with the maximum likelihood estimates obtained using the non-parametric estimation procedure of Bowman and George (1995). No dose response function is assumed for the non-parametric estimates and these estimates of mean response in Table 2 follow the observed response rates closely, including the downturn in the second dose group following control. An inherent assumption in the Weibull-type model is that response increases with dose. This accounts for the modest differences between the Weibull-type model and the non-parametric procedure in estimates of mean response. The non-parametric estimates for P2 differ from those obtained using the Weibull model as well. The standard errors, given in parenthesis, of the estimates of correlation for both the Weibull model and the non-parametric estimates are larger than the estimated correlations in all dose groups including control. Fig. 1 shows a graph of the estimated dose response curve for the Weibull type model in Eq. (11). The Weibull-model seems to provide reasonable fit to the observed data, shown on Fig. 1 as the discrete dots. A goodness-of-fit test for the Weibull model resulted in a chi-square value of 122.99 with 1 degree of freedom, indicating no lack of fit.
D. Bowman I Statistics & Probability Letters 41 (1999) 1-7 Table 1 Frequency distribution of female mice according to littersize and the number of dead fetuses following exposure to diethylene glycol dimethyl ether (DYME) Price et al. (1987) Dose (g/kg)
No. deaths
Littersize 3
0
7
8
9
10
11
0
1
1
1
1
1
2 3 0.0625
1 1
1
1
2 0.125
0
1
1
15
16
17
4
2
2
1
1
l
2
2
1 1
1
3
1
2
2
1
1
2
1
2
1
2
3
1
2 1
1
1
2
1
1
2 1
2 3 4 7 8 1 2 4 5 6 7 9 11 13
1
1
2
0 1
0.5
14
1
2 3 0.250
13
1
0 1
12
2
1
3
2
1
1 1
1
1 1
2
1 1
1 1 1 1 1 1
1 1
1
2
1 1 1
1 2 1 1 1 1
1 1
2 1
A likelihood ratio test is constructed as an overall test for independence in the data. A Weibull-type model with 2 k = e x p ( - ( 0 1 + 0 2 d i ) * k 1)
(12)
is fit to the data. The model in Eq. (12) is appropriate if the data are independent. In the D Y M E data a likelihood ratio test based on the models in Eqs. (11 ) and (12) was constructed to test for overall independence in the data. The test statistic was 12.097 with 2 degrees o f freedom, resulting in a p - v a l u e o f 0.0024. Thus the L R T for i n d e p e n d e n c e is rejected in the D Y M E data. The estimated values o f al + azdi for the four dose groups were all less than 1, indicating positive correlation in these dose groups. In the control group the estimate was 1.004 which m a y indicate independence a m o n g responses there. A test for independence versus positive correlation within each dose group was derived. Specifically, the null hypothesis, H0 : al + a2 * di = 1 versus H1 : al + a2 * di < 1 is tested for each dose group. The z-scores for each dose group, respectively, were found to be 0 . 0 0 9 , - 0 . 0 6 8 , - 0 . 1 4 4 , - 0 . 2 9 2 , a n d - 0 . 5 7 2 4 , n o n e o f which were significant. For this data, the hypothesis o f independence is not rejected in any dose group, although
6
D. Bowman / Statistics & Probability Letters 41 (1999) 1 - 7
Table 2 Estimates of mean response and correlation in the DYME data i
Observed rates
1 2 3 4 5
0.0473 0.0744 0.0705 0.1271 0.5053
Weibull model
0.0410 0.0590 0.0832 0.1567 0.4689
Non-parametric
(0.0172) (0.0256) (0.0375) (0.0769) (0.2704)
0.0488 0.0877 0.0705 0.1202 0.5041
(0.7468) (0.5820) (0.4397) (0.2731) (1.2234)
0.0179 -0.0564 -0.0124 0.1125 0.1418
p2
(0.0142) (0.0188) (0.0135) (0.0287) (0.0481)
P2
-0.0007 0.0072 0.0196 0.0600 0.2001
(0.0950) (0.1689) (0.0832) (0.1335) (0.3433)
/
ml(d Yi
I --
cl, xi I
model I observed
Fig. 1.
the overall test rejected independence. These results correspond with the estimates of intralitter correlation, P2 given in Table 2 which clearly tend to increase with increasing dose although none are significantly different from zero. This indicates the presence of a dose effect on correlation as was shown by the significance of the estimate of a2. Other likelihood ratio tests may be constructed to test for dose effects. An overall test for dose effect is constructed by comparing the likelihood of a model of the form 2k z e -O~*(k)~
for all dose groups with the full model specified as Eq. (11 ). The results of this test for the DYME data gave a test statistic of 108.38 with 2 degrees of freedom. The null hypothesis of no dose effect is overwhelmingly rejected.
D. Bowman / Statistics & Probability Letters 41 (1999) 1 - 7
7
5. Discussion The model developed here for exchangeable binary observations has the desirable property of allowing for a parametric test for independence. A likelihood ratio test for overall independence across all dose groups was compared with an independence test within each dose groups. Using this method it was shown that for the DYME, an overall test of independence across all dose groups rejected the independence hypothesis, although tests performed separately within each dose group failed to reject the hypothesis of independence in any dose group. This effect may be caused by the significance of ~2 in the first model indicating correlation between littermates is influenced by dose. References Bowman, D., George, E.O., 1995. A saturated model for analyzing exchangeable binary data. J. Amer. Statist. Assoc. 90, 871-879. de Finetti, B., 1974. Theory of Probability. Wiley, New York.
George, E.O., Bowman, D., 1995. A full likelihood procedure for analyzing exchangeable binary data. Biometrics 51, 512-523. George, E.O, Kodell, R., 1996. An independence test for exchangeable data, submitted. Nelder, J.A., Mead, R., 1965. A simplex method for function minimization. Comput. J. 7, 308-313. Price, C.J., Kimmel, C.A., George, J.D., Mart, M.C., 1987, The developmental toxicity of diethylene glycol dimethyl ether in mice. Fundam. Appl. Toxicol. 8, 115-126. Ryan, L., 1992. Quantitative risk assessment for developmental toxicity. Biometrics 48, 163-174. Williams, D.A., 1975. The analysis of binary responses from toxicological experiments involving reproduction and teratogenieity. Biometrics 31, 949-952.
Statistics & Probability Letters 41 (1999) 1 7
A parametric independence test for clustered binary data Dale Bowman * Department of Mathematics, University of Mississippi, University, MS 38677, USA Received February 1997; received in revised form January 1998
Abstract
This paper proposes an independence test for a set of clustered binary observations, such as might be encountered in developmental toxicity studies. An exchangeable binary model is employed, under an assumption of exchangeability among cluster elements, to model probability of positive response. With a Weibull form assumed for response, the independence test is equivalent to testing whether a parameter value is unity. This Weibull form also allows for parametric tests for a covariate-response relationship and for covariate effects on correlation. The procedure is illustrated using data obtained from a developmental toxicity study. (~) 1998 Elsevier Science B.V. All rights reserved Keywords:
Exchangeability; Developmental study; Correlated data
1. Introduction
In many experimental situations binary responses may be elicited in clusters. It is often noted that responses among individuals within the same cluster are dependent while responses across clusters are independent. Such an assumption is usually based on a priori belief and on observable overdispersion relative to a model based upon an assumption o f independence. George and Kodell (1996) develop a likelihood ratio test for independence based on an assumption of exchangeability, which is realistic in many situations. Using the definition of exchangeability given by de Finetti (1974), a set of random variables, X1,X2 . . . . . is exchangeable if for any n, and any vector ( x l . . . . . x n ) , P(X~(1) : X l , . . .
,Xrc(n) =
Xn ) = P ( X , = x l , . . . , X n = xn ),
(1)
for any permutation n(1) . . . . . rffn) of 1. . . . . n. Bowman and George (1995) develop maximum likelihood estimation procedures for response rates for clustered binary data under an assumption of exchangeability. Specifically, suppose X = ( X 1 . . . . . X n ) is a
* Tel.: +1 601 232 7409; e-mail: [email protected]. 0167-7152/99/$ - see front matter (~) 1998 Elsevier Science B.V. All rights reserved PII: S01 6 7 - 7 1 5 2 ( 9 8 ) 0 0 0 8 6 - 8
2
D. Bowman I Statistics & Probability Letters 41 (1999) 1 - 7
cluster of n exchangeable binary random variables. Define 2k =P(X1 = 1..... Xk = 1).
(2)
Note that since the random variables are exchangeable, 2k will be the same for any subset of k of the random variables. Note also that 21 is the marginal response rate. Define R to be the sum of the X/'s. Using an inclusion, exclusion argument, it can be shown that
P(R=r)=
(7)
Z(_l)k
n
k=0
r 2r+k.
(3)
Using this model, Bowman and George (1995) develop maximum likelihood estimates of the parameters 2k and George and Kodell (1996) construct a likelihood ratio test of independence among cluster elements based on the fact that if the observations are independent the probability distribution given above reduces to the ordinary binomial distribution. An alternative approach to testing for independence is derived here by assuming a Weibull-type model for the marginal probabilities 2k. This allows for a parametric test for independence. The estimation of the parameters of the Weibull form is based on the procedure developed in George and Bowman (1995) where a folded logistic form is assumed for 2k. The Weibull-type model is developed and hypothesis tests are derived in Section 2. Section 3 discusses the estimation procedure and Section 4 illustrates the model and testing procedures on data from a developmental study.
2. Parametric model and a test for independence
Let Xl ..... Xn be a cluster of exchangeable binary random variables and let 2k be as before, where, due to the assumption of exchangeability, 2k has the same value for any subset of the random variables of size k. Then it can be shown that the probability that ~ i n l X, = r is given by (-1)k k=O
() n-r k
2r+k.
(4)
George and Bowman (1995) term this the exchangeable binary model. In this model the marginal response rate of individual cluster elements is 21. Correlations between cluster elements may be obtained using this model. In many instances primary interest is in correlations of order two, but correlations of all order may be obtained. Defining Pk as (var(X/))k/2pk = E (0(1 - #)(X2 - p) .... (Ark -- g)),
(5)
where E(X,)= 21 = #, an expression for Pk is given by E j = o ( - 1 ) k-' (~)2~-sJ.y Pk = [21(1 - - /~1)] k/2 k
•
"
(6)
Specifically, the correlation of order 2, often called intra-cluster correlation, P2 is given by (22 - 22) P2 - - ~ 1 ( 1 _ 21 )
(7)
D. Bowman I Statistics & Probability Letters 41 (1999) 1-7
3
In order to estimate a response curve using the exchangeable binary model, George and Bowman (1995) assume a folded logistic model for 2k. That is ,~ =
2 1 + efl~g(k)+fl2zg(k)'
(8)
where fll and t2 are parameters, z is a covariate, and 9(k) is a function of k specified prior to estimation. George and Bowman selected g(k) to be log(k + 1). The resulting expression for 2k is a completely monotonic function when the exponent fllO(k)+ fl2zg(k) is non-negative. The restriction of complete monotonicity on the form assumed for 2k is necessary to ensure that all estimated response probabilities are non-negative. The parameters of the folded logistic model control the effect of covariates on both mean response, 21, and higher-order correlation structures. In some cases this may not be flexible enough to capture all effects of covariates. A more flexible model in this respect is the Weibull-type model where 2k = e -°k',
(9)
where 0 and 7 are functions of covariates. Using the Weibull-type model, mean response is given by 21 = e -°, demonstrating that only 0 influences mean response, while both 0 and a influence the correlation structure. Based on the definition of 2k it is clear that the exchangeable random variables are independent if and only if 2k = (21)k. Under the assumption of independence, the exchangeable model reduces to the binomial distribution. Using the Weibull form, 2k = (21)k if and only if a = 1. Hence a test for independence of the exchangeable observations reduces to a hypothesis test of H0 : a = 1 versus H1 : a ~ 1. It can also be shown that positive correlation among cluster elements corresponds to the case in which a < 1, allowing for a test for positive correlation to be constructed as well.
3. Maximum likelihood estimation Let nl ..... nm be the cluster sizes of m clusters of exchangeable binary observations. Further let x~..... Xm be the number of positive responses within each exchangeable cluster. Further let zi and vi be vectors of covariates effecting mean response and correlation, respectively. Further let 0 and a be vectors of parameters. Then the likelihood function of 0 and a may be written
Z(O,a)= i= 1
ni xi
Z(_I)j
ni~xi
e_O,z(xi+j),,'~
.
(10)
j=0
Maximum likelihood estimates of 0 and a are obtained using an algorithm based on the Nelder Mead simplex method (Nelder and Mead, 1965). The simplex method obtains the minimum of a function of n parameters (here the function is the negative log-likelihood) by comparing the function values at n + 1 vertices of a general simplex. The vertex associated with the highest function value is replaced by a new point in each iteration by contraction, extension or reflection, ultimately achieving the function minimum. The advantage of the simplex method over Newton-Rhapson-type procedures is that no evaluation of gradients or Hessians is required, making the estimation process less sensitive to ill-conditioning. Estimates of the covariance matrix for the parameter estimates are obtained using a fit quadratic surface as suggested by Nelder and Mead (1965). Using a transformation of the form Cov = D r V D where D is a matrix of partial derivatives and V is a covariance matrix, estimates of standard errors of 21 and p may be obtained. Estimates of 0 and a must be such that the Weibull-type model is completely monotonic. At the minimum, this requires that both O'z and a'v be non-negative. More precise requirements for complete monotonicity depend on the covariates z and v. Constrained optimization may be required to force the parameter estimates to remain in the feasible region. With the simplex method the estimates can be constrained by assigning
4
D. Bowman I Statistics & Probability Letters 41 (1999) 1-7
arbitrarily large function values to points outside the feasible region. The simplex method then replaces the non-feasible point within that iteration.
4. Example The procedure is illustrated using a data set from a developmental toxicity study designed to measure the teratogenic effect of diethylene glycol dimethyl ether (DYME) in CD-1 mice (Price et al., 1987). In this experiment five dose groups of pregnant females were administered dose levels of 0,0.0625,0.125,0.25, and 0.5 g/kg of body weight. A dose scaling of log(dose + 1) was used to facilitate the estimation procedure. Each litter represents a cluster and it is recognized (e.g. Williams, 1975; Ryan 1992) that littermates do not respond independently. However, George and Bowman (1995) assume exchangeability of littermate response, in which case the exchangeable binary model is appropriate. The binary endpoint under consideration in this analysis is fetal death. Table 1 gives the frequency distribution of litter size and number of deaths within each litter. The model considered for 2 is
~i)= exp(-(Ol ÷
(11)
02di) * (k)a'+az*di),
where di is the dose level in the ith dose group, 01 + 02di>~O, and al +a2di>>-Ofor i = 1..... 5. This model includes a dose term in both the mean response covariate vector and the correlation covariate vector. This allows for a test of dose effect on both effects as follows. If 02 = 0 then there is no dose effect on mean response rate. If 02 = 0 and a2 is non zero then there is an effect of dose on correlation only. It is also possible with this model to test for independence within each dose group by testing H0:al + a2di = 1 for i = 1. . . . . 5. The simplex estimation procedure described in Section 4 resulted in the following estimates of the parameters in Eq. (11 ), &~= (3.1947, -6.0114, 1.0038, -0.5247) with covariance matrix
A
Cov(~) =
0.1755
-0.0874
0.1167
-0.0061
-0.0874
0.5243
0.1439
-0.1408 /
0.1167
0.1439
0.1749
-0.0607 /
-0.0061
-0.1408
-0.0607
0.0461 /
for o~t=(Ol,O2,al,a2). Corresponding z-scores for the parameter estimates are 7.62,-8.30,2.40 and -2.44, respectively. Note that all estimates are significant indicating a dose effect on both mean and correlation. In Table 2, the corresponding estimates of mean response, 2~ and correlation, P2 are given along with observed response rates. These estimates are compared with the maximum likelihood estimates obtained using the non-parametric estimation procedure of Bowman and George (1995). No dose response function is assumed for the non-parametric estimates and these estimates of mean response in Table 2 follow the observed response rates closely, including the downturn in the second dose group following control. An inherent assumption in the Weibull-type model is that response increases with dose. This accounts for the modest differences between the Weibull-type model and the non-parametric procedure in estimates of mean response. The non-parametric estimates for P2 differ from those obtained using the Weibull model as well. The standard errors, given in parenthesis, of the estimates of correlation for both the Weibull model and the non-parametric estimates are larger than the estimated correlations in all dose groups including control. Fig. 1 shows a graph of the estimated dose response curve for the Weibull type model in Eq. (11). The Weibull-model seems to provide reasonable fit to the observed data, shown on Fig. 1 as the discrete dots. A goodness-of-fit test for the Weibull model resulted in a chi-square value of 122.99 with 1 degree of freedom, indicating no lack of fit.
D. Bowman I Statistics & Probability Letters 41 (1999) 1-7 Table 1 Frequency distribution of female mice according to littersize and the number of dead fetuses following exposure to diethylene glycol dimethyl ether (DYME) Price et al. (1987) Dose (g/kg)
No. deaths
Littersize 3
0
7
8
9
10
11
0
1
1
1
1
1
2 3 0.0625
1 1
1
1
2 0.125
0
1
1
15
16
17
4
2
2
1
1
l
2
2
1 1
1
3
1
2
2
1
1
2
1
2
1
2
3
1
2 1
1
1
2
1
1
2 1
2 3 4 7 8 1 2 4 5 6 7 9 11 13
1
1
2
0 1
0.5
14
1
2 3 0.250
13
1
0 1
12
2
1
3
2
1
1 1
1
1 1
2
1 1
1 1 1 1 1 1
1 1
1
2
1 1 1
1 2 1 1 1 1
1 1
2 1
A likelihood ratio test is constructed as an overall test for independence in the data. A Weibull-type model with 2 k = e x p ( - ( 0 1 + 0 2 d i ) * k 1)
(12)
is fit to the data. The model in Eq. (12) is appropriate if the data are independent. In the D Y M E data a likelihood ratio test based on the models in Eqs. (11 ) and (12) was constructed to test for overall independence in the data. The test statistic was 12.097 with 2 degrees o f freedom, resulting in a p - v a l u e o f 0.0024. Thus the L R T for i n d e p e n d e n c e is rejected in the D Y M E data. The estimated values o f al + azdi for the four dose groups were all less than 1, indicating positive correlation in these dose groups. In the control group the estimate was 1.004 which m a y indicate independence a m o n g responses there. A test for independence versus positive correlation within each dose group was derived. Specifically, the null hypothesis, H0 : al + a2 * di = 1 versus H1 : al + a2 * di < 1 is tested for each dose group. The z-scores for each dose group, respectively, were found to be 0 . 0 0 9 , - 0 . 0 6 8 , - 0 . 1 4 4 , - 0 . 2 9 2 , a n d - 0 . 5 7 2 4 , n o n e o f which were significant. For this data, the hypothesis o f independence is not rejected in any dose group, although
6
D. Bowman / Statistics & Probability Letters 41 (1999) 1 - 7
Table 2 Estimates of mean response and correlation in the DYME data i
Observed rates
1 2 3 4 5
0.0473 0.0744 0.0705 0.1271 0.5053
Weibull model
0.0410 0.0590 0.0832 0.1567 0.4689
Non-parametric
(0.0172) (0.0256) (0.0375) (0.0769) (0.2704)
0.0488 0.0877 0.0705 0.1202 0.5041
(0.7468) (0.5820) (0.4397) (0.2731) (1.2234)
0.0179 -0.0564 -0.0124 0.1125 0.1418
p2
(0.0142) (0.0188) (0.0135) (0.0287) (0.0481)
P2
-0.0007 0.0072 0.0196 0.0600 0.2001
(0.0950) (0.1689) (0.0832) (0.1335) (0.3433)
/
ml(d Yi
I --
cl, xi I
model I observed
Fig. 1.
the overall test rejected independence. These results correspond with the estimates of intralitter correlation, P2 given in Table 2 which clearly tend to increase with increasing dose although none are significantly different from zero. This indicates the presence of a dose effect on correlation as was shown by the significance of the estimate of a2. Other likelihood ratio tests may be constructed to test for dose effects. An overall test for dose effect is constructed by comparing the likelihood of a model of the form 2k z e -O~*(k)~
for all dose groups with the full model specified as Eq. (11 ). The results of this test for the DYME data gave a test statistic of 108.38 with 2 degrees of freedom. The null hypothesis of no dose effect is overwhelmingly rejected.
D. Bowman / Statistics & Probability Letters 41 (1999) 1 - 7
7
5. Discussion The model developed here for exchangeable binary observations has the desirable property of allowing for a parametric test for independence. A likelihood ratio test for overall independence across all dose groups was compared with an independence test within each dose groups. Using this method it was shown that for the DYME, an overall test of independence across all dose groups rejected the independence hypothesis, although tests performed separately within each dose group failed to reject the hypothesis of independence in any dose group. This effect may be caused by the significance of ~2 in the first model indicating correlation between littermates is influenced by dose. References Bowman, D., George, E.O., 1995. A saturated model for analyzing exchangeable binary data. J. Amer. Statist. Assoc. 90, 871-879. de Finetti, B., 1974. Theory of Probability. Wiley, New York.
George, E.O., Bowman, D., 1995. A full likelihood procedure for analyzing exchangeable binary data. Biometrics 51, 512-523. George, E.O, Kodell, R., 1996. An independence test for exchangeable data, submitted. Nelder, J.A., Mead, R., 1965. A simplex method for function minimization. Comput. J. 7, 308-313. Price, C.J., Kimmel, C.A., George, J.D., Mart, M.C., 1987, The developmental toxicity of diethylene glycol dimethyl ether in mice. Fundam. Appl. Toxicol. 8, 115-126. Ryan, L., 1992. Quantitative risk assessment for developmental toxicity. Biometrics 48, 163-174. Williams, D.A., 1975. The analysis of binary responses from toxicological experiments involving reproduction and teratogenieity. Biometrics 31, 949-952.