METHODS

Previously constructed demographic and exposure history files were obtained for the 2342 members of the diatomaceous earth cohort comprising all employees with at least 1 year of employment, and active anytime during 1942–94 at the diatomaceous earth plant.⁴ The cohort was exclusively white men about 25% of whom were Hispanic in origin. The mean duration of employment was 7.4 years. Follow up was from the date on which 1 year of employment was reached or 1 January 1942, whichever was the later, until the earlier date of 31 December 1994 or the date of death. The work history file provided a complete sequence of job assignments with dates and corresponding estimated exposures for total dust, silica (mainly as cristobalite), and asbestos. Smoking information, available for half of the study population, was used in some analyses as a three level stratification variable (smokers, non-smokers, unknown).

We selected for study all deaths due to LDOC also excluding pneumonia and infectious diseases (international classification of diseases, 9th revision (ICD-9): 490–519 and equivalent ranges in ICD-5 to ICD-8). We made this choice because (a) there were too few deaths from silicosis for analysis (n=10), (b) other respiratory diseases are related to silica,^17,18 and (c) some silicosis is very likely to be misclassified as other respiratory diseases.^19,23 Two of the 67 deaths from LDOC available for analysis were coded as asbestosis. Of four deaths from tuberculosis (not included in LDOC) none were silicotuberculosis.⁴ For analysis of silicosis, 70 cases with positive radiographic evidence of silicosis (International Labour Organisation (ILO): ≥1/0 or greater) were incident in the period of follow up.

Exposure history

The diatomaceous earth mining and mill workers were employed at a California plant in operation since 1902. A detailed retrospective exposure assessment has been conducted for this population.^5,24 Measured total respirable dust concentrations at this factory averaged 3.55 mg/m³ (SD 1.25) before 1949 but declined more than 10-fold after 1973. Crystalline silica content of the diatomaceous earth dusts varied from 1% to 25% (bulk samples) depending on process location. Thus, exposures to crystalline silica were low in mining (quarrying) activities but high in mill work in the form of cristobalite, a product of the high temperature calcining process.²⁴ The mean concentration of respirable crystalline silica averaged over the years of employment of the cohort was 0.29 mg/m³.

Cumulative exposures to respirable crystalline silica dust were estimated for each worker by linking the historical exposure data with detailed work history files.⁴ Cumulative exposures as time weighted sums of exposure concentrations were calculated annually for total and crystalline silica respirable dusts during the period 1942–94 for each worker.²² The mean final cumulative exposures to total respirable dust and respirable crystalline silica dust were 7.31 and 2.16 mg/m³.years respectively (maximums 168.84 and 62.52, respectively).⁴ Estimates of asbestos exposure (fibres/ml) for all workers in the cohort were derived from historical exposure monitoring data for jobs with direct handling of asbestos, production records, and documented quantities of asbestos contained in various mixed products.²⁵ Cumulative exposure to respirable asbestos (fibres/ml.years) was calculated in the same manner as the exposure to silica.

Mortality analyses

A life table program (Epicure: DATAB, version 2.06)²⁶ was used to generate files for Poisson regression modelling of LDOC exposure-response. In creating analytical files, cumulative exposure to silica was classified in 50 levels but was entered into regression models as a continuous variable defined by the mean values of the classification strata. Exposure lag periods from 2 to 15 years were considered.

Multiplicative relative rate and additive excess rate models for mortalities from LDOC were evaluated as in previous analyses of lung cancer.²² The general form for relative rate models was:where, exp(â₀) represents the background rate and E is the cumulative exposure to silica. The following variants of f(E), which is the rate ratio (RR), were evaluated:

where the â, â₁, and â₂ are parameters to be estimated. The additive excess rate model evaluated was of the form:Models with the largest decrease in deviance (decrease in −2 log(likelihood)) with addition of exposure terms were considered the best fitting. Relative or excess rates were modelled with stratification on age (13 levels: <25, 25–29, 30–34, 35–39...,≥80), calendar time (nine levels:1942–54, 1955–59,...1990–94), ethnicity (Hispanic v non-Hispanic), and time since first observed (three levels: <10, 10–19, ≥20 years). Alternatively, in relative rate models, external standardisation on age and calendar time was achieved by using United States mortalities from LDOC during 1940–1994²⁷ as a multiplier of person-years thereby producing models of standardised RRs. This external standardisation, which accommodated a general departure of the workers' risk from United States rates as would be expected with a healthy worker effect, resulted in a considerably better fit than in models with exclusively internal standardisation (stratification).

Analytical files were also constructed for Cox's proportional hazards models with continuous cumulative exposures and risk sets defined on age with a FORTRAN program as described elsewhere.²² Results using Cox's proportional hazard analyses were similar to the Poisson analyses (increasingly so with finer stratification of exposure) and are not shown.

Unlagged models seemed to provide the best fit to the data in Poisson analyses although lagged models performed almost as well. For this reason and on grounds of biological plausibility (recent exposures could contribute to an adverse effect), unlagged analyses are presented. In the case of models with external standardisation, plausibility of intercepts (representing standardised mortality ratios for the unexposed, or baseline risk) was also considered. In addition, models were graphically evaluated by comparing fit to categorical models with five exposure strata, and cubic spline smoothing models. The spline models, providing a flexible basis for detecting smooth departures from a specified model form, were fitted as generalised additive models with three degrees of freedom for the effect of exposure²⁸ using S-Plus, version 4.5.²⁹ A final functional form for modelling the relation between exposure to crystalline silica and LDOC was chosen for estimating excess lifetime risk.

Radiographic analyses

Chest radiographs had been generated for the diatomaceous earth cohort since the 1930s in a surveillance programme.²⁰ Eighty five per cent of workers participated in the programme, 77% among those hired before 1940, but only a small proportion (4%) chose to continue after the end of employment. As described in detail elsewhere,²⁰ a worker's most recent readable film was rated for pneumoconiosis by three B readers according to the ILO 1980 classification system. If two or more readers considered the film positive for opacities (≥1/0 for small opacities or with large opacities) they then reviewed the films immediately before until a negative film was found and recorded the dates of first positive and last negative films. The median dates for first positive film were used (or the mean date if only two readers considered the films positive).²⁰ Poisson regression and proportional hazards analyses were performed as for LDOC with the (median) date of first positive film used to define the time of onset of silicosis. Follow up was ended on the date of the last readable film or the (median) date of the first positive film for silicosis, whichever occurred earlier. Workers who had onset of silicosis or last film before 1942 or within 1 year of the start of employment were excluded (reducing cases of silicosis from 84 to 70). Although the previous investigation²⁰ had included all cases of silicosis identified through surveillance, we limited the study to those occurring in the defined cohort to achieve balance with the mortality analysis. It was assumed throughout, that if a worker's last readable film was negative, then all earlier films were also negative.

Estimation of working lifetime risks

Excess lifetime risk of death from LDOC (excluding pneumonia) was estimated for a range of crystalline silica dust concentrations by an actuarial method developed for a risk analysis of radon.³⁰ The procedure is not specific to any type of disease but depends on the model of risk specified; it accounts for competing risks and assumes that workers were exposed to a constant silica concentration for 45 years (or until death) between the ages of 20 and 65. Annual excess risks were accumulated up to the age of 85. Excess risk occurring among the small proportion of those surviving past the age of 85 was not calculated because rates beyond that age are generally unstable. Rate ratios for mortality from LDOC corresponding to work at various concentrations of silica were derived from the final externally adjusted Poisson linear relative rate model. Age specific competing cause of death rates came from a life table for the United States population.³¹ Excess lifetime risk for radiographic silicosis was calculated with the same actuarial procedure but using the estimated background rate (in unexposed people) for new onset, positive radiographic silicosis derived from the model.

RESULTS

Mortality from LDOC

There were 67 deaths in the diatomaceous earth cohort for which the underlying cause was lung disease other than cancer (excluding pneumonia and infectious diseases). Among Hispanic workers the estimated mortality from LDOC from Poisson regression, based on 11 deaths and without controlling for smoking, was about half that of other workers (data not shown), probably a reflection of lower smoking rates among Hispanic workers and is consistent with cancer rates related to lower smoking in the California Hispanic population.³²

The results from fitting Poisson regression models with continuous or categorical unlagged exposures and controlling for smoking are shown in figure 1. Estimates with categorical exposure concentrations (with different cut off points) were consistent with those reported by the original investigators.⁴ Because there were no deaths from LDOC with cumulative exposures greater than 32 mg/m³.years (maximum 62.5) the exposure response diminished at high cumulative exposures as showed in the log linear spline plot which showed substantial downward deviation above 25 mg/m³.years (fig 1). The improvement in fit of this spline over the log linear model (1a) was highly significant (14.9 reduction in deviance for about two added degrees of freedom, p<0.0001) whereas the improvement over the power model (1d) was not significant (change=2.07, p=0.4). This decline in exposure-response with increasing cumulative exposure, and a similar finding with radiographic silicosis, we think most likely resulted from some form of survivor selection. Based on it, we decided to exclude from analysis observation time with cumulative exposures over 10 mg/m³.years (representing 3% of follow up, 16% of deaths from LDOC (56 v 67)). A cumulative exposure of 10 mg/m³.years corresponds to 45 years of exposure at more than four times the current permissible respirable exposure limit for this form of crystalline silica.

The best fitting model with truncated cumulative exposure was the linear relative rate model (1e) although the power model (1d) performed almost as well among single variable models (table 1). The log linear model (1a) (which fitted very poorly without excluding follow up above 10 mg/m³.years, fig 1), fitted almost as well. For models that used external adjustment for age and calendar time and controlled for smoking, the changes in deviance upon adding the exposure term were 16.49, 16.36, 15.87, and 14.98 respectively in the linear relative rate, power, log-square root, and log linear Poisson models (table 1). The log quadratic model (form 1c) fitted somewhat better (change=17.72) but not enough to justify adding another term. Addition of a variable for the exponent of exposure in the linear relative rate model (shape model, form 1f) produced no improvement (change=0.03). Among the relative rate models, external compared with internal adjustment for age and calendar time provided considerable improvement in fit for the exposure terms, contributing reductions in deviance of about 5.0 (linear relative rate, power, log square root models) and 4.5 (log linear model). As a result, the externally adjusted silica effect estimates were slightly larger (and more precise) than those with internal adjustment. Compared with relative rate models, the additive excess rate model (form 2a) gave a decidedly inferior fit (data not shown).

Only in the poorly fitting log linear relative rate model, when high cumulative exposure follow up was not excluded, did average silica intensity (cumulative exposure divided by duration of exposure to silica at the time of follow up) produce a superior fit compared with cumulative exposure to silica. Duration of exposure to silica alone performed much less well (change=2.4). With the power, square root, and linear relative rate models, cumulative exposure to silica out performed models with average silica concentration or duration terms.

The linear relative rate Poisson model with external standardisation was chosen as the basis for estimating excess risk. For this model, the RR estimate was 4.2 (p<0.0001) at the mean cumulative exposure experienced by the cases with LDOC (5.8 mg/m³.years); among the other better fitting models it ranged from 4.0 to 4.7 (table 1). Extrapolating to the maximum cumulative exposure of the cases with LDOC (31.8 mg/m³.years), the RRs ranged from 16.0 to 38.4. Intercepts for the Poisson externally adjusted log linear, log square root, power, and linear relative rate models (without adjustment for smoking or time since first observed) corresponded to baseline (unexposed) standardised mortality ratios (SMRs) of 1.27, 0.87, 0.93, and 1.02, respectively, for white men.

Morbidity from LDOC: analyses of radiographic opacities

The 70 cases of silicosis incident during 1942–94 had a highly skewed distribution in calendar time, 51 or 73% of cases occurring during the first 13 of 53 years (25%) of follow up (1942–54). In Poisson regression models, radiographic silicosis incidence in the 1942–54 period, controlling for cumulative exposure to silica, was 13.3 times higher than in subsequent years, which was highly significant. The mean cumulative exposures of the 51 incident cases in 1942–54 (7.1 mg/m³.years) was somewhat greater than for the five cases incident in the subsequent 10 years (5.4 mg/m³.years), or those in later decades.

Without excluding high cumulative exposures, the shape model (1f) with Poisson regression provided the best fit (change in −2 ln(likelihood)= 30.01) to the radiographic data, the log linear model provided the worst (change in −log(likelihood)=3.02), and the linear relative rate and power models provided intermediate fits to the data. Deviations from the fitted models were apparent above 10 mg/m³.years with a spline model (fig 2). As with the LDOC outcome, analyses were then restricted to observation with up to 10 mg/m³.years unlagged cumulative exposure (excluding 2.5% of follow up, and 10 (14%), of cases with silicosis, table 2). The best fitting, linear relative rate model predicted an RR of 25.6 (p<0.0001) for cases of silicosis at the mean exposure, based on an estimated background rate of opacities in the very few cases among low or unexposed workers.

In the population at risk both for mortality and for incidence of silicosis (excluding those not participating in radiographic surveillance or who had a positive film before 1942 or before the start of follow up), there were subsequently 38 deaths due to LDOC, 14 of which were incident silicotis. Having an incident positive film (in most cases before the end of employment) was a strong predictor of death from LDOC (unadjusted OR=6.9, p<0.0001, table 3). However, most workers (96%) ended participation in radiographic surveillance after the end of employment. Poisson regression analysis for the combined outcome of a positive x ray film or death from LDOC was carried out with follow up until a first positive film or death even though in most cases, workers were no longer participating in surveillance after leaving employment. The resulting RR estimates were greater than for LDOC alone but smaller than for radiographic silicosis alone, due to the incomplete follow up, and, accordingly are not presented.

Estimates of excess lifetime risk

Estimates of excess lifetime mortality from LDOC resulting from up to 45 years of exposure to silica at concentrations of 0.001–0.20 mg/m³ were calculated for the analyses with and without exclusion of observations above 10 mg/m³.years (table 4). Forty five years of exposure at 0.20 mg/m³ corresponds to a cumulative exposure of 9 mg/m³.years. At an exposure of 0.05 mg/m³, which is roughly the Occupational Safety and Health Administration (OSHA) standard for the cristobalite form of silica, the lifetime excess risk for mortality from LDOC from the model excluding higher cumulative exposures was 54/1000 workers (95% CI 17 to 150); without the exclusion, the excess risk was 50/1000. At higher concentrations, exclusion made more of a difference: at 0.20 mg/m³ the excess lifetime risk estimates were 190/1000 with the exclusion and 140/1000 without it (table 4).

The excess lifetime risk for radiographic silicosis was calculated from a simplified Poisson regression model (controlling confounding by age in three levels and by calendar year in two levels) to estimate the background rate as well as relative rates. At the OSHA limit with analyses excluding higher cumulative exposures, the lifetime excess risk for radiographic silicosis was 75/1000 workers, or 7.5%, and without the exclusion it was 68/1000 (table 4). At lower exposure concentrations—for example, 0.005 mg/m³—exclusion resulted in substantially smaller excess lifetime risk estimates whereas at higher concentrations—for example, 0.20 mg/m³—risks were considerably greater.

DISCUSSION

We have estimated that 54 out of 1000 workers would die from silica related lung diseases other than cancer or pneumonia when exposed to respirable crystalline silica for their working lifetime at the current OSHA standard (about 0.05 mg/m³ for cristobalite dust) or the current NIOSH recommended exposure level (REL, 0.05 mg/m³ for silica). At the OSHA limit for respirable quartz dust (about 0.10 mg/m³), 10% of workers would experience excess mortality from LDOC. This estimate was essentially the same or increased with various models and lag periods. Based on surveillance of x ray films that largely ended with employment, about 75/1000 would develop silicosis at the 0.05 mg/m³ exposure concentration and 140 at the 0.10 concentration. These risks far exceed levels generally considered acceptable by OSHA based on judicial review of the OSHA benzene standard—for example, <1/1000).³³ Furthermore, the risk of silicosis has undoubtedly been underestimated due to the end of surveillance with the end of employment (discussed later) and because radiography has a relatively low detection sensitivity for silicosis compared with necropsy assessment.²¹

For the diatomaceous earth workers themselves, who were exposed to mean concentrations of crystalline silica of about 0.29 mg/m³, there would have been more than 200 lifetime excess deaths from LDOC /1000 workers, or 20%, and the lifetime excess of radiographic cases of silicosis (based on incomplete surveillance) would have been considerably more than 200/1000 workers.

We also found that an initial radiographic finding of ILO ≥1/0 or any large opacities among workers exposed to silica at the end of employment was a strong predictor of subsequent death from LDOC (RR=6.9). However, we were unable with this dataset to evaluate the extent to which an early end of exposure of a worker with a positive x ray film (ILO ≥1/0) would diminish subsequent risk of attributable death from LDOC because the first positive film had not been coded in the ILO classification and because worker participation in surveillance was variable over time and very incomplete after the end of employment.

Today, populations exposed to silica exist not only in traditional industries where continuous exposures occur—for example, mining, tunnelling, ceramics, foundry—but also, for example, in construction where very high exposures can result from dry wall sanding and the casual cutting of materials containing silica with powered cutting tools such as diamond or tungsten carbide bladed saws. Exposures for some workers in these industrial sectors³⁴ and in some specific tasks³⁵ have been measured at concentrations of 10 to 50 times the OSHA PEL of today.

Model choice

Model choice was not critical because several model forms and a range of lag periods gave similar findings. Model choices favoured on deviance grounds in the externally adjusted mortality analyses produced estimates of background (unexposed) mortalities from LDOC (from model intercepts) that were high (SMR 0.93 and 1.02, respectively, from the power and linear relative rate models for non-Hispanic white men). LDOC is one of the outcomes most affected by healthy worker selection. In predominantly white men industrial cohorts largely free of respiratory hazards, SMRs for LDOC typically range from 0.5 to 0.8.³⁶ Finding a background SMR of 0.93 or 1.03 for the best fitting models suggests that some work related LDOC has not been accounted for. It is against this background rate that excess deaths are calculated suggesting that excess lifetime risks have been underestimated. Although not fitting quite as well as the linear relative rate model, the log square root model had an intercept corresponding to an SMR of 0.87 and a larger RR for exposure to silica (4.7 for log square root v 4.2 for linear relative rate models). Basing the excess lifetime risk on the square root model would have appreciably increased the estimates of lifetime risk.

The relative rate model for silicosis assumes a non-zero background rate. A survey of the scientific literature on prevalence of positive radiographs (small opacities ≥1/0) within populations presumably not exposed occupationally to dusts showed a wide range of prevalence³⁷: 0.2%–12.0%. In our analysis, in which there were essentially no unexposed workers, the regression model produced a low (but imprecise) estimate of background, implying that almost all cases were related to exposure to silica from work.

Observing a slower rate of increase (and possibly decrease) of risk of LDOC or silicosis at higher cumulative exposures (square root, power, and spline models) prompts several explanations. Accelerated attrition from the cohort of smokers³⁸ or others with compromised respiration unrelated to workplace exposures is one explanation, although control for time since first observed and smoking (ever, never) would have partially corrected that bias. Other contributing factors might be (a) depletion of high susceptibility subpopulations, (b) exposure misclassification, particularly for the early years where sampling data were lacking, and (c) the inherent nature of relative risk models where both background rates and cumulative exposures increase with age. Restricting observation to lower cumulative exposures (<10 mg/m³.years) that still covered a wide range in excess of what could accrue at currently allowed concentrations of crystalline silica resulted in improved fit and plausible linear effects for both mortality from LDOC and onset of silicosis.

Evaluation of confounding and other issues

Data on smoking, collected since the 1960s in the company's radiographic screening programme, were available for 1171 of the subjects (50%). However, smoking habits were unknown for 45 of the 67 workers that died from LDOC (67%). Our Poisson regression analyses for LDOC, stratified on smoking, have partially rectified the confounding by smoking issue. Furthermore, analyses performed without control for smoking produced slightly smaller and less precise estimates of the effects of silica, suggesting that smoking is a negative confounder. In their analysis of this cohort, Checkoway et al⁴ applied the method of Axelson³⁹ concluding that it was very unlikely that cigarette smoking could account for the association found between mortality from LDOC and cumulative exposure to silica.

A possible confounding role for asbestos was a concern in earlier analyses of this diatomaceous earth cohort. Based on a retrospective exposure assessment for asbestos, Checkoway et al⁴ found no evidence that exposure to asbestos accounted for the observed association between mortality from LDOC and cumulative exposure to silica in the diatomaceous earth cohort. We too found no confounding by asbestos in Poisson regression models and omitted asbestos exposure in our final models.

Multiple cause of death coding, not available in this study, would have increased associations by including additional cases in which silicosis and other respiratory diseases were contributing but not underlying causes of death.

Radiographic surveillance

The apparent higher incidence of radiographic silicosis during 1942–54 after controlling for exposure to silica suggests that surveillance underwent a major expansion some time after 1942, identifying many cases that were present earlier but were undetected. This scenario is consistent with the history of this employer, which, in 1940, denied to the State of California that there had been any cases of silicosis,¹⁹ and which denied the State access to medial records until after a public expose of conditions in 1952.¹⁹ Late detection would identify cases at more advanced stages of disease. Those detected before 1954 did have a somewhat higher cumulative exposure at the time of first positive film. Because the first positive films had not been ILO rated, the late detection hypothesis could not be assessed from severity. Bias from late detection of silicosis would result in underestimation of the exposure-response.

Other estimates of risk from silica

We are not aware of other estimations of excess lifetime mortality from LDOC in populations exposed to silica, however, several investigators have estimated prevalence, cumulative risk, or excess lifetime risk for onset of silicosis in populations primarily exposed to quartz (table 5).^40–46 Kreiss et al,⁴³ Greaves,²³ and Finkelstein⁴⁷ propose that the wide differences in cumulative risk reported from the various studies probably reflect differing periods of follow up after employment. For example, the Muir study that found the lowest risk (1.2% for 40 years at 0.010 mg/m³) had no follow up after employment.⁴⁰ With the linear relative risk model for onset of positive films, we estimate that 45 years of exposure to 0.05 mg/m³ respirable silica (cristobalite) corresponded to an excess lifetime risk (for small opacities ≥1/0) of about 7.5%, also at the low end of reported cumulative risks (table 5). However, only 4% of diatomaceous earth workers continued participation in surveillance after the end of employment (for an average of 3.9 years).

For mortality from LDOC in the diatomaceous earth cohort, follow up after employment was largely complete (91%), but the estimated excess lifetime risk for mortality from LDOC, 5.4%, was far below estimates from studies of the onset of silicosis. The relatively low risks might be explained by (a) use of mortality from LDOC (excluding pneumonia and infectious disease) as an end point, (b) use of underlying as opposed to multiple cause of death coding, and, (c) less than full accounting for the exposure-response by models with higher than plausible background rates (intercepts). These findings do not support a special risk for cristobalite compared with quartz dust, the subject of most studies, although lack of follow up after employment for silicosis makes this conclusion uncertain. At the OSHA PEL for quartz (about 0.10 mg/m³, or twice the cristobalite PEL), the excess lifetime risk for mortality from LDOC and silicosis based on the diatomaceous earth cohort would be 10% and 14%, respectively.

Extrapolation to ambient exposure to silicas

If the linear exposure-response for LDOC extends to concentrations of exposure to silica in the range 0.001–0.005 mg/m³, which includes ambient exposure concentrations, our findings would suggest considerable risks to the population due to air pollution. Although this may be the case for lung cancer,²² for LDOC and silicosis low dose linearity has not been established and the data from this study are inadequate to test that hypothesis. Moreover, the cristobalite form of silica studied here is not generally present in ambient exposures which tend to be quartz.⁴⁸ Detecting population excesses of LDOC at ambient concentrations of silica—for example, predicted here to be about 15 deaths per 1000 (0.005 mg/m³ for 24 hours/day for 45 years) roughly corresponding to an all cause SMR of 1.015, and an LDOC SMR of about 1.25—would require a large epidemiological investigation that would have to characterise in detail various exposures and other risk factors as well as silica.

CONCLUSION

We found the most appropriate and best fitting model for the effects of exposure to crystalline silica on mortality from LDOC to be the linear relative rate model, although the power model with log transformed cumulative exposure to silica performed almost as well. Based on the linear relative rate model, the excess lifetime risk at the current silica standard (about 0.05 mg/m³ for the cristobalite form of silica) was estimated to be over 5%. Even exposures at 0.01 mg/m³ may confer an unacceptable risk according to OSHA precedent (>1/1000). This risk of LDOC adds to that of lung cancer which seems also to be excessive for workers exposed below the current standard to this form of silica.²² This study adds weight to the already persuasive evidence^{6–16,22,23,41–47} that workplace exposure to silica at present limits are associated with risks that are generally considered by OSHA and other agencies to be unacceptably high.