When considering breast and pubic hair development separately, age at occurrence may be treated as a survival time, perhaps following a univariate normal distribution. To characterize the relationship between ages at the two events, a natural extension is the bivariate normal distribution. Parametric survival analyses allow the estimation of the mean (and standard deviation) age at event. However, the correlation between ages at occurrence of the two events is not easily obtained. Most extensions of survival analysis for multiple outcomes use stratification or competing risk scenarios, in which the correlation is either accounted for but not estimated or only one outcome is observed per individual. A further difficulty is introduced by interval censoring of event times. To estimate the correlation between ages at occurrence of two distinct and observed events, a maximum likelihood approach will be used, with modification to accommodate the interval censored nature of the available data.

A pair of event times, defined as

*X* for the first event and

*Y* for the second, follow some joint probability density distribution

*f*_{X,Y}(

*x*,

*y*). The corresponding cumulative distribution function of the paired variables,

*F*_{X,Y}(

*x*,

*y*) =

*P*(

*X*≤

*x*,

*Y*≤

*y*). If

*X* and

*Y* follow a bivariate normal distribution,

*f*_{X}_{,Y}(

*x*,

*y*) takes the following form (equation

1):

whereas

*F*(

*x*,

*y*) takes the following form (equation

2):

Because there is no closed expression for *F*_{X,Y}(*x*, *y*), the double integrand is approximated with the use of iterative numerical techniques.

In the case of interval censoring,

*X* and

*Y* are not observed directly, but are known to fall within a certain interval; in the case of left and right censoring, the end point of the interval is negative or positive infinity, respectively. These intervals are indicated as (L

_{x}, R

_{x}) for the interval containing x, and (

*L*_{y},

*R*_{y}) for the interval containing

*y*. Given these intervals, a region of support can be obtained which contains the pair of event times (

*x*,

*y*) as follows: [

*F*_{X,Y}(

*R*_{x},

*R*_{y}) –

*F*_{X,Y}(

*R*_{x},

*L*_{y}) –

*F*_{X,Y}(

*L*_{x},

*R*_{y}) +

*F*_{X,Y}(

*L*_{x},

*L*_{y})]

(5–7). Evaluating

*F*_{X}_{,Y} over this region approximates the likelihood contribution for each individual. To obtain a likelihood function for the cohort, the individual likelihood contributions are multiplied together and the resulting function maximized with respect to ρ. To increase efficiency without loss of specificity, an alternative is to minimize the negative of the natural log of the likelihood.

This method was used to estimate the correlation between age at entry into Tanner stage 2 (or greater, if stage 2 is not observed) of breast and of pubic hair development among girls participating in the ALSPAC

(8). The physical development of this cohort has been assessed by the use of annual questionnaires that ascertain Tanner stage of breast and pubic hair development (self-assessed at time of questionnaire completion), mailed to participants from the ages of 8 to 14. Because data were collected annually, entry into breast and pubic hair development is known to occur within intervals defined by timing of questionnaire completion. The median age at entry for each marker and accompanying standard deviation were calculated by the use of the LIFEREG procedure in SAS (SAS Institute, Cary, NC), with interval censoring and the normal distribution specified. This median age is equivalent to the mean, attributable to the normal distribution assumption and the fact that the event occurs within a very specific timeframe during which there are no outliers or skewness. Then, the likelihood procedure was implemented in SAS by use of the PROBBNRM function to estimate the cumulative distribution function of the bivariate normal distribution.

To determine whether the correlation between ages at entry into breast and pubic hair development is a function of some underlying variable, the analysis was repeated for selected maternal and child characteristics that are not intermediate between breast and pubic hair development and are known to be associated with pubertal development. These included mother's prepregnancy body mass index (BMI), birth order, child's BMI at age 8, and child's race. Human subject protection was assessed and approved by the ALSPAC Law and Ethics Committee, the Local Research Ethics Committees, and the Centers for Disease Control and Prevention Institutional Review Board.