Article Text

## Abstract

**Background and objectives** We aimed to estimate credibility intervals for the British occupational cancer burden to account for bias uncertainty, using a method adapted from Greenland’s Monte Carlo sensitivity analysis.

**Methods** The attributable fraction (AF) methodology used for our cancer burden estimates requires risk estimates and population proportions exposed for each agent/cancer pair. Sources of bias operating on AF estimator components include non-portability of risk estimates, inadequate models, inaccurate data including unknown cancer latency and employment turnover and compromises in using the available estimators. Each source of bias operates on a component of the AF estimator. Independent prior distributions were estimated for each bias, or graphical sensitivity analysis was used to identify plausible distribution ranges for the component variables, with AF recalculated following Monte Carlo repeated sampling from these distributions. The methods are illustrated using the example of lung cancer due to occupational exposure to respirable crystalline silica in men.

**Results** Results are presented graphically for a hierarchy of biases contributing to an overall credibility interval for lung cancer and respirable crystalline silica exposure. An overall credibility interval of 2.0% to 16.2% was estimated for an AF of 3.9% in men. Choice of relative risk and employment turnover were shown to contribute most to overall estimate uncertainty. Bias from using an incorrect estimator makes a much lower contribution.

**Conclusions** The method illustrates the use of credibility intervals to indicate relative contributions of important sources of uncertainty and identifies important data gaps; results depend greatly on the priors chosen.

- Attributable fraction
- bias
- MCSA
- occupational cancer burden
- uncertainty

## Statistics from Altmetric.com

### What this paper adds

Basing CIs on random variability only is insufficient to represent the real uncertainty around burden of disease estimates.

We have used Monte Carlo sensitivity analysis to estimate credibility intervals that take account of bias as well as statistical uncertainty for the UK burden of occupational cancer study.

A new graphical presentation of credibility intervals enables the relative importance of possible sources of bias to be illustrated in a hierarchical structure.

Choice of relative risk and employment turnover were shown to contribute most to overall estimate uncertainty; biases from using an incorrect estimator make much lower contributions.

## Estimating the burden of occupational cancer: assessing bias and uncertainty

Estimates of the burden of disease, which extrapolate from specific data sources and epidemiological studies to population-level measures, are subject to a range of uncertainty and bias. These come from multiple data sources and value choices arising from both the input data and the methods used for modelling and extrapolation.1 2 Sensitivity analysis to investigate this is therefore an important part of burden of disease estimation and allows identification and quantification of the sources of uncertainty that have the largest impact on the results.

Uncertainty may be due to random error, which represents a lack of precision, or bias, which occurs when the variable to be estimated does not have an expectation equal to the (unknown) true population value.

In Britain (GB), we have estimated the current burden of cancer attributable to occupation for all carcinogenic agents and occupations classified by the International Agency for Research on Cancer as group 1 (established) or group 2A (probable) carcinogens.3 In this paper, we present our approach to investigating the uncertainties in the occupationally related cancer burden by establishing credibility intervals for these estimates using multiple bias analysis. We describe two complementary methods: (1) graphical sensitivity analysis (GSA), which systematically displays the effect of varying the parameters contributing to the estimates, and (2) estimation of credibility intervals using Monte Carlo sensitivity analysis (MCSA), in which bias parameters (based on selected priors) are sampled and the bias model is then inverted to provide a distribution of bias-corrected estimates.4 5 We illustrate these methods for several key biases contributing to the estimation of burden using the example of lung cancer and exposure to respirable crystalline silica (RCS).

## Background

The methodology for estimating the current burden of occupational cancer is described elsewhere.6 Briefly, as risk estimates mostly from industry-based studies were used along with national estimates of exposed numbers, Levin’s formula was used to estimate the attributable fraction (AF), that is, the proportion of cases due to occupational exposure. This requires (i) an estimate of the risk of disease such as relative risk (RR), obtained from published literature, and (ii) the proportion of the population ever exposed over a risk exposure period (REP), defined as the relevant exposure period accounting for cancer latency (from 10 to 50 years before the appearance of a cancer for solid tumours, 0–20 years for blood cancers).

UK data sources were used to provide a point estimate of numbers exposed to each carcinogen, with CAREX7 giving numbers by industry sector and the Labour Force Survey (LFS) giving numbers for occupational circumstances such as ‘work as a welder’ and ‘work as a painter’.8 Numbers ever exposed over the REP were then calculated by extrapolating over the REP accounting for employment turnover (turnover factors were estimated using length of employment data from the LFS) and life expectancy and adjusted for employment trends over the period.6 The proportion ever exposed over the REP was estimated by dividing the numbers ever exposed by the total population who had ever been of working age over the same period.

As exposure–response risk estimates and proportions exposed at different levels are not widely available, RRs were selected for exposure levels identified as generally ‘high’ and generally ‘low’ with industries allocated accordingly. Estimated AFs were applied to current deaths and registrations to give attributable cancer numbers.

### Key sources of bias in AF estimation

Table 1 gives the equations used to estimate the AF and attributable cancer numbers (in column 2). Column 3 gives the key sources of bias and uncertainty for each equation, and column 4 gives the model parameters used to represent these biases in MCSA. The biases may arise from (i) inadequate or inaccurate data, which may affect the estimates directly, for example, the point estimate of numbers exposed, or indirectly through data-based parameter values such as employment turnover or the length of the REP; (ii) inherent problems in the risk estimate study designs, including general epidemiological issues such as unmatched exposures and/or confounders, exposure misclassification and publication bias and specific occupational study issues such as healthy worker bias; or (iii) statistical and modelling issues, for example, Levin’s equation used with risk estimates adjusted for confounding, or AFs not estimated separately by age-group (addressed elsewhere9). To identify appropriate risk estimates, we sought published studies that were either of British workers or of populations that represented a reasonably good match and therefore were thought to be ‘portable’ to the British situation, in terms of the industrial processes, technology, time period and distribution of exposures and relevant confounders. Differences regarding these issues between the study chosen and the British situation, plus exposure misclassification and a healthy worker population effect, may have occurred in the published study. There was a notable absence of information on the risks of exposure to occupational carcinogens for women and also on cancer incidence risk so that very often risk estimates for men/mortality were used. Bias associated with the ages of workforce recruitment and retirement have been described elsewhere,6 and timing of the REP is addressed in this paper using GSA only.

Point estimates were available for numbers of workers exposed to specific carcinogens from CAREX and from LFS data for occupations; for these occupations, estimates of the proportion exposed p_{wk}(E) (equation 1) were not available. The number ever exposed over the REP was estimated using equation 2 accounting for staff turnover; this equation is modified in practice to account for the timing of the REP with respect to the target year, which determines the survival of workers to that year, as estimated from life expectancy data.6 Equation 3 gives the proportion ever exposed over the REP. Levin’s equation (equation 4) uses this with RR to estimate the AF. Estimated AFs applied to current deaths and registrations give attributable cancer numbers (equation 5).

## Methods

### Graphical sensitivity analysis

GSA illustrates the impact that different values of data and parameters have on final estimates, which is useful for identifying or testing a range of plausible values for priors of MCSA and to indicate sensitivity to biases.

As a starting point for identifying ranges of plausible values for the main AF variables, we identified all the values used in the British occupational cancer burden study across all the cancer–carcinogen combinations for CAREX and LFS estimates of numbers exposed (n_{0}), the proportion ever exposed over the REP (p(E)), the RRs and the AF.10 The maximum values for the variables n_{0} and therefore p(E) occurred when estimating the burden of non-melanoma skin cancer for mineral oils, and the minimum data points were found for burden estimation of liver cancer from vinyl chloride monomer exposure. Online supplementary figure 1 shows how varying annual staff turnover from 5% to 25% a year and cancer latency from 10 to 60 years affects the proportion of the population ever exposed over the REP. Online supplementary figure 2 shows the effect of the resulting change on the proportion ever exposed on the AFs for RRs up to 5; these ranges were used in the credibility estimates presented below.

### Supplementary Material

### Monte Carlo sensitivity analysis

To obtain an MCSA credibility interval for each AF, a prior distribution source was identified for each of the variables and parameters in the estimation equations (table 1); for example, a range of plausible values for variables used for estimating the proportion ever exposed was explored as above using GSA, or an alternative RR was chosen from another published source. These prior estimates were used to estimate a bias parameter (column 4, table 1) for each bias source. Equations for AF_{i}(i indicates iterations) were built by re-writing equations 2 and 4 in table 1 to recalculate AF, incorporating the bias parameters as multipliers of the equation variables (see online supplementary statistical appendix 1). For example, to estimate a credibility interval for an RR adjusted for confounding rather than an unadjusted RR, equation 4 is re-written as

where RR_{ai} is the RR adjusted for confounding (RR_{ai} varies randomly in i iterations), and the priors C_{ui}=RR_{ai}/RR_{ui} and C_{ri} = RR_{ai}/RR_{t} (RR_{ui}=unknown unadjusted RR, RR_{t}=unknown true target RR) represent the biases due to using an adjusted estimate of RR and additional confounding bias, respectively.

The calculation of AF was replicated 10 000 times using random values from suitable distributions for each bias parameter (see table 2 for those used for the example of lung cancer from exposure to RCS). The 97.5 and 2.5 percentiles give upper and lower 95% credibility interval limits. The influence of each source of bias can also be estimated separately and illustrated relative to other sources graphically.

## Effect of biases on the AF

### Biases due to inadequate or inaccurate data

#### Proportion exposed

Equations 1–3 in table 1 are the proportion ever exposed over the REP and rely on the adequacy of national data sources used. CAREX gives numbers exposed to a limited list of carcinogens at a single point in time (1990–1993) and only by broad industry. For occupational circumstances, all those recorded by occupation from the LFS have been assumed to be exposed. These sources do not allow accurate estimation of the proportions exposed at different levels of exposure. Also, little information is available on latency of cancers associated with occupational carcinogens, leading to uncertainty as to the most appropriate REP; if exposure agents are cancer promoters as opposed to cancer initiators, as assumed here, or there is recovery after exposure ceases, the REP will be more recent.

Online supplementary figure 1 shows the effect of different latencies for solid tumours and staff turnover on proportion ever exposed over the REP. This proportion, p(E), increases approximately in proportion to the increase in staff turnover and to the length of the REP from 10–20 to 10–60 years. For example, for an annual staff turnover of 10%, p(E) increases by about 0.4% per additional 10 years added to the REP for a point estimate of numbers exposed of 100 000 (n_{0} in equation 2). If the true REP is 10–40 years instead of 10–50 years, the latter overestimates p(E) by 19% (54% for a 10–30 REP; 125% for 10–20 years). These percentage errors are constant across the prevalence estimates (n_{0}) but increase for higher staff turnover estimates (to 23%, 69% and 179%, respectively, for 20% TO). This impacts substantially on the estimate of AF, but as the same latency and staff turnover estimates have been used across all cancer sites and carcinogens, the bias affects all AF estimates equally. There is little difference between estimates of proportions exposed when shorter, and therefore more recent latencies are assumed (results not shown).6

### Risk estimate

As discussed earlier, there can potentially be several different sources of bias and uncertainty that may occur when selecting a risk estimate to use in Levin’s equation (equation 4). These include non-portability between the studies used for the risk estimate and the British situation (ie, dissimilarity between the published study populations and the GB situation as regards the nature and levels of exposures, the study exposure period and the distribution of known confounders11), differences in exposure scenarios and levels and definitions of exposed and non-exposed. Good quality studies for some carcinogens, particularly with dose–response data and estimates of low-level exposure risks, can be scarce. Our allocation of industry sectors to higher or lower exposure categories may have caused over- or under-estimation. Where no risk estimate could be identified for a low exposure for specific carcinogens, we took the harmonic mean of the ratios of ‘high’ to ‘low’ RR estimates for all cancer–exposure pairs for which data were available and applied this average ratio to the high-level estimate of the carcinogens to obtain a low-level RR. If these were less than 1, RR was set to 1. If such a ratio is applied to a particularly high-level risk estimate, over-estimation of the low-level risk could have occurred; more often however (for 15 of the 21 harmonic mean estimates), this resulted in ‘low’ exposed groups being allocated a zero excess risk. In fact, we allocated a zero excess risk (ie, RR=1) to the low exposure category for 23 out of the 42 cancer/carcinogen pairs for which low exposure estimation was carried out, eight from published data, a potentially major source of underestimation. Online supplementary figure 2 demonstrates the effect on the AF of changes in RR, for the range of population proportions exposed and RRs encountered in this study.12 Doubling a small excess risk of 0.1 for low exposure nearly doubles the AF; a similar absolute reduction in a high exposed risk estimate (eg, RR=2 to 1.9) reduces AF by about 10%.

### Biases due to choice of model equation

#### Levin’s equation for the AF

Bias is known to result from using Levin’s equation in the presence of confounding.13 The magnitude of bias resulting from using an adjusted estimate of RR (RR_{a}) rather than an unadjusted one (RR_{u}) in Levin’s equation14 15 increases in proportion to the departure of the ratio of RR_{a}/RR_{u} from 1.16 This is illustrated in online supplementary figure 3 for a range of values of RR_{a} and ratios of RR_{a}/RR_{u} and two contrasting values of the proportions exposed. The ‘biased’ AF is calculated using equation 9 in online supplementary statistical sppendix 2 in the online Supplementary material with RR_{a} substituted for RR_{u}, and the ‘corrected’ (unbiased) AF is calculated using equation 10. When the adjusted RR is less than the unadjusted one (RR_{a}<RR_{u}), this gives a bias downwards from the true AF and vice versa. The magnitude of bias in absolute terms increases as the adjusted RR increases and as the proportion ever exposed (p(E)) increases (shown in online supplementary figure 3 for the range of values, up to 0.30, presented in our study). However, the magnitude of bias as a percentage of true AF is greatest for lower values of RR and p(E) and rises linearly in proportion to the ratio RR_{a}/RR_{u} as this departs up or down from one .

Figure 1 shows the effect on the AF of (i) changing RRa/RRu from 0.5 (representing confounding due to smoking for lung cancer, see below) through to 1.2 (representing a 20% positive bias due to age and sex confounding), illustrated for (ii) the proportion ever exposed from 0 to 0.20, which represents most of the range of exposed proportions encountered in the GB study, for (iii) the maximum (RR=4.3), minimum (RR=1.01) and midrange (RR=2) relative risks encountered in our GB burden estimates. At most, a bias of −0.13 to +0.035 for AF results from choosing the prior range modelled in figure 1.

## Example: lung cancer from exposure to RCS in men

Exposure to RCS occurs in industries such as pottery, foundries, concrete, brick making, stonework, mining and construction. For our published estimate of AF for lung cancer and RCS, we used an overall meta-RR of 1.32 (95% CI 1.24 to 1.41) from Kurihara and Wada17 for a high exposed group and an RR of 1.0 (0.85–1.3)18 for a low exposed group.19 A point estimate of numbers exposed to RCS was obtained from CAREX with industries such as those listed above allocated to the high category, accounting for 97% of an estimated 550 000 exposed in 1990–1993. From these data, it was estimated that 12.7% of men had been exposed to RCS at a high level at some point during a 40 year REP, with an additional 0.3% exposed at a low level; the AF for lung cancer due to RCS was 3.9% for men.

### Bias affecting the equation variables

Table 2 gives the priors used to estimate bias parameters for estimating credibility intervals for RCS and lung cancer. The priors chosen for the variables contributing to p(E) representing the length of the REP and staff turnover are based on a full range of possible values tested for plausibility using GSA (see online supplementary figure 1), but a narrower range of values (±20%) was chosen to represent our uncertainty about the CAREX exposed number estimates, with the effect tested using GSA.

Figure 2 illustrates how the data priors were screened graphically for plausibility and shows the effect on the proportion ever exposed, p(E) and therefore on the AF of using the range of 20% around the point estimate of numbers of exposed workers from CAREX.

To assess the contribution of uncertainty regarding the risk estimates, we have used three priors to cover the biases due to (i) choosing an RR from a study not well matched to the GB situation as regards confounders (bias source 7 in table 1 represented by C_{r}) or (ii) not well matched as regards the pattern of exposure levels (8 in table 1 represented by C_{m}) and (iii) modelling issues (9, C_{u}), respectively. Principally, to address C_{m}, we first split our original high exposure group so that high-level exposure was restricted to only the construction and pottery industries, and the remaining industry sectors were redefined as a medium exposure group. The original low exposure group of industries remained the same.

To address C_{r}, we used a higher risk estimate of 1.66 (95% CI 1.14 to 2.41) adjusted for smoking from a case–control study nested in the UK potteries cohort20 as a prior for the pottery/construction high exposure group to account for possible non-portability issues regarding confounders in the meta-analysis of Kurihara and Wada (original RR 1.32). To address C_{m}, a lower risk estimate of RR=1.17 (95% CI 1.12 to 1.22)21 was used as a prior for the newly defined medium exposure group; this group now accounts for 17.7% of the CAREX exposed in 1990–1993 (1.9% of men ever exposed during the 40 year period of the REP (p(E)=0.019)). Although we have used multiple priors for RR, as the AF is estimated separately for each exposure level, and these RR priors are specific to one exposure level, these ‘partial’ credibility intervals are estimated independently of one another. The effect on AF of using these two priors for the reallocated industry sectors versus the original RR (1.32) is illustrated in figure 2.

Using RR=1.3217 for the redefined medium exposed group has a greater impact on the AF (up 31%–62% compared with using RR=1.17 for this group, measured on the AF% axis in figure 2) than the (arbitrary) effect of changing the point estimate of numbers exposed (which changes AF by between ±19% for either risk estimate). Similarly, using RR=1.32 rather than the higher RR=1.66 prior for construction/potteries only also has a greater impact on the AF for this group (from 139% below to 43% above the original estimate).

### Levin’s equation for the AF

For our original AF estimate for RCS and lung cancer, the RR we selected (RR=1.32) was adjusted for age, sex, calendar period and smoking.17 To account for the bias introduced by using Levin’s equation in this situation, a third prior for the ratio of adjusted to unadjusted RR (C_{u}=RR_{a}/RR_{u}) was chosen using the following information. Kurihara and Wada suggest that the risk of lung cancer in silicotics is doubled among smokers versus non-smokers, with RRs of 4.47 and 2.24, respectively. A smaller bias due to confounding by age and sex has been suggested, for example, 17% compared with the effect of smoking.22 We have therefore assumed that an RR unadjusted for smoking would be double the adjusted estimate and that there would be a 20% increase in the adjusted estimate versus the unadjusted for age/sex confounding or effect modification.

Thus, using a lower limit of RR_{a}/RR_{u}=0.5 to represent confounding due to smoking, for lung cancer, and an upper limit of RR_{a}/RR_{u}=1.2 representing 20% positive bias due to age and sex confounding gives a 16.8% positive bias in the AF estimate for men highly exposed to silica (with RR_H=1.32, and using Greenland’s formulation,4 equation 3 (see also online supplementary statistical appendix 2 in the supplementary material) for estimating an unbiased AF), which is similar to the estimate of Flegal *et al*
22 of 17% bias due to age/sex confounding.

Figure 1 illustrates the bias range for lung cancer and RCS due to using this AF estimator with adjusted estimates of RR (from +0.003 for RR_{a}/RR_{u}=0.5 to −0.025 for RR_{a}/RR_{u}=1.2) using the prior assumptions described above for the range of ratios of RR_{a}/RR_{u}. This represents a 95% credibility interval of 3.6% to 6.4% around the central AF estimate of 3.9%, contributing relatively little to an all bias source credibility interval of 2.0% to 16.2% (figure 3).

### Random error only

Published CIs around the original RRs from Kurihara and Wada17 and Steenland *et al*
18 were used to estimate a random error only CI for the AF of 3.1% to 4.8%.

### Credibility intervals for lung cancer and RCS

The credibility intervals for each of the separate biases described above, including the CI calculated on RR random error alone, are given in table 2 column 5 and illustrated in figure 3. The bias in the RR estimate makes the largest contribution to the overall credibility interval, mainly due to confounder bias. Bias from estimates of numbers ever exposed, mostly attributable to the range of staff turnover estimates tested, makes a smaller contribution, and the model bias (Levin’s estimator with adjusted RRs) is demonstrated to be less important.

## Discussion

We have presented an approach for establishing credibility intervals for burden of disease estimates, for sensitivity analyses on data uncertainties and methodological assumptions and for identifying data gaps; the relative widths of specific credibility intervals aids quantification of the specific source of uncertainty or bias. Our lung cancer and silica example demonstrates the generally increased width of credibility intervals compared with CIs based on random error of the RRs only. We show that uncertainties in the risk estimate were more important than uncertainties in the estimates of proportion exposed and bias due to using Levin’s equation with adjusted risk estimates.

The credibility intervals depend on the choice of priors. We have demonstrated the use of GSA to select plausible data ranges for these, for example, by graphing exposed proportions and then AF against the contributing variables (see online supplementary figures 1 and 2). Choices of annual staff turnover (0.5–2.5) and REPs (10–50 years in length) may have given overly wide credibility intervals; similarly, use of point estimates from CAREX ±20% may also overestimate the bias.

Identification of priors that address all the possible sources of error for RRs can be problematical. In the lung cancer example, we identified RRs for an additional medium exposure level, and we have used a UK study to address portability concerns. Although our original estimates were generally chosen from meta-analyses or high-quality key studies, most of the risk estimates were related to some estimate of cumulative exposure and rarely a dose–response relationship, which we then assign to ‘higher’ and ‘lower’ exposure categories to which we have allocated the CAREX industry groups, requiring implicit assumptions regarding the similarity of durations and intensities of exposure. Our final choices of priors and alternative exposure categories were based on review of available data or eliciting expert opinion.23

There are other sources of bias that have not been addressed in our example of an estimated credibility interval. We did not take account of bias introduced by not estimating AFs by age; an extension to take account of age was used for the global burden of disease estimation of occupational cancer.24 We have also not addressed bias from assumptions about ages of workforce recruitment and retirement, and other issues such as the assumption of a single exposure cancer initiation model, and risk reversible if exposure ceases, and thus a wrongly timed REP. In addition, the use of risk estimates derived from studies of men for women and mortality risk estimates for incidence may have biased the AFs. Epidemiological studies of occupational groups often also result in a ‘healthy worker effect’, that is, a reduced overall risk estimate compared with the general population. The underestimation of the true effect and thus of the resulting cancer burden together with potential misclassification of exposure in epidemiological studies has also not been addressed explicitly.

The choice of distribution function for the priors can also influence the results. Turner *et al*
23 reported the results from prior distributions elicited from several experts, and both Greenland5 and Greenland and Kheifets25 assess sensitivity of results to disturbance in some of the parameters. With large numbers of parameters, this can lead to large numbers of scenarios to explore; however, disagreement over prior distributions may be less important than disagreement over which point value the unknown parameters should take.26–28

We have chosen a single example of a disease/exposure pair to illustrate how a credibility interval can be estimated for each pair contributing to an estimate of national burden derived using Levin’s equations. Additional information from the studies contributing to the selected RRs (including meta-analysis RRs) could have been used to derive a more appropriate CI for this particular AF example, but with loss of generalisation of the method to other disease/exposure pairs. If all occupational exposures and associated diseases could be identified from a single, very large and comprehensive national population-based study, then CIs could be properly estimated for a combined exposure AF, the ultimate aim of our approach, from the proportions of cases exposed and RRs from within the study, avoiding most of the bias problems described above. However, such a study is rarely available and unlikely to be able to address the rarer occupational exposures and diseases.

One of the most important uses of credibility intervals in disease burden estimation is in providing a tool with which to compare the influence of different sources of bias. For example, it is clear from figure 3 that bias from the incorrect use of Levin’s estimator with adjusted RR estimates contributes less to overall uncertainty than the data-based sources of bias, even for the lung cancer and silica example where adjustment of between −20% (for the possible effect of age-sex confounding) and +50% (for smoking) has been allowed for.

Our example demonstrated that uncertainty in the risk estimates may be more important than data limitations as regards numbers exposed. The most important data gap for occupational AF estimation is probably therefore the exposure level information needed to allocate workers to categories to which appropriate risk estimates can be applied.

## Conclusions

The approach developed for estimating credibility intervals for disease burden estimation can be used to quantify as well as indicate the direction of bias, and the widths of the intervals are useful for summarising and simplifying the presentation of results from sensitivity analyses. In particular, the methods presented here allow direct comparisons to be made between alternative sources of bias as regards their contribution to overall uncertainty of an estimate.

## References

## Footnotes

Contributors SH developed the statistical methodology. LR contributed to the design of the original study and the writing of the paper.

Funding This work was supported by the UK Health and Safety Executive (grant no. JN3117).

Competing interests None declared.

Provenance and peer review Not commissioned; externally peer reviewed.